Skip to content


  • Research article
  • Open Access

Alternated selection mechanisms maintain adaptive diversity in different demographic scenarios of a large carnivore

  • 1,
  • 1,
  • 2,
  • 1,
  • 3,
  • 1,
  • 1, 4 and
  • 1, 4Email author
BMC Evolutionary Biology201919:90

  • Received: 15 January 2019
  • Accepted: 4 April 2019
  • Published:



Different population trajectories are expected to impact the signature of neutral and adaptive processes at multiple levels, challenging the assessment of the relative roles of different microevolutionary forces. Here, we integrate adaptive and neutral variability patterns to disentangle how adaptive diversity is driven under different demographic scenarios within the Iberian wolf (Canis lupus) range. We studied the persistent, the expanding and a small, isolated group within the Iberian wolf population, using 3 MHC class II genes (DRB1, DQA1, and DQB1), which diversity was compared with 39 microsatellite loci.


Both the persistent and the expanding groups show evidence of balancing selection, revealed by a significant departure from neutrality at MHC loci, significant higher observed and expected heterozygosity and lower differentiation at MHC than at neutral loci, and signs of positive selection. However, despite exhibiting a significantly higher genetic diversity than the isolated group, the persistent group did not show significant excess of MHC heterozygotes. The expanding group, while showing a similar level of genetic diversity than the persistent group, displays by contrast a significant excess of MHC heterozygotes, which is compatible with the heterozygote advantage mechanism. Results are not clear regarding the role of drift and selection in the isolated group due to the small size of this population. Although diversity indices of MHC loci correspond to neutral expectations in the isolated group, accelerated MHC divergence, revealed by a higher differentiation at MHC than neutral loci, may indicate diversifying selection.


Different selective pressures were observed in the three different demographic scenarios, which are possibly driven by different selection mechanisms to maintain adaptive diversity.


  • Canis lupus
  • Demography
  • Major histocompatibility complex
  • Microsatellites
  • Selection mechanisms


Genetic drift and natural selection are often seen as two counteracting evolutionary forces across populations [1]. In large and stable populations, where the random effects of genetic drift are limited, balancing selection maintains high levels of adaptive variation affecting population fitness and evolutionary potential [2, 3]. By contrast, in small and isolated populations, or within the context of demographic events such as bottlenecks, genetic drift is predominant, masking possible effects of natural selection [4, 5]. Neutral molecular markers allow assessing interactions between demographic history and genetic diversity [6], but do not provide direct information on selective processes [7, 8]. Thus, evidence of selection and adaptive potential is best approached by contrasting neutral and functional molecular markers (e.g., [9, 10]).

In vertebrates, the major histocompatibility complex (MHC) is the prime candidate for pathogen resistance genes and contains some of the most polymorphic functional loci [11, 12]. The exceptional level of MHC polymorphism is believed to be driven by the antagonistic coevolution with pathogens and occurs through pathogen-mediated balancing selection [11, 12]. Beyond its clear significance in modulating pathogen resistance [2, 13, 14], MHC genes have been shown to influence other biological traits such as maternal-fetal interactions, kin recognition, life-time reproductive success and mate choice [1518]. MHC is intimately linked with factors likely to affect individual fitness, population viability and evolutionary potential in changing environments. Thus, patterns of MHC diversity have been repeatedly used in a conservation context in populations of particular interest [1922].

In Europe, abundance and distribution of large carnivores (bear, wolf and lynx) have been dramatically shaped by humans during the past few centuries. Although persecuted and driven to or close to extirpation in several countries, most of the remnant populations have stabilized or even expanded recently [23]. One example is the Iberian wolf (Canis lupus), which erstwhile ranged over most of the Iberian Peninsula, but after decades of severe human persecution became confined to small and fragmented populations, representing about 1/5 of the former population range [2325]. In the 1970s, the Iberian wolf population was estimated to be reduced to ca. 700 individuals, mostly restricted to the northwestern region [23, 26, 27]. After the 1970s, the population started to grow and expanded southwards and eastwards (Fig. 1) [28]. Currently, the Iberian wolf is considered under a stable or positive demographic trend [23]. It occurs in a core area where the species have always persisted in northwestern Iberia [26, 27, 29], and the adjacent re-colonized area (the expansion front), alongside two small and isolated populations [3032]. Such a complex demographic history have translated into a cryptic population structure in the Iberian wolf at small spatial scales, with moderate level of genetic differentiation, including the differentiation of the re-colonized area [33].
Fig. 1
Fig. 1

Map showing location of samples corresponding to three demographic groups (persistent in green, expanding in blue and isolated in yellow). Pie charts represent the relative frequency distribution of seven three-locus MHC haplotypes per demographic group. Each color of the pie chart represents a three-locus MHC haplotype, and size is proportional to haplotype frequency within demographic group (Additional file 2: Table S2). Filled gray polygon represents the estimated wolf distribution in the Iberian Peninsula in 2005 [24]

Remarkable MHC class II variability has been shown in wolf populations [3437]. However, results vary regarding the predominant role of selection over neutral forces. Some studies provided evidence that wolf MHC diversity is maintained by balancing selection, including studies on populations under demographic decline [34, 36], where parasite resistance was suggested as the possible driving force [36]. In contrast, MHC diversity in the bottlenecked wolf population of Scandinavia was shown to be compatible with neutral evolution [38]. Although MHC variability patterns of small and isolated populations may differ from that in large and outbred populations, the power of balancing selection acting on MHC can be outweighed by demographic events, such as bottlenecks and fragmentation with consequent genetic drift [38, 39]. Hence the signature of selection and/or population demographic effects is expected to vary across populations under different demographic histories, challenging the assessment of the relative role of different microevolutionary forces [40].

In this study, we integrate neutral and adaptive diversity patterns to disentangle how adaptive diversity is driven under three different demographic within the Iberian wolf range. Taking into account the genetic population structure of this population based on neutral markers [33], we used the permanent wolf range (i.e., persistent group) as the baseline for demographic stability and compared it to the population dominated by the re-colonization process after the 1970s bottleneck [27] (i.e., expanding group) and to the remaining small and isolated population in the South of Douro River in Portugal (i.e., isolated group). For the persistent group we expect to find natural selection acting to maintain adaptive variation at a stronger intensity than genetic drift; and therefore MHC diversity is expected to exceed neutral diversity. In this scenario, we can either detect balancing or diversifying selection. We may observe balancing selection when MHC exhibits lower population differentiation than neutral loci, due to enhanced effective migration rate and similar allele frequencies in MHC loci [11]. Alternatively, we can detect diversifying selection when MHC exhibits accelerated divergence as a response to spatial or temporal heterogeneity in parasite abundance and diversity [41]. Contrastingly, for expanding and isolated groups, genetic drift is expected to outweigh selection; and therefore similar diversity patterns between MHC and neutral markers are expected.


MHC diversity

We found seven DRB1, four DQA1, and six DQB1 alleles across the 113 Iberian wolf samples (Additional file 1: Table S1), all previously reported in other wolf populations (Additional file 2: Table S2, Additional file 3: Table S3). Diversity at DRB1 and DQB1 was higher than diversity at DQA1 (Additional file 4: Table S4). Observed heterozygosity was higher than expected at all three loci considering the whole dataset, with significant departure from Hardy-Weinberg Equilibrium (Additional file 4: Table S4). Nucleotide and amino acid distances were higher for DRB1 and DQB1 than for DQA1, with amino acid distances higher than the nucleotide distances (Additional file 4: Table S4).

When clustering the three loci (DRB1/DQA1/DQB1), we observed seven haplotypes (Fig. 1, Table 1). Each DRB1 allele is found in haplotypic association with a different DQ pair, with exception of alleles DRB1*05401 and DRB1*05501, which share the same DQ pair, DQA1*00301/ DQB1*00401.
Table 1

Sequence diversity and neutrality tests of the three demographic groups (persistent, expanding and isolated) and the whole Iberian wolf range for the three-locus (DRB1/DQA1/DQB1) MHC haplotypes








Tajima’s D

Fu & Li D*





































Information in the table includes sample size (n), number of haplotypes (h), number of segregating sites (S), nucleotide diversity (π), number of mutations (η), Watterson’s mutation parameter (θW), neutrality tests of Tajima’s D and Fu &Li D*. Statistical significance: *P < 0.05, **P < 0.02

Five haplotypes are shared with other European wolf populations (Additional file 3: Table S3). We cannot confirm the presence of two haplotypes in other European wolf populations due to the absence of haplotype reconstruction in [37], where all alleles that compose these two haplotypes were reported. The three most frequent three-locus haplotypes are shared with Croatia, Italy and Finland (Fig. 2, Additional file 2: Table S2, Additional file 3 Table S3). Iberian haplotypes are distributed throughout the NJ tree (Fig. 2). Seven, six and four three-locus haplotypes were found in the persistent, expanding and isolated groups, respectively (Fig. 1, Table 1). The most frequent haplotype was different for each demographic group (Fig. 1, Additional file 2: Table S2). Sequence diversity as the number of segregating sites, nucleotide diversity, number of mutations and Watterson’s mutation parameter were similar among the three different demographic groups (Table 1).
Fig. 2
Fig. 2

Neighbor joining tree showing relationships between three-locus haplotypes of Iberian and other European wolf populations using p-distances. Colored circles indicate populations where haplotypes were found. Iberian haplotypes are in bold. Haplotype nomenclature is given in Additional file 1: Table S1

Positive selection

All MHC loci revealed positive and significant departure from neutrality (Additional file 4: Table S4), indicating an excess of high-frequency segregating sites. Persistent and expanding groups revealed positive and significant departure from neutrality for Tajima’s D and Fu & Li D* tests (Table 1). The isolated group did not show significant departure from neutrality for Tajima’s D and Fu & Li D* tests (Table 1). The dN/dS ratio was significantly different from neutral expectations in all demographic groups both at all sites and at the peptide binding region for DRB1 loci, and only at all sites for DQB1 loci (Table 2). For DQA1, the dN/dS ratio was significantly different from neutral expectations at all sites in isolated group and at the peptide binding region in persistent and expanding group (Table 2).
Table 2

Z-tests of selection on all sites, peptide binding region (PBR) inferred from [42], and the non-peptide binding region (non-PBR) of the three demographic groups (persistent, expanding and isolated) of Iberian wolves for the three MHC loci (DRB1, DQA1, and DQB1)












































Statistical significance: *P < 0.05, **P < 0.01, ***P < 0.001

Additionally, we found several MHC codons under positive selection in the three demographic groups using OmegaMap, a few having high posterior probability (PP > 95%, Additional file 5: Table S5). DRB1 and DQB1 loci had higher proportions of selected codons than DQA1 (Additional file 5: Table S5). We matched our alignment with the human orthologue [42] and observed that almost all selected codons with high posterior probability (DRB1: 6 in 7, DQA1: 1 in 2, DQB1: 4 in 4 codons) match the same peptide-binding sites described by these authors (Additional file 5: Table S5).

Genetic differentiation between demographic groups

We found a single model with support from the data to explain allelic richness and expected heterozygosity, with a deviance explained of more than 65% (Additional file 6: Table S6). Allelic richness is explained by a combination of demographic group and locus type, while expected heterozygosity is explained by demographic group, both including the effect of locus (Additional file 6: Table S6). For both the number of alleles and observed heterozygosity we found three models with support, presenting a deviance explained between 20 and 43% (Additional file 6: Table S6). Number of alleles is explained by demographic group or a combination of demographic group and locus type (Additional file 6: Table S6). Observed heterozygosity is explained by the null model or by locus type (Additional file 6: Table S6).

According to the selected models, allelic richness, expected heterozygosity and the number of alleles are significantly higher in the persistent than in the isolated group (Additional file 7: Table S7). Number of alleles is also significantly higher in the persistent than in the expanding group (Additional file 7: Table S7). Allelic richness and observed heterozygosity are significantly higher in MHC than in microsatellites.

Comparing MHC and neutral diversity

For both MHC and microsatellite loci, significant deviations from Hardy-Weinberg Equilibrium were observed for the persistent and expanding groups (Table 3). However this deviation is in opposite directions (Table 3). The persistent group showed significant higher expected than observed heterozygosity in both molecular markers, while the expanding group showed significant higher observed than expected heterozygosity in MHC (Table 3). Observed heterozygosity was significantly higher at MHC than at microsatellite loci for the persistent and expanding groups (p < 0.05, t-test, Additional file 8: Table S8). For the expanding group, the mean allelic richness and expected heterozygosity at MHC were also significantly higher than the ones at microsatellite loci (p < 0.05, t-test, Additional file 8: Table S8). For the isolated group, none of the diversity indices was significantly different between MHC and microsatellite loci (p > 0.05, t-test, Additional file 8: Table S8).
Table 3

Genetic diversity of the three demographic groups (persistent, expanding and isolated) and the whole Iberian wolf range for microsatellite loci and for the three-locus (DRB1/DQA1/DQB1) MHC haplotypes





















5.5 (0.3)

3.2 (0.2)


0.578 (0.027)

0.639 (0.028)

0.10* (0.015)


5.7 (0.9)

4.1 (0.7)


0.675 (0.017)

0.742 (0.048)

0.085* (0.039)



4.2 (0.2)

2.8 (0.1)


0.561 (0.033)

0.584 (0.03)

0.04* (0.028)


5.0 (0.6)

3.7 (0.5)


0.813 (0.035)

0.724 (0.034)

−0.124* (0.008)



3.3 (0.2)

2.3 (0.1)


0.546 (0.36)

0.523 (0.024)

−0.05 (0.047)


4.0 (0.6)

2.8 (0.8)


0.733 (0.167)

0.613 (0.089)

−0.160 (0.131)



4.3 (0.2)

2.8 (0.1)






4.9 (0.4)

3.6 (0.3)


0.741 (0.053)

0.693 (0.036)

−0.066 (0.055)

Information in the table includes sample size (n), mean number of alleles per locus (Na), mean number of effective alleles (Ne), mean allelic richness (AR), mean observed and expected heterozygosities (Ho and He), mean inbreeding coefficient (FIS). Standard error is given between brackets. Statistical significance: *P < 0.02

The first two components of PCA based on MHC loci did not differentiate among the three demographic groups (Fig. 3). Pairwise FST values using MHC data support the lack of genetic differentiation between persistent and expanding groups, while moderately significant differentiation is observed between all groups using Jost’s D differentiation (Table 4). Contrastingly, PCA based on microsatellite data showed clear genetic differentiation between all demographic groups, especially for the isolated group, while a partial overlap is observed between persistent and expanding groups (Fig. 3). All FST and Jost’s D pairwise comparisons between groups were significant for genetic differentiation using microsatellite loci (Table 4). Overall, pairwise genetic differentiation using FST was lower in MHC than in microsatellite loci, whereas for Jost’s D differentiation index, the isolated group exhibited higher values for MHC than for microsatellite loci (Table 4). The global co-inertia RV coefficient observed was high (0.98), indicating strong correlation between MHC and microsatellite matrices, but this coefficient is not different from what is expected by chance (p-value = 0.17).
Fig. 3
Fig. 3

Co-inertia analysis (CoA) between neutral and adaptive data for all demographic groups. a. Principal Component Analysis (PCA) for MHC loci; b. PCA for microsatellite loci; dots represent different individuals and circles around dots with different colors represent demographic groups. c. CoA plot indicating the relative position of each demographic group on the factorial plane for the first two CoA eigenvalues. Dots and arrows represent the projected co-ordinates of each dataset (MHC and microsatellite loci, respectively), which are joined by a vector, where the length of the vector is proportional to the divergence between the datasets. d. Canonical weights for MHC alleles. e. Canonical weights for microsatellite alleles

Table 4

Pairwise differentiation values for microsatellite and MHC loci (X/X, respectively), considering three demographic groups of Iberian wolves, estimated for FST and for Jost’s differentiation index (D)







Persistent /Isolated






Statistical significance: *P < 0.01, **P < 0.001


This is the first report of MHC class II diversity in the Iberian wolf. Overall, Iberian wolves exhibit lower MHC diversity than their European counterparts, both in number of alleles and number of three-locus haplotypes (i.e., 6–7 alleles and 7 haplotypes in Iberia vs 6–13 alleles and 13–14 haplotypes in other European populations) [3436]. The exception is the isolated Scandinavian wolf population, which went through a drastic bottleneck recovering from only three individuals [38]. The reduced MHC diversity in the Iberian wolf has also been observed for other genomic regions [31, 4345]. Long-term isolation and past bottlenecks in the Iberian population [27, 45] have certainly reduced neutral diversity and may also be at the origin of this depleted MHC diversity.

All MHC alleles have been previously identified in other canids [37, 46, 47], a pattern that is better explained by the trans-species polymorphism described in MHC canid phylogenies [37, 48] than by hybridization which is not common in the Iberian wolf [49]. The observed 17 alleles were confined in just seven three-locus haplotypes of which five have previously been reported [3436, 48]. The remaining two haplotypes were not unequivocal confirmed to be present in other wolf populations due to the absence of haplotype reconstruction in [37]. The occurrence of conserved three-locus MHC haplotypes in Iberian wolves confirms the tight linkage of these gene loci, supporting a strong selective pressure maintaining haplotype combinations. The association observed between each DRB1 allele and a specific DQ pair supports the preferential association between alleles at these loci [46]. A similar trend is observed in other European wolf population [3436], though not as extreme as in Iberia, where seven DRB1 alleles yielded only seven haplotypes. However, similarly to other European populations, the seven haplotypes are widely distributed in the NJ tree, supporting maximal MHC diversity in the Iberian wolf.

Within Iberia, we found contrasting MHC diversity patterns according to different demographic trajectories. The persistent group exhibited significant higher genetic diversity than the isolated group, but did not show significant excess of MHC heterozygotes. However, we found evidence of balancing selection acting to maintain MHC adaptive variation in this group, compensating the effects of the severe population decline in the 19th and 20th centuries [50, 51]. This is revealed by a significant departure from neutrality at MHC loci, significant higher observed and expected heterozygosity and lower differentiation at MHC than neutral loci, and signs of positive selection with higher proportions of selected codons lying in peptide-binding regions, as inferred from humans. Similar results have been reported for Finnish wolves where, despite a strong genetic signal of population decline [52], the population still reveals a predominant role of balancing selection maintaining MHC diversity [36]. Although selective advantage has been associated with a heterozygote advantage mechanism, MHC diversity can also be promoted by negative frequency-dependent selection [2]. This model proposes that the emergence of a novel pathogen can promote the frequency of a rare or novel allele/haplotype [2, 3, 53]. The frequency shift of an adaptive allele/haplotype by such process can occur in very few generations [53]. We can therefore not exclude that negative frequency-dependent selection may have contributed to maintain high allelic richness possibly compensating for the low number of heterozygotes observed.

The expanding group exhibited similar genetic diversity than the persistent group but show significant excess of MHC heterozygotes. Reduced genetic diversity is expected for marginal populations [10, 5456], but increased heterozygosity at MHC supports a model of heterozygote advantage for this demographic group. The positive and significant result for neutrality statistics provides further evidence for balancing selection acting in this demographic group. The overdominance model postulates that heterozygotes have higher fitness than homozygotes due to the wider spectrum of MHC receptors able to induce parasite resistance [2]. MHC heterozygotes have indeed been associated with resistance to several infections in wolves [36]. Thus, the high MHC genetic diversity exhibited by the expanding group is expected to maintain its high adaptive potential.

The isolated group showed the lowest neutral and adaptive diversity and no significant excess of MHC heterozygotes, a likely result of the dramatic population decline that culminated in isolation in the early twenty-first century [26, 30]. The observed depletion of MHC diversity potentially increases population susceptibility to disease, and adds to current concern for its survival [12, 57]. Indeed, in a wolf serologic survey the isolated group showed no canine parvovirus antibodies, possibly related to high case-fatality rate [58]. No significant departure from neutrality neither significant difference between MHC and neutral loci were observed, suggesting that MHC loci behave according to neutral expectations in the isolated group. Other examples are known for small and isolated populations where genetic drift outweighs the strength of selection [3, 38, 39]. However, one interesting result comparing the isolated group with the others may suggest selective forces acting in this group. The higher Jost’s D differentiation index observed for MHC relative to neutral loci may be explained by diversifying selection as consequence of spatial variation in pathogen-mediated responses [40] and/or due to stronger divergence due to short-term positive selection on beneficial alleles. Interestingly, the isolated wolf population shares several pathogens with wild and domestic species [5961], but holds also new pathogens [60], which could be associated with the observed differentiation in MHC. Inconclusive results are probably due to the small sample size of this group and should be interpreted with caution [62].


By empirically testing adaptive diversity under three different demographic scenarios within the Iberian wolf range, we were able to show that these different demographic scenarios are associated with different selective pressures, which are possibly driven by different selection mechanisms to maintain adaptive diversity. Future studies are warranted to investigate the link between the Iberian wolf immunogenetic diversity and the pathogen community structure to better understand the mechanisms underlying adaptation in each demographic group.


Sampling and DNA extraction

A total of 113 wolf tissue samples were collected across the range of wolves in the Iberian Peninsula in the last two decades. Tissue samples were collected from several Portuguese and Spanish entities (please see Ethics section). To avoid close familiar relationships and the presence of hybrids between wolf and dog in our dataset we selected samples according criteria defined by [33]. Samples were classified into one of three demographic groups: i) persistent group (n = 78), represented by samples from the northwest Iberian Peninsula where the wolf has been continuously present in the last decades, ii) expanding group (n = 25), represented by samples originally from the area where wolves were absent during the recent low of the population in the 1970s, representing a well-defined genetic cluster within Iberia (Castilla y León cluster) [33], and iii) isolated group (n = 10), represented by samples from the small and isolated wolf area south of the Douro river in Portugal (Fig. 1). Samples were stored in ethanol and DNA was further extracted using QIAGEN DNeasy Blood & Tissue Kit.

MHC amplification and sequencing

DLA class II genes DRB1, DQA1 and DQB1 were amplified and sequenced using intronic and locus specific primers described in [34]. Polymerase chain reactions (PCR) were carried out in a total volume of 10 μl containing approximately 10 ng of DNA, 0.4 μM of each primer and 5 μl 2× MyTaq HS Mix (Bioline) (see Additional file 9: Table S9 for PCR conditions). PCR products were sequenced in both directions using BigDye® Terminator v3.1 Cycle Sequencing Kit (Applied Biosystems) according to the manufacturer’s protocol and sequencing products were separated on a 3130xl Genetic Analyzer (Applied Biosystems).

MHC alleles and haplotype assignment

DRB1, DQA1 and DQB1 loci sequences of each of the 113 samples were analyzed and aligned with CodonCode Aligner (CodonCode Corporation, and BioEdit [63].

MHC sequences were phased in DnaSP v. 5.10 [64] using PHASE [65] with the “recombination” model (−MR0) and 1000 iterations after 100 burn-ins. The threshold for accepting true alleles was set to 90%. This methodology has been shown to be accurate for MHC haplotype reconstruction [66], retrieving very low false positives [67]. MHC alleles obtained for each locus were compared with available information in public datasets (NCBI: for the genus Canis. As the three genes studied are chromosomally clustered in tandem, we analyzed each locus individually or as three-locus haplotypes. Three-locus haplotypes (DRB1/DQA1/DQB1) were reconstructed firstly in individuals homozygous for each locus, and then for heterozygous individuals basing inferences on haplotypes observed in homozygous state. Haplotype reconstruction was confirmed with available information for other canid populations.

MHC diversity

MHC diversity was analysed for each locus and for the three-locus haplotypes. For both we estimated number of alleles (NA), mean number of effective alleles (Ne), mean observed (Ho) and expected (He) heterozygosities and inbreeding coefficient (FIS) using GENALEX 6.5 [68]. Test of deviation from Hardy-Weinberg equilibrium was calculated using GENEPOP 4.6 [69]. Allelic richness (AR) was estimated using FSTAT [70], which performs rarefaction compensating sampling disparity, reducing all populations to a minimum sample size of n = 10. Sequence-level diversity values, such as number of haplotypes (h), number of segregating sites (S), nucleotide (π) and haplotype (HD) diversities, and Watterson’s mutation parameter (θw) were estimated either in each locus or each group (with the three-locus haplotypes) using DnaSP 5.10.

Molecular divergence analysis over all sequence pairs, also either in each locus or each group (with the three-locus haplotypes) was conducted using Mega 7 [71] with a 1000 replicate bootstrap procedure. The best nucleotide and amino acid substitution models were chosen according to the Bayesian Information Criterion and used to compute nucleotide and amino acid evolutionary distances by the maximum-likelihood method through 1000 bootstrap replicates. The relationship between three-locus haplotypes of Iberian and other European wolf populations was reconstructed using the Neighbor Joining (NJ) method and p-distances in Mega 7. Three-locus haplotypes of other European wolf populations were retrieved from previous studies [3436, 38].

Testing positive selection

To detect evidence of selection on MHC loci at the population-level, we estimated Tajima’s D [72] and Fu & Li D* [73] statistics using DnaSP 5.10. Under neutrality, D and D* are not different from 0, and significant negative deviations may indicate a recent selective sweep or population expansion while positive deviations indicate balancing selection or population contraction. Neutrality tests were computed for each MHC locus and for the three-locus haplotypes in the whole Iberian wolf range and in each demographic group. Deviations from neutrality were also assessed by a Z-test (non-synonymous substitution ≠ synonymous substitution; dN ≠ dS) implemented in MEGA 7, using Nei-Gojobori method with Jukes-Cantor correction for multiple substitutions and 1000 bootstrap replicates. Z-tests were computed for all sites, antigen-binding sites inferred from the human orthologue [42], and non-antigen binding sites for each demographic group. Additionally, we also used a “coalescent with recombination” Bayesian approximation implemented in OmegaMap [74] to find signatures of selection in each demographic group. OmegaMap allows both selection parameter (ω) and recombination rate (ρ) to vary along the sequence [75]. Codons with excess of non-synonymous polymorphism compared to synonymous polymorphism are assumed under positive selection. We followed [36] to set distribution and parameter values of priors to represent neutrality. We performed two runs for each locus and each demographic group with 250,000 repeats and a burn-in of 20,000. Convergence of runs was check through summary function implemented in R 3.5.0 (R Core Team, 2018). Codons are considered as positively selected with posterior probabilities (PP) above 95%.

Microsatellite diversity

We selected a set of 39 microsatellite loci previously genotyped by [33] for all but three samples (two from the persistent and one from the isolated groups) sequenced for MHC to estimate neutral genetic diversity in each demographic group. This set included only loci with less than 25% missing data and with no significant departure from Hardy–Weinberg equilibrium after Bonferroni correction estimated using GENALEX 6.5. Genetic diversity, including number of alleles (NA), mean number of effective alleles (Ne), mean allelic richness (AR), mean observed (Ho) and expected (He) heterozygosities and inbreeding coefficient (FIS) were estimated using GENALEX 6.5 and FSTAT For allelic richness all populations were reduced to a minimum sample size of n = 5. Test of deviation from Hardy-Weinberg equilibrium was calculated using GENEPOP 4.6.

Testing for genetic differentiation between demographic groups

To test for statistically significant differences of genetic diversity between the demographic groups, we performed linear mixed effects models (e.g.,[76, 77]) using lme4 [78] and lmerTest [79] implemented in R platform (R Core Team, 2018). We used number of alleles, allelic richness, observed and expected heterozygosities as dependent variables, demographic group and locus type (MHC and microsatellites) as explanatory factors, and locus as random variable. Ten competing models were tested, including i) null, ii) demographic group-dependent, iii) locus type-dependent, iv) dependent on both demographic group and locus type, and v) dependent on both demographic group and locus type plus the association between these two factors. For each of these options we run two models accounting or not for locus random effect (for model syntaxes see Additional file 6: Table S6). We performed independent model testing for each diversity measure. We used Akaike Information Criterion corrected for small sample sizes, AICc, to determine the best model. Models with delta AICc (ΔAICc) lower than 2.0 were assumed as equally plausible. For each diversity measure we present estimates, standard errors and probability inferred from t-value, indicating the precision of the estimates.

Comparing neutral and functional diversity

Expectations that functional diversity exceeds neutral diversity were tested by a one-tailed t-test (p-value), following [35]. The relationship between microsatellite and MHC loci was estimated using a co-inertia analysis (CoA) following [10]. CoA is a multivariate method that estimates the covariance structure between two tables with the same observations (e.g., same individuals). For each genetic matrix (microsatellite and MHC loci) we performed a Principal Component Analysis (PCA) followed by factorial PCA using demographic groups as explanatory variable. For this analysis we used the adegenet and ade4 R packages [80, 81]. We then performed a CoA using the two most important PCA components. In a CoA plot, an arrow and a dot define each demographic group, where the length of the vector connecting arrow and dot is proportional to the divergence between the types of markers. The correlation coefficient RV was used to evaluate the strength of the relationship between the two matrices. Values of RV close to 1 indicate strong correlation between matrices. The significance of the RV coefficient was tested using 1000 permutations, where null hypothesis is that the two matrices are related by chance. Canonical weights for microsatellite and MHC alleles were used to disclose the contribution of each allele to the divergence of demographic groups.

MHC and neutral differentiation between groups was estimated through FST [82] and Jost’s D [83] using Arlequin [84] and the DEMEtics R package [85], respectively. FST and Jost’s D differentiation measures were both used because FST may be affected by highly variable markers [39].



Akaike Information Criterion corrected for small sample sizes


Mean allelic richness


Co-inertia analysis


Non-synonymous substitution


Synonymous substitution


Inbreeding coefficient


number of haplotypes


Haplotype diversity


Expected heterozygosity


Observed heterozygosity


Major histocompatibility complex


Number of alleles


Mean number of effective alleles


Neighbor Joining


Principal Component Analysis


Polymerase chain reaction


Posterior probability


Correlation coefficient


Number of segregating sites


Delta AICc


Watterson’s mutation parameter


Nucleotide diversity


Recombination rate


Selection parameter



We acknowledge the Portuguese Institute for Nature Conservation and Forestry (SMLM/ICNF) for providing wolf tissue samples in Portugal. For collection of wolf samples in Spain, we thank Consejería de Medio Ambiente del Pricipado de Asturias, Consellería de Medio Ambiente de la Xunta de Galicia, Consejería de Medio Ambiente de la Junta de Castilla y León, Gobierno de Cantabria and Parque Nacional Picos de Europa. We are grateful to the Subject Editor, Dr. Helmut Schaschl, and to Dr. Franz Suchentrunk and another anonymous reviewer for valuable comments in a previous version of this manuscript. We thank Pedro Monterroso for helping implementing the linear mixed effects models. This is scientific paper number 23 of the Iberian Wolf Research Team (IWRT).


RGR was supported by a post-doctoral fellowship (project IF/00564/2012/CP0159/CT0003, funded by the Portuguese Foundation for Science and Technology, FCT), RG by the research contract IF/00564/2012 from FCT, and JVLB was supported by a Ramon & Cajal research contract (RYC-2015-18932) from the Spanish Ministry of Economy, Industry and Competitiveness. This work was partially supported by FCT (project PTDC/BIA-EVF/2460/2014 to RG) and by Norte Portugal Regional Operational Program (NORTE2020), under the PORTUGAL 2020 Partnership Agreement, through the European Regional Development Fund (ERDF) (NORTE-01-0145-FEDER-000007). These funding sources had no role in the design of the study and collection, analysis, and interpretation of data and in writing the manuscript.

Availability of data and materials

Microsatellite and MHC genotype data were uploaded to the Open Science Framework repository under the identifier

Authors’ contributions

RG, JVLB and PJE conceived and designed the study; VM conducted laboratory analyses; RG, FA and LL contributed with acquisition of data; RGR, WvdL, VM and RG performed data analyses; RGR, WvdL, VM and RG wrote the first draft of the manuscript. All authors actively contributed and approved the final version of the manuscript.

Ethics approval and consent to participate

The Portuguese Institute for Nature Conservation and Forestry (SMLM/ICNF), Consejería de Medio Ambiente del Pricipado de Asturias, Consellería de Medio Ambiente de la Xunta de Galicia, and Consejería de Medio Ambiente de la Junta de Castilla y León, Gobierno de Cantabria were responsible for the collection of samples and for permits to use them in this study.

Consent for publication

Not applicable.

Competing interests

The authors declare that they have no competing interests.

Publisher’s Note

Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Open Access This article is distributed under the terms of the Creative Commons Attribution 4.0 International License (, which permits unrestricted use, distribution, and reproduction in any medium, provided you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made. The Creative Commons Public Domain Dedication waiver ( applies to the data made available in this article, unless otherwise stated.

Authors’ Affiliations

CIBIO/InBio - Centro de Investigação em Biodiversidade e Recursos Genéticos, Universidade do Porto, Campus de Vairão, 4485-661 Vairão, Portugal
Research Unit of Biodiversity (UO/CSIC/PA), University of Oviedo, 33600 Mieres, Spain
A.RE.NA, S.L. Asesores en Recursos Naturales S.L., 27003 Lugo, Spain
Departamento de Biologia, Faculdade de Ciências, Universidade do Porto, 4169-007 Porto, Portugal


  1. Hartl DL, Clark AG. Principles of population genetics. Sunderland: Sinauer Associates; 1997.Google Scholar
  2. Sommer S. The importance of immune gene variability (MHC) in evolutionary ecology and conservation. Front Zool. 2005;2:1–18.View ArticleGoogle Scholar
  3. Sutton JT, Nakagawa S, Robertson BC, Jamieson IG. Disentangling the roles of natural selection and genetic drift in shaping variation at MHC immunity genes. Mol Ecol. 2011;20:4408–20.View ArticleGoogle Scholar
  4. Ejsmond MJ, Radwan J. MHC diversity in bottlenecked populations: a simulation model. Conserv Genet. 2011;12:129–37.View ArticleGoogle Scholar
  5. Eimes JA, Bollmer JL, Whittingham LA, Johnson JA, van Oosterhout C, Dumm PO. Rapid loss of MHC class II variation in a bottlenecked population is explained by drift and loss of copy number variation. J Evol Biol. 2011;24:1847–56.View ArticleGoogle Scholar
  6. Holderegger R, Kamm U. F G. adaptive vs. neutral genetic diversity: implications for landscape genetics. Landsc Ecol. 2006;21:797–807.View ArticleGoogle Scholar
  7. Brumfield RT, Beerli P, Nickerson DA, Edwards SV. The utility of single nucleotide polymorphisms in inferences of population history. Trends Ecol Evol. 2003;18:249–56. ArticleGoogle Scholar
  8. van Tienderen PH, de Haan AA, van der Linden CG, Vosman B. Biodiversity assessment using markers for ecologically important traits. Trends Ecol Evol. 2002;17:577–82. ArticleGoogle Scholar
  9. Agudo R, Alcaide M, Rico C, Lemus JA, Blanco G, Hiraldo F, et al. Major histocompatibility complex variation in insular populations of the Egyptian vulture: inferences about the roles of genetic drift and selection. Mol Ecol. 2011;20:2329–40.View ArticleGoogle Scholar
  10. Rico Y, Ethier DM, Davy CM, Sayers J, Weir RD, Swanson BJ, et al. Spatial patterns of immunogenetic and neutral variation underscore the conservation value of small, isolated American badger populations. Evol Appl. 2016;9:1271–84.View ArticleGoogle Scholar
  11. Bernatchez L, Landry C. MHC studies in nonmodel vertebrates: what have we learned about natural selection in 15 years? J Evol Biol. 2003;16:363–77.View ArticleGoogle Scholar
  12. Ujvari B, Belov K. Major histocompatibility complex (MHC) markers in conservation biology. Int J Mol Sci. 2011;12:5168–86.View ArticleGoogle Scholar
  13. Hedrick PW, Lee RN, Buchanan C. Canine parvovirus enteritis, canine distemper, and major histocompatibility complex genetic variation in Mexican wolves. J Wildl Dis. 2003;39:909–13. ArticlePubMedGoogle Scholar
  14. Schwensow N, Fietz J, Dausmann KH, Sommer S. Neutral versus adaptive genetic variation in parasite resistance: importance of major histocompatibility complex supertypes in a free-ranging primate. Heredity. 2007;99:265–77. ArticleGoogle Scholar
  15. Edwards SV, Hedrick PW. Evolution and ecology of MHC molecules: from genomics to sexual selection. Trends Ecol Evol. 1998;13:305–11. ArticlePubMedGoogle Scholar
  16. Piertney SB, Oliver MK. The evolutionary ecology of the major histocompatibility complex. Heredity. 2005;96:7–21. ArticleGoogle Scholar
  17. Kalbe M, Eizaguirre C, Dankert I, Reusch TBH, Sommerfeld RD, Wegner KM, et al. Lifetime reproductive success is maximized with optimal major histocompatibility complex diversity. Proc R Soc B Biol Sci. 2009;276:925–34. ArticleGoogle Scholar
  18. Sin YW, Annavi G, Newman C, Buesching C, Burke T, Macdonald DW, et al. MHC class II-assortative mate choice in European badgers (Meles meles). Mol Ecol. 2015;24:3138–50.View ArticleGoogle Scholar
  19. Aguilar A, Roemer G, Debenham S, Binns M, Garcelon D, Wayne RK. High MHC diversity maintained by balancing selection in an otherwise genetically monomorphic mammal. Proc Natl Acad Sci. 2004;101:3490–4. ArticlePubMedGoogle Scholar
  20. Marsden CD, Woodroffe R, Mills MGL, McNutt JW, Creel S, Groom R, et al. Spatial and temporal patterns of neutral and adaptive genetic variation in the endangered African wild dog (Lycaon pictus). Mol Ecol. 2012;21:1379–93. ArticleGoogle Scholar
  21. Palomares F, Godoy JA, Lopez-Bao JV, Rodriguez A, Roques S, Casas-Marce M, et al. Possible extinction vortex for a population of Iberian lynx on the verge of extirpation. Conserv Biol. 2012;26:689–97.View ArticleGoogle Scholar
  22. Ishibashi Y, Oi T, Arimoto I, Fujii T, Mamiya K, Nishi N, et al. Loss of allelic diversity in the MHC class II DQB gene in western populations of the Japanese black bear Ursus thibetanus japonicus. Conserv Genet. 2017;18:247–60. ArticleGoogle Scholar
  23. Chapron G, Kaczensky P, Linnell JDC, von Arx M, Huber D, Andren H, et al. Recovery of large carnivores in Europe’s modern human-dominated landscapes. Science. 2014;346:1517–9.Google Scholar
  24. Álvares FB, Blanco JC, Correia J, Cortés Y, Costa G, Llaneza L, et al. Wolf status and conservation in the Iberian Peninsula. In: Frontiers of wolf recovery. Colorado Springs: International Wolf Center; 2005. p. 66–7.Google Scholar
  25. Núñez-Quirós P, Garcia-Lavandera R, Llaneza L. Analysis of historical wolf (Canis lupus) distributions in Galicia: 1850, 1960 and 2003. Ecologia. 2007;21:195–205.Google Scholar
  26. Petrucci-Fonseca F. O lobo (Canis lupus signatus Cabrera, 1907) em Portugal. Problemática da sua conservação: Universidade de Lisboa; 1990.Google Scholar
  27. Valverde JA. El lobo español. Montes. 1971;159:229–41.Google Scholar
  28. Blanco JC, Cortés Y. Dispersal patterns, social structure and mortality of wolves living in agricultural habitats in Spain. J Zool. 2007;273:114–24. ArticleGoogle Scholar
  29. Sommer R, Benecke N. Late-Pleistocene and early Holocene history of the canid fauna of Europe (Canidae). Mamm Biol. 2005;70:227–41. ArticleGoogle Scholar
  30. Bencatel J, Ferreira CC, Barbosa AM, Rosalino LM, Álvares F. Research trends and geographical distribution of mammalian carnivores in Portugal (SW Europe). PLoS One. 2018;13:e0207866. ArticlePubMedPubMed CentralGoogle Scholar
  31. Hindrikson M, Remm J, Pilot M, Godinho R, Stronen AV, Baltrunaite L, et al. Wolf population genetics in Europe: a systematic review, meta-analysis and suggestions for conservation and management. Biol Rev Camb Philos Soc. 2017;92:1601–29.View ArticleGoogle Scholar
  32. López-Bao JV, Blanco JC, Rodríguez A, Godinho R, Sazatornil V, Alvares F, et al. Toothless wildlife protection laws. Biodivers Conserv. 2015;24:2105–8.View ArticleGoogle Scholar
  33. Silva P, López-Bao JV, Llaneza L, Álvares F, Lopes S. Blanco JC, et al. Sci Rep: Cryptic population structure reveals low dispersal in Iberian wolves; 2018.Google Scholar
  34. Arbanasic H, Huber D, Kusak J, Gomercic T, Hrenovic J, Galov A. Extensive polymorphism and evidence of selection pressure on major histocompatibility complex DLA-DRB1, DQA1 and DQB1 class II genes in Croatian grey wolves. Tissue Antigens. 2013;81:19–27.View ArticleGoogle Scholar
  35. Galaverni M, Caniglia R, Fabbri E, Lapalombella S, Randi E. MHC variability in an isolated wolf population in Italy. J Hered. 2013;104:601–12.View ArticleGoogle Scholar
  36. Niskanen AK, Kennedy LJ, Ruokonen M, Kojola I, Lohi H, Isomursu M, et al. Balancing selection and heterozygote advantage in major histocompatibility complex loci of the bottlenecked Finnish wolf population. Mol Ecol. 2014;23:875–89.View ArticleGoogle Scholar
  37. Seddon JM, Ellegren H. MHC class II genes in European wolves: a comparison with dogs. Immunogenetics. 2002;54:490–500.View ArticleGoogle Scholar
  38. Seddon JM, Ellegren H. A temporal analysis shows major histocompatibility complex loci in the Scandinavian wolf population are consistent with neutral evolution. Proc R Soc B Biol Sci. 2004;271:2283–91. ArticleGoogle Scholar
  39. Strand TM, Segelbacher G, Quintela M, Xiao L, Axelsson T, Höglund J. Can balancing selection on MHC loci counteract genetic drift in small fragmented populations of black grouse? Ecol Evol. 2012;2:341–53.View ArticleGoogle Scholar
  40. Alcaide M. On the relative roles of selection and genetic drift in shaping MHC variation. Mol Ecol. 2010;19:3842–4.View ArticleGoogle Scholar
  41. Loiseau C, Richard M, Garnier S, Chastel O, Julliard R, Zoorob R, et al. Diversifying selection on MHC class I in the house sparrow (Passer domesticus). Mol Ecol. 2009;18:1331–40.View ArticleGoogle Scholar
  42. Brown JH, Jardetzky T, Saper MA, Samraoui B, Bjorkman PJ, Wiley DC. A hypothetical model of the foreign antigen binding site of class II histocompatibility molecules. Nature. 1988;332:845–50.View ArticleGoogle Scholar
  43. Pilot M, Branicki W, Jedrzejewski W, Goszczyński J, Jedrzejewska B, Dykyy I, et al. Phylogeographic history of grey wolves in Europe. BMC Evol Biol. 2010;10:104. ArticlePubMedPubMed CentralGoogle Scholar
  44. Randi E. Genetics and conservation of wolves Canis lupus in Europe. Mamm Rev. 2011;41:99–111. ArticleGoogle Scholar
  45. Montana L, Caniglia R, Galaverni M, Fabbri E, Ahmed A, Bolfíková BČ, et al. Combining phylogenetic and demographic inferences to assess the origin of the genetic diversity in an isolated wolf population. PLoS One. 2017;12:1–19. ArticleGoogle Scholar
  46. Kennedy LJ, Barnes A, Short A, Brown JJ, Lester S, Seddon J, et al. Canine DLA diversity: 1. New alleles and haplotypes. Tissue Antigens. 2007;69(SUPPL. 1):272–88.View ArticleGoogle Scholar
  47. Kennedy LJ, Barnes A, Happ GM, Quinnell RJ, Bennett D, Angles JM, et al. Extensive interbreed, but minimal intrabreed, variation of DLA class II alleles and haplotypes in dogs. Tissue Antigens. 2002;59:194–204.View ArticleGoogle Scholar
  48. Kennedy LJ, Angles JM, Barnes A, Carmichael LE, Radford AD, Ollier WER, et al. DLA-DRB1, DQA1, and DQB1 alleles and haplotypes in north American gray wolves. J Hered. 2007;98:491–9.View ArticleGoogle Scholar
  49. Godinho R, Llaneza L, Blanco JC, Lopes S, Álvares F, García EJ, et al. Genetic evidence for multiple events of hybridization between wolves and domestic dogs in the Iberian Peninsula. Mol Ecol. 2011;20:5154–66.View ArticleGoogle Scholar
  50. Fernández JM, Ruiz de Azua N. Historical dynamics of a declining wolf population: persecution vs. prey reduction. Eur J Wildl Res. 2009;56:169–79.View ArticleGoogle Scholar
  51. Rico M, Torrente JP. Caza y rarifacación del lobo en España: investigación histórica y conclusiones biológicas. Galemys. 2000;12:163–79.Google Scholar
  52. Aspi J, Roininen E, Ruokonen M, Kojola I, Vilà C. Genetic diversity, population structure, effective population size and demographic history of the Finnish wolf population. Mol Ecol. 2006;15:1561–76. ArticlePubMedGoogle Scholar
  53. Eizaguirre C, Lenz TL, Kalbe M, Milinski M. Rapid and adaptive evolution of MHC genes under parasite selection in experimental vertebrate populations. Nat Commun. 2012;3:621–6. ArticlePubMedPubMed CentralGoogle Scholar
  54. Ciborowski K, Jordan WC, Garcia De Leaniz C, Consuegra S. Temporal and spatial instability in neutral and adaptive (MHC) genetic variation in marginal salmon populations. Sci Rep. 2017;7 February:1–11.
  55. Eckert CG, Samis KE, Lougheed SC. Genetic variation across species’ geographical ranges: the central-marginal hypothesis and beyond. Mol Ecol. 2008;17:1170–88.View ArticleGoogle Scholar
  56. Volis S, Ormanbekova D, Shulgina I. Role of selection and gene flow in population differentiation at the edge vs. interior of the species range differing in climatic conditions. Mol Ecol. 2016;25:1449–64.View ArticleGoogle Scholar
  57. O’Brien SJ, Evermann JF. Interactive influence of infectious disease and genetic diversity in natural populations. Trends Ecol Evol. 1988;3:254–9.View ArticleGoogle Scholar
  58. Santos N, Almendra C, Tavares L. Serologic survey for canine distemper virus and canine parvovirus in free-ranging wild carnivores from Portugal. J Wildl Dis. 2009;45:221–6.View ArticleGoogle Scholar
  59. Figueiredo A, Oliveira L, Madeira de Carvalho L, Fonseca C, Torres RT. Parasite species of the endangered Iberian wolf (Canis lupus signatus) and a sympatric widespread carnivore. Int J Parasitol Parasites Wildl. 2016;5:164–7. Scholar
  60. Conceição-Neto N, Godinho R, Álvares F, Yinda CK, Deboutte W, Zeller M, et al. Viral gut metagenomics of sympatric wild and domestic canids, and monitoring of viruses: insights from an endangered wolf population. Ecol Evol. 2017;7:4135–46.View ArticleGoogle Scholar
  61. Müller A, Silva E, Santos N, Thompson G. Domestic dog origin of canine distemper virus in free-ranging wolves in Portugal as revealed by hemagglutinin gene characterization. J Wildl Dis. 2011;47:725–9.View ArticleGoogle Scholar
  62. Zhai W, Nielsen R, Slatkin M. An investigation of the statistical power of neutrality tests based on comparative and population genetic data. Mol Biol Evol. 2009;26:273–83. ArticlePubMedGoogle Scholar
  63. Hall TA. BioEdit: a user-friendly biological sequence alignment editor and analysis program for Windows 95/98/NT. Nucleic Acids Symposium Series. 1999;41:95–98. doi:citeulike-article-id:691774.Google Scholar
  64. Librado P, Rozas J. DnaSP v5: a software for comprehensive analysis of DNA polymorphism data. Bioinformatics. 2009;25:1451–2.View ArticleGoogle Scholar
  65. Stephens M, Donnelly P. A comparison of Bayesian methods for haplotype reconstruction from population genotype data. Am J Hum Genet. 2003;73:1162–9.View ArticleGoogle Scholar
  66. Bos DH, Gopurenko D, Williams RN, DeWoody JA. Inferring population history and demography using microsatellites, mitochondrial DNA, and major histocompatibility complex (MHC) genes. Evolution. 2008;62:1458–68. ArticleGoogle Scholar
  67. Garrick RC, Sunnucks P, Dyer RJ. Nuclear gene phylogeography using PHASE: dealing with unresolved genotypes, lost alleles, and systematic bias in parameter estimation. BMC Evol Biol. 2010;10:118. ArticlePubMedPubMed CentralGoogle Scholar
  68. Peakall R, Smouse PE. GenAlEx 6.5: genetic analysis in excel. Population genetic software for teaching and research--an update. Bioinformatics. 2012;28:2537–9.View ArticleGoogle Scholar
  69. Rousset F. genepop’007: a complete re-implementation of the genepop software for windows and Linux. Mol Ecol Resour. 2008;8:103–6. ArticlePubMedGoogle Scholar
  70. Goudet J. FSTAT (version 1.2): a computer program to calculate F-statistics. J Hered. 1995;86:485–6. ArticleGoogle Scholar
  71. Kumar S, Stecher G, Tamura K. MEGA7: molecular evolutionary genetics analysis version 7.0. Mol Biol Evol. 2016;33:1870–4.View ArticleGoogle Scholar
  72. Tajima F. Statistical method for testing the neutral mutation hypothesis by DNA polymorphism. Genetics. 1989;123:585–95.PubMedPubMed CentralGoogle Scholar
  73. Fu YX, Li WH. Statistical tests of neutrality of mutations. Genetics. 1993;133:693–709.PubMedPubMed CentralGoogle Scholar
  74. Wilson DJ, McVean G. Estimating diversifying selection and functional constraint in the presence of recombination. Genetics. 2006;172:1411–25. ArticleGoogle Scholar
  75. Wilson DJ, McVean G. Estimating diversifying selection and functional constraint in the presence of recombination. Genetics. 2006;172:1411–25.View ArticleGoogle Scholar
  76. Awadi A, Ben Slimen H, Smith S, Knauer F, Makni M, Suchentrunk F. Positive selection and climatic effects on MHC class II gene diversity in hares (Lepus capensis) from a steep ecological gradient. Sci Rep. 2018;8:1–14.Google Scholar
  77. Demirbaş Y, Albayrak İ, Özkan Koca A, Stefanović M, Knauer F, Suchentrunk F. Spatial genetics of brown hares (Lepus europaeus Pallas, 1778) from Turkey: different gene pool architecture on either side of the Bosphorus? Mamm Biol. 2019;94:77–85.View ArticleGoogle Scholar
  78. Bates D, Mächler M, Bolker B, Walker S. Fitting linear mixed-effects models using lme4. J Stat Softw. 2014;67:1–48.
  79. Kuznetsova A, Brockhoff PB, Christensen RHB. lmerTest Package: Tests in Linear Mixed Effects Models. J Stat Softw. 2017;82:1–26. ArticleGoogle Scholar
  80. Dray S, Dufour A-B. The ade4 package: implementing the duality diagram for ecologists. J Stat Softw. 2007;22:1–20.
  81. Jombart T. Adegenet: a R package for the multivariate analysis of genetic markers. Bioinformatics. 2008;24:1403–5.View ArticleGoogle Scholar
  82. Slatkin M. A measure of population subdivision based on microsatellite allele frequencies. Genetics. 1995;139:457–62.PubMedPubMed CentralGoogle Scholar
  83. Jost L. GST and its relatives do not measure differentiation. Mol Ecol. 2008;17:4015–26. ArticlePubMedGoogle Scholar
  84. Excoffier L, Lischer HEL. Arlequin suite ver 3.5: a new series of programs to perform population genetics analyses under Linux and windows. Mol Ecol Resour. 2010;10:564–7.View ArticleGoogle Scholar
  85. Gerlach G, Jueterbock A, Kraemer P, Deppermenn J, Harmand P. Calculations of population differentiation based on GST and D: forget GST but not all of statistics! Mol Ecol. 2010;19:3845–52. ArticlePubMedGoogle Scholar


© The Author(s). 2019