Two decades of non-invasive genetic monitoring of the grey wolves recolonizing the Alps support very limited dog introgression

Potential hybridization between wolves and dogs has fueled the sensitive conservation and political debate underlying the recovery of the grey wolf throughout Europe. Here we provide the first genetic analysis of wolf-dog admixture in an area entirely recolonized, the northwestern Alps. As part of a long-term monitoring program, we performed genetic screening of thousands of non-invasive samples collected in Switzerland and adjacent territories since the return of the wolf in the mid-1990s. We identified a total of 115 individuals, only 2 of them showing significant signs of admixture stemming from past interbreeding with dogs, followed by backcrossing. This low rate of introgression (<2% accounting for all wolves ever detected over 1998–2017) parallels those from other European populations, especially in Western Europe (<7%). Despite potential hybridization with stray dogs, few founders and strong anthropogenic pressures, the genetic integrity of the Alpine population has remained intact throughout the entire recolonization process. In a context of widespread misinformation, this finding should reduce conflicts among the different actors involved and facilitate wolf conservation. Real-time genetic monitoring will be necessary to identify potential hybrids and support an effective management of this emblematic population.

In the recent years, a wide array of studies have comprehensively monitored the rate and modality of wolf-dog hybridization in expanding European populations (summarized in Table 1), which appears to depend on their level of disturbance (see below), in link with the abundance of feral (returned to the wild state) and stray (free-ranging) dogs. Accordingly, admixture remains low in Western Europe (i.e. 0 to 6.5% across the Italian Apennines and Iberia, Table 1), being locally higher in anthropogenic environments where stray dogs are common 19,20 , and where residing packs are introgressed 21 . The figures are expectedly higher in eastern countries (10-14% in Georgia, Bulgaria and nearby Greece), where free-ranging large-size guarding dogs are an issue 22 . All studies, however, have focused on recovering populations that never faced full eradication. In contrast, the risk of crossbreeding with dogs is expected to be significantly higher in newly recolonized areas. This is firstly due to the low number of founders involved: hybridization in the initial stages of colonization could lead to widespread introgression as the population expands. Second, the areas where the wolf formerly vanished are usually human-dominated, with low acceptance by locals. These conditions, and notably the potentially high rate of poaching, may disrupt their social structure, in turn favoring crossbreeding with stray dogs. For instance, in Croatia and Italy, many hybrids are found near human settlements with a high rate of human-caused mortality 19,20 .
The Alpine wolf population remains a major gap in the European wolf-dog hybridization literature. Exterminated in the French and Swiss Alps in the late 19 th century, and in the Italian Alps during the early 20 th century, wolves have been re-establishing across this range from relict Italian Apennine pockets ever since the early 1990s [23][24][25][26] . Their formerly high mitochondrial DNA diversity has now essentially been replaced by a single Control Region (CR) haplotype (thereafter the Italian haplotype), diagnostic of these regions and never found in dogs nor wolves elsewhere 24,[27][28][29][30] . This new population thus makes an ideal system to test whether hybridization with dogs had stronger genetic consequences in recolonized areas. Listed as "endangered" (based on data from 2007 31 ,), the Alpine wolf population is constantly increasing, forming about 65 packs in the last update (2015-2016), mostly located in the Western parts between France and Italy 32 . In addition to the lack of suitable habitat and persecution, potential hybridization with dogs is considered as a major threat 33 . In Switzerland, only three established packs occur, but tens of itinerant individuals migrate throughout the country and neighboring regions. Wolf presence is highly controversial in this country: despite its strictly protected status, the species is under very stringent management through selective removal of stock-raiding individuals (when predation impact is considered too high), but is facing strong pressure by livestock breeders, hunters and local  34 , which conveys public rumors such as that wolf-dog hybrids are widespread. Hence, inferring the rate of wolf-dog introgressive hybridization in the Alpine population is also crucial for its future management and social perception.
In this study, we conducted a genetic survey of all the wolves detected in the Swiss and adjacent territories, as part of an ongoing non-invasive genetic monitoring since the early stages of the recolonization process started more than 20 years ago 24,25 . We genotyped thousands of non-invasive DNA samples and performed admixture analyses of multilocus genotypes of putative wolves and reference dogs, in order to estimate the rate of hybridization and admixture in the newly colonized Alpine population.

Results
Genetic screening. We sequenced the left domain of the hypervariable mitochondrial (mtDNA) Control Region (CR) in a total of 3,463 non-invasive unidentified samples and 23 tissues from dead wolves (359-360 bp or 235-236 bp, see Methods). A total of 1,645 samples could be assigned to the private Italian wolf haplotype (47.2%), 409 to typical dog haplotypes (C. l. familiaris; 11.7%), 709 to red foxes (Vulpes vulpes; 20.3%), 361 to other mammals (10.4%) and 362 (10.4%) could not be either PCR-amplified or sequenced (due to contamination by several taxa, e.g. saliva samples on preys). The 68 reference dogs did not carry the Italian wolf CR haplotype, but only typical dog ones (Table S2).

Population genetic and admixture analyses. The Swiss wolves met Hardy-Weinberg Equilibrium
(HWE) and appear panmictic among the three mainly occupied management compartments defined for Switzerland (n = 73, non-significant pairwise F st ); note that this analysis excluded presumptive cubs from the three packs, known from field observations and confirmed by their genetic profiles. We found similar levels of diversity between these compartments (observed heterozygosity, H e = 0.62-0.64; Allelic richness, A r ~ 2.4, scaled on two individuals). This lack of genetic structure is also supported by the frequent movement of individuals between Swiss cantons (i.e. states; Table S1). The full wolf dataset had a probability of identity (PID) of 2.2 × 10 -8 and a PID sibs of 3.8 × 10 −4 .
Bayesian clustering analyses of the 115 putative wolves and 68 reference dogs (obtained from Switzerland) by STRUCTURE 35 unambiguously discriminated between the two groups (K = 2, ΔK = 1485.0) ( Fig. 2A). Principal Component Analysis (PCA) on individual genotypes also recover the two corresponding clusters along the first axis ( Fig. 3). At least four individuals, however, were not perfectly assigned to their corresponding group, i.e. STRUCTURE admixture coefficients (Q) below 1.0 and intermediate positions on the first axis of the PCA (Fig. 3). One challenge is to interpret whether such pattern arises from recent gene flow or imperfect clustering due to shared polymorphism (i.e. false positive), which is bounded by the number and informativeness of the genetic marker used 19,36 . To address this issue, we simulated different classes of hybrids, using genotypes that were correctly assigned to their clusters ("pure" wolves and dogs, Q > 0.99) as parental references, on which we conducted similar analyses (Fig. S1). F1s, F2s and wolf first generation backcrosses received average wolf ancestry Q w values of 0.48 (95% CI: 0.28-0.72), 0.50 (95% CI: 0.26-0.79), 0.73 (95% CI: 0.50-0.94) respectively, and thus marginally encompass those of pure individuals (95% CI: Q w = 0.90-1.0; Fig. S1). However, the distinction becomes difficult with second-generation backcrosses (average Q w = 0.89, 95% CI: 0.61-0.96; Fig. S1). Following the methodology of previous studies 19, 36 , we used the minimum Q w value of pure individuals (0.86) as an ad-hoc threshold to consider individuals as admixed.
Moreover, Bayesian assignment of empirical genotypes to hybrid classes (F1, F2, first and second backcross generations in each direction) with NewHybrids 37 confirmed that almost all Alpine wolves had the best probabilities to be purebred animals, except for M51 and F16, which shared an equal probability to be first or second generation backcrosses (Fig. 2B). Five additional individuals had cumulated probabilities above 0.1 to be backcrosses, although there were always primarily assigned to the parental wolf class, rather than to any other class.

Discussion
Our long-term genetic monitoring supports very limited wolf-dog hybridization and introgression in the Swiss Alps and surroundings, confirming the genetic integrity of this recently established wolf population throughout time. Over the past two decades, no F1 hybrid was detected and only two wolves (out of 115) featured significant signs of introgression by dogs, likely stemming from backcrossing. There were accordingly assigned as first or second generation wolf backcrosses. We cannot exclude that other individuals descend from older hybridization events, with the little dog ancestry they hypothetically retained confounding with the standard variation shared by the two parental genomes at our limited set of loci. Because these crossbreeding events are very rare, the wolf gene pool is rapidly restored due to gene flow from pure individuals. While our resolution seems sufficient to  17,19 . Here, the two individuals M80 and M83, with Q w < 0.95, were among the very few with > 15% of missing data, and hence unreliable estimates (Table S1 and Methods). The anecdotal rate of admixture in the Alpine wolves ever since their very first establishment (<2% considering all individuals since 1998) is in line with other European populations, although these were mainly studied over shorter timeframes (Table 1). It is also supported from a recent assessment in the Italian Alps, which did not identify any hybrid (www.lifewolfalps.eu). Feral and stray dogs are supposedly absent in Switzerland, and wolves usually avoid the immediate vicinity of human settlements. All dogs detected here are doubtlessly pets, hunting or guarding dogs, the Swiss mountains being extensively hiked and exploited for pastoralism. Although local crossbreeding cannot be ruled out, these rare hybridization events may have rather occurred in the Apennines, the source population of Alpine wolves 24,25 . About a million of free-ranging dogs has been reported over Italy by perhaps outdated estimates (10% of feral individuals 38,39 ). While the overall dog introgression rate of Italian wolves still remains low (<7%, Table 1), some localized areas host introgressed packs 21 . Among them, western Tuscany, at the northwestern edge of the Apennine population, appears to be a hotspot of wolf × dog hybridization 20 , near the corridor connecting the Apennines and the Alps. Alpine wolves thus seems to have avoided gene flow with dogs and admixed individuals, which is remarkable given the low numbers of effective founders that arrived in this human-dominated region (~15 25 ). The high mortality rate in Switzerland (21 out of the 115 individuals analyzed) could have also hampered the rapid establishment of stable social structures, promoting mating with dogs instead (as shown in wolf-coyote hybridization 40 ); the first Swiss pack settled >16 years after the species' arrival. Therefore, and despite these suboptimal initial conditions, the now established Alpine population has so far resisted the risks of hybridization. Note that because we identified putative wolves from non-invasive samples based on maternally transmitted mtDNA (i.e. the diagnostic CR Italian haplotype), we only tested admixture resulting from crossbreeding between female wolves and male dogs. Accordingly, all wolf genetic surveys from Italy and the Alps published to date (featuring thousands of samples) failed to detect the reverse (dog mitochondrial introgression in wolf populations; reviewed in 41 ). Potential reasons for this asymmetry are multiple, e.g. less aggressiveness of female wolves towards dogs (while male wolves prey on them), longer partner-seeking activity by female wolves, continuous physiological availability of male dogs for mating, greater survival chances of pups raised by wolf mothers in the wild ( 12 and references therein). Accordingly, all tissues from dead wolves identified here featured the Italian CR haplotype (n = 23). As populations are reconnecting, other wolf haplotypes have recently been found in the eastern Alps (originating from the Balkans and Central Europe 42,43 ) and even southern Italy 44,45 , but so far have not reached our study area 30,41 . The absence of nuclear differentiation between Swiss wolves also confirms that the recolonization of the Alpine arc stems from a single source, as previously suggested 24,25 . Moreover, it indicates that the environmental conditions in Switzerland are favorable enough that itinerant and resident wolves remain panmictic, implying a good connectivity, as also suggested by the migration patterns of individuals between cantons (Table S1).  The idea sometimes propagated that Alpine wolves are "hybrids", largely based on unconvincing and unpublished evidence, is thus not supported by our genetic data. This misconception is confusing the public debate between the various actors involved, e.g. livestock breeders, hunters and institutional stakeholders. Wolves and dogs obviously share a common evolutionary history, intertwined by recurrent crossbreeding ever since the earliest stages of domestication, between 35,000 and 11,000 years ago 46 . As a result, ancient dog breeds feature traces of admixture with wolves and, reciprocally, most Eurasian wolf populations acquired dog alleles through past introgression events, which can confound with ancestral polymorphisms maintained across these recently diverged taxa [47][48][49] . While this should not be confused with contemporary hybridization, it conditions why detecting wolf-dog hybrids for management purposes is far from trivial, given the limitation of resources (i.e. number and nature of analyzable genetic loci) imposed by genetic profiling of wild individuals from low-quality non-invasive samples. New methodologies and technologies (e.g. high-throughput microsatellite genotyping 50 ; SNPs 51,52 ) will allow better resolution to detect hybrid backcrosses from pure individuals in future genetic screening.
In Switzerland, the recolonizing wolves have thus retained their genetic integrity. The two confirmed backcrossed individuals are no longer in the country. Female F16 arrived and settled in the Central Alps from June 2014 to February 2017, when it was ultimately poached (Table S1). Male M51 was only detected at nine instances from February to August 2015 in the eastern and southern parts of the country. None of the parents, and by extension their offspring, of the three resident wolf packs feature signs of dog introgression. In order to protect this population, we stress the need to prevent dog vagrancy in the species' expansion range, as well as to legally remove any F1 hybrids as soon as they are detected by real-time genetic screening and/or suspicion from morphological characters. In the presumed absence of dogs, these F1 hybrids, potentially migrating from Italy, are the proximate cause of subsequent dog introgression into the wolf gene pool, and should be the main target of legal regulations. In contrast, managing admixed individuals beyond F1s seems neither relevant nor efficient, since the wolf genome is being restored, and because effectively tracing such level of admixture remains challenging. External criteria are indeed unreliable to identify backcrosses, even from dead animals: despite its dog introgression, no question was raised regarding female F16. Reciprocally, this questions the validity of so-called "reference" wolves identified by morphology; putative wolves should only be selected based on molecular analyses (e.g. 21,36,53,54 , this study). Finally, we emphasize the crucial role played by such genetic monitoring, which provides an empirical basis to move the public debate forward, and hopefully improve the tense relationship between the many actors discussing the fate of large predators claiming back their former ranges.

Methods
Study area and sample collection. The study was conducted in Switzerland and neighboring territories (nearby France and Italy), covering an area of about 45,000 km 2 , as a part of a long-term non-invasive genetic monitoring of the species since its recolonization of the Alpine range in the 1990s. A total of 3,486 samples were analyzed. These include 3,463 unidentified samples, comprising 1,191 scats (34.4%), 2,007 saliva swabs collected on preys (58.0%), 106 hair (3.1%), 119 urine samples (3.4%), 7 regurgitates (0.2%) and 33 blood samples (0.9%), which were non-invasively collected in the field from November 1998 to December 2017 by trained rangers and wardens from Swiss regional wildlife offices and researchers. In addition, tissues from 23 dead wolves (accidentally, illegally or legally killed) were also analyzed. The geographic origin of samples and their number of detections are given in Table S1. Sites for non-invasive sampling were chosen opportunistically based on known or presumed wolf presence, documented prey kills or random direct observations. Tissue, scat and hair samples were stored in 80-95% ethanol at +4 or −20 °C. Swab samples were stored dried and immediately processed after arrival in the laboratory. Regurgitate and urine samples were stored at −20 °C.
No ethics approval was necessary to work with non-invasive samples or tissues from dead animals. Fieldwork procedures were specifically approved by regional wildlife offices and the Federal Office for the Environment (FOEN) as a part of national wolf monitoring activities. DNA from urine samples collected in snow was extracted according to 55 . Species identification was assessed by sequencing a 359-360 bp portion of the left domain of the mtDNA Control Region (CR; primers and methods in 24 ), except for saliva samples for which we used a newly-designed ungulate-unspecific internal primer (HW3: 5′-GCCCTTATTGGACTAAGTG-3′; 235-236 bp amplified) in order to avoid co-amplification of undesirable prey DNA. DNA sequencing was performed on an ABI3100 platform (Applied Biosystems), after purification with the QIAquick PCR Purification kit (Qiagen) or Wizard ® SV Gel and PCR Clean-Up System (Promega). All DNA extractions and pre-PCR setups occurred in a physically-separated laboratory exclusively devoted to the analysis of low copy number DNA samples. Negative controls were employed during all extraction and amplification experiments to monitor contaminations.
Eleven microsatellite (i.e. STR) loci (FH2054, FH2140, FH2161, FH2096, FH2137, FH2088, FH2001, FH2010 56 , PEZ17, PEZ1 57 , CPH5 58 ) and a Y-chromosome sex marker (Sryw f: 5′-GCCGAGTCCTCTCCTGTA-3′/Sryw r: 5′-TTGTATGAACCATCATTGTGA-3′) were amplified, following the multiplex preamplification method 59 , with some modifications. The principle of this two-step PCR approach is that an initial large-volume PCR is carried out including primers for all 12 markers to be genotyped. A post-amplification aliquot from this PCR is then used as template in separate amplifications for each marker to genotype individuals. This procedure was repeated eight times (four times for tissues) to obtain a consensus genotype, following the multi-tube approach 60  The second PCR step was performed in 20 µl reactions containing 2.5-5.0 mM MgCl 2 , 0.1 mM of each dNTP, 0.2-0.5 µmol of each primer, 0.2 mg/ml BSA, 1 x PCR Buffer and 1 U of AmpliTaq ® Gold DNA polymerase (Applied Biosystems), and 2.0 µl of extract. We used an annealing temperature of 50-57 °C and 40 PCR cycles. One primer of each pair was synthesized with a 5′-end fluorescent dye (ATTO532, ATTO550, Dyomics 630, FAM and HEX) to allow detection and sizing of fragments on an ABI Prism 3100 DNA sequencer (Applied Biosystems). Alleles were scored using GeneMapper 4.0 (Applied Biosystems). Consensus genotypes, allelic dropout and rates of false alleles were computed manually as follows: heterozygous at one locus were accepted if both alleles appeared at least four times among the eight replicates, and homozygous were accepted if the genotype was observed in five independent PCRs. If neither of those cases occurred, the alleles were treated as missing data (NA). In addition, we genotyped saliva samples of dogs (n = 68) from 35 different breeds, obtained in Switzerland through veterinary practice, to be used as a reference (Table S2). The final dataset (n = 183) included 131 individuals with no NA (72%), 41 with 0-10% NA (22%), 7 with 10-20% NA (4%), and 4 with 20-45% NA (2%) (Details in Tables S1 and S2).

Standard Genetic Analyses.
We assigned individuals to the management compartments defined by the FOEN for this species (Fig. 1), based on the area they spent most of their time according to our genetic detections. Individuals sampled in neighboring France and Italy where assigned to the closest compartment. For the three main occupied compartments (n = 7, 35 and 31 for compartments III, IV and V, respectively, not-including cubs from packs) we computed observed heterozygosity (H e ), allelic richness (A r ), pairwise F st (with significance tested by 1,000 permutations), and performed test for Hardy-Weinberg Equilibrium (HWE) in FSTAT 61 . Probabilities of identity (PID) were estimated with GenAlEx 6.5 62 .
Admixture analyses. We assigned putative wolves (n = 115) and dogs (n = 68) to groups with the Bayesian clustering algorithm of STRUCTURE 35 . We used the admixture model and ran 10 replicates for K from 1 to 11, each including 100′000 iterations after a burnin of 10′000. Replicates were combined with CLUMP 63 and admixture proportions Q were visualized with DISTRUCT 64 . The ΔK statistic 65 was computed to infer the number of K best explaining the data (STRUCTURE HARVESTER 66 ). In parallel, we conducted a PCA to assess the genetic variance between individuals (ade4 and adegenet R package 67 ).
In order to infer the nature of individuals that featured intermediate STRUCTURE Q values, we used the software NewHybrids 37 to compute their probability of belonging to the six following hybrid classes: F1, F2, wolf first generation backcross, wolf second generation backcross, dog first generation backcross and dog second generation backcross. We pre-assigned pure individuals (Q > 0.99 in the STRUCTURE analyses) as parental references (n = 99 wolves and 55 dogs). The chain was ran for 100′000 iterations, which was sufficient to reach stationarity.
We further implemented a simulation approach to evaluate the power to discriminate between pure and introgressed individuals from our dataset. To this end, we generated F1, F2 and backcrossed hybrids (100 for each category) from the reference individuals selected above, using the hybridize function of adegenet. These simulated and pure genotypes were then analyzed by STRUCTURE (K = 2) and PCA, in the same way as the empirical dataset (see above). This allowed to obtained the distribution of Q w (the wolf ancestry of individual) for each hybrid class, and calculate the threshold below which individuals can be considered admixed without confusion with pure ones 19,36 .