Balancing selection, genetic drift, and human‐mediated introgression interplay to shape MHC (functional) diversity in Mediterranean brown trout

Abstract The extraordinary polymorphism of major histocompatibility complex (MHC) genes is considered a paradigm of pathogen‐mediated balancing selection, although empirical evidence is still scarce. Furthermore, the relative contribution of balancing selection to shape MHC population structure and diversity, compared to that of neutral forces, as well as its interaction with other evolutionary processes such as hybridization, remains largely unclear. To investigate these issues, we analyzed adaptive (MHC‐DAB gene) and neutral (11 microsatellite loci) variation in 156 brown trout (Salmo trutta complex) from six wild populations in central Italy exposed to introgression from domestic hatchery lineages (assessed with the LDH gene). MHC diversity and structuring correlated with those at microsatellites, indicating the substantial role of neutral forces. However, individuals carrying locally rare MHC alleles/supertypes were in better body condition (a proxy of individual fitness/parasite load) regardless of the zygosity status and degree of sequence dissimilarity of MHC, hence supporting balancing selection under rare allele advantage, but not heterozygote advantage or divergent allele advantage. The association between specific MHC supertypes and body condition confirmed in part this finding. Across populations, MHC allelic richness increased with increasing admixture between native and domestic lineages, indicating introgression as a source of MHC variation. Furthermore, introgression across populations appeared more pronounced for MHC than microsatellites, possibly because initially rare MHC variants are expected to introgress more readily under rare allele advantage. Providing evidence for the complex interplay among neutral evolutionary forces, balancing selection, and human‐mediated introgression in shaping the pattern of MHC (functional) variation, our findings contribute to a deeper understanding of the evolution of MHC genes in wild populations exposed to anthropogenic disturbance.


| INTRODUC TI ON
The major histocompatibility complex (MHC) is an extraordinarily polymorphic, gene-rich family involved in the immune response of vertebrates. Proteins encoded by classical MHC genes recognize and bind antigens derived from parasites and pathogens initiating the adaptive immune cascade (Kaufman, 2018). The range of antigens that a MHC protein can recognize is determined by the amino acid composition of particular highly variable sites, referred to as antigenbinding sites, encoded by the exon 2 in classical MHC class II genes.
Positive selection for nonsynonymous mutations, mainly in antigenbinding sites, and recombination contribute to generating the characteristic high levels of allele and sequence polymorphism at MHC loci . Notably, given the adaptive significance and the long-time persistence of MHC polymorphism, MHC has become a paradigm for investigating pathogen-mediated balancing selection. Two mutually non-exclusive forms of balancing selection are commonly postulated: the rare allele advantage (RAA) and the heterozygote advantage (HA) Spurgin & Richardson, 2010).
The RAA, also called negative frequency-dependent selection, is thought to mainly contribute to maintaining MHC polymorphism by preventing the loss of rare variants. RAA arises from host-parasite arms race (Red Queen dynamics) in which pathogens rapidly adapt to escape recognition by the commonest host MHC alleles, ultimately resulting in cyclical fluctuations of MHC allele frequencies.
Computational studies supported the RAA theory (e.g., Borghans et al., 2004) and experimental trials conducted in semi-natural condition demonstrated that fish carrying novel MHC (functional) variants experienced a reduction of parasite infection . As an interesting consequence, the RAA may facilitate the introduction and spread of novel, initially rare alleles into populations via introgression. For instance, Dudek et al. (2019) demonstrated that introgression coupled with balancing selection caused a massive allelic exchange and overall higher MHC variation in a newt hybrid zone. Under the HA model, MHC heterozygous individuals are predicted to better cope with pathogens because of their ability to recognize a wider range of antigens compared to homozygotes. By extension, heterozygotes carrying genetically divergent MHC alleles should recognize a still broader spectrum of antigens than heterozygotes carrying similar alleles, that is, the divergent allele advantage (DAA) hypothesis . Evidence for HA, such as a significant positive association between MHC heterozygosity and individual fitness-related traits (e.g., parasite diversity, infection intensity, body condition, survival), has been detected sometimes (e.g., Evans & Neff, 2009;Niskanen et al., 2014;Savage et al., 2019).
As mentioned above, contemporary pathogen-mediated balancing selection can be assessed at the individual level by testing the association between fitness-related traits and components of MHC diversity, for example, sequence dissimilarity, frequency of variants, zygosity, and number of alleles (Lenz et al., 2009(Lenz et al., , 2013Seifertová et al., 2016;Trachtenberg et al., 2003). On the other hand, exploring selection over the history of populations requires comparing levels of diversity and/or genetic structuring between MHC and putative neutral-evolving loci, such as microsatellites (Garrigan & Hedrick, 2003). Similar levels of structuring at MHC and neutral loci are expected when demographic processes overcome selection in shaping MHC diversity across populations. Conversely, a weaker or stronger MHC structuring compared to that of neutral loci may indicate either uniform balancing selection or fluctuating selection/local adaptation, respectively (Garrigan & Hedrick, 2003;Spurgin & Richardson, 2010). In addition, mixed and/or scaledependent patterns are possible (Miller et al., 2001;Talarico, Babik, Marta, Pietrocini, et al., 2019). Therefore, which balancing selection processes underpin MHC diversity is still an open question, and the relative contribution of each process remains unclear, especially in wild populations under natural or semi-natural condition.
Another innovative perspective concerns the functionality of MHC diversity, with MHC studies increasingly focusing on MHC supertypes. Supertypes are classes of sequences with similar biochemical properties at those amino acid positions determining the specificity of antigen binding (Schwensow et al., 2007;Sette & Sidney, 1998). The supertyping approach has the advantage of accounting for the (actual) functional diversity of MHC variants, which is biologically relevant especially in ecological studies where phenotypic effects of MHC alleles should be examined (Naugler & Liwski, 2008;Trachtenberg et al., 2003). Moreover, supertyping implies reducing the usually high number of MHC variants that often represents a challenge for statistical analyses (Schwensow et al., 2007;Sommer, 2005).
For instance, associations between specific supertypes and pathogen resistance, lifetime reproductive success, and survival were found in the wild great tit (Sepil et al., 2012(Sepil et al., , 2013, and an advantage for rare MHC supertypes in HIV disease progression was detected in human (Trachtenberg et al., 2003). On an evolutionary perspective, Lighten et al. (2017) innovatively proposed a different role of selection on MHC supertypes compared to alleles. According to this hypothesis, balancing selection would maintain functional (supertype) diversity, while positive selection and/or genetic drift would erode allelic diversity, giving rise to a population structuring more pronounced for MHC alleles than supertypes. Nevertheless, Ejsmond et al. (2018) questioned empirical support to this hypothesis, and a successive study failed to confirm such a dual role of selection in the same target species (Herdegen-Radwan et al., 2020).
Brown trout is one of the most genetically and phenotypically variable species, referred to as Salmo trutta complex, and among the most broadly distributed salmonids worldwide (Lobón-Cerviá & Sanz, 2018).
The species complex includes several genetic lineages, which extensively hybridize with hatchery-reared brown trout (typically domestic lineages of Atlantic origin), released into the wild for restocking purposes for more than one century (reviewed in Lobón-Cerviá & Sanz, 2018).
Various native lineages occur in Italy and inhabit inland freshwaters, mostly upper mountain streams and, occasionally, lowland cold waters.
Previous studies targeting brown trout revealed high polymorphisms and historical positive selection on the MHC class II gene (the single-copy MHC-DAB) and examined multiple aspects of MHC variation and evolution. Campos et al. (2006) (Forsberg et al., 2007;Jacob et al., 2010).
Here, we characterize putative adaptive (MHC) and neutral (microsatellite loci) variation in genetically unrelated wild populations of brown trout from central Italy. The study area was previously recognized as a hotspot of genetic diversity for native brown trout lineages, although almost all populations were variously subject to stocking with domestic brown trout (Fabiani et al., 2018;Rossi et al., 2019). Such activities may ultimately favor introgressionwhich refers to the consequence of hybridization between lineages rather than distinct species in this case-and maintenance of novel genetic variation into wild native populations through balancing selection mechanisms, especially at putatively adaptive genes such as MHC. For these reasons, we used brown trout as a suitable model to investigate and disentangle evolutionary forces driving both neutral and adaptive diversity in a context of genetic disturbance caused by human activities. Specifically, the purposes of this study are twofold.
Firstly, we aim to assess the relative contribution of selection and nonselective forces, i.e., genetic drift and (human-mediated) gene flow, to drive putatively adaptive MHC variation in the recent past, by comparing genetic structuring and diversity at MHC and microsatellite loci, while controlling for the extent of admixture between native and domestic lineages. Secondly, we look for signatures of current balancing selection by testing the association between individual body condition (a putative proxy of parasite load) and frequency of MHC variants carried by an individual (i.e., testing for TA B L E 1 Genetic variation at the MHC, when considering alleles (MHC-DAB) or supertypes (MHC-ST), and 11 microsatellite loci (STR) for six wild brown trout populations from central Italy

MHC-DAB alleles
rare allele advantage), MHC zygosity (i.e., testing for heterozygote advantage), and within-individual MHC dissimilarity (i.e., testing for divergent allele advantage). We also explored associations between specific MHC functional variants (supertypes) and body condition.

| Sample collection and DNA extraction
During summer 2019, 156 wild brown trout were sampled by electrofishing from six putatively unrelated populations in central Italy (Latium) (Table 1; Figure 1). For each fish, we collected a fin clip and recorded weight and standard length (i.e., the total length excluding the caudal fin) before releasing in the field. Fin clips were stored in ethanol at −20°C until total DNA extraction, which was carried out following a salt-extraction protocol (Aljanabi & Martinez, 1997). The "Direzione agricoltura, promozione della filiera e della cultura del cibo, caccia e pesca" and "Direzione Regionale Politiche Ambientali e Ciclo dei Rifiuti" of the Latium Region approved all sampling procedures (authorization no. G10101-27/7/2019 and Q169/2019).  (Sebastian et al., 2016). Both laboratory procedures and bioinformatics are extensively described in Appendix S1.

| Amplification and genotyping of MHC, STR loci, and the LDH gene
Individuals were genotyped at 11 polymorphic microsatellite loci (hereafter STR) to estimate neutral diversity (details on amplification and genotyping procedures are given in Appendix S2). Additionally, we followed the RFLP-based method described in McMeel et al. (2001) to obtain individual profiles of the lactate dehydrogenase (LDH) C1 gene, which is routinely used to assess the admixture between native and domestic lineages: the LDH*90 allele is fixed within European hatchery stocks of domestic brown trout, while the *100 allele naturally occurs in wild, Mediterranean-native brown trout populations.
Two approaches were adopted to look for signatures of historical positive selection on MHC-DAB sequences. Firstly, the one-tailed Z test performed in MEGA7 was used to test for ω > 1 (i.e., the ratio between nonsynonymous and synonymous substitution rates, according to the Nei-Gojobori method of pairwise comparison with the Jukes-Cantor correction; Nei & Gojobori, 1986) in three alignment partitions: all sites; only positively selected sites (PSSs) detected in brown trout (see below); and only non-PSSs. Secondly, we choose the best-fitting model, based on the minimum Akaike information criterion (AIC) value, among three codon-based evolution models as performed in PAML v4.8 (Yang, 2007): M0 assuming a single ω value for all codons; M7 considering purifying selection and neutrality (ω ≤ 1); and M8 allowing evolution by positive selection (ω > 1) in a proportion of sites. We run the Single BreakPoint recombination algorithm (SBP; Kosakovsky Pond et al., 2006) on the Datamonkey server (Delport et al., 2010) to detect putative recombination events.
To account for functional MHC diversity, MHC alleles were grouped into theoretical supertypes following a widely adopted statistical approach (Doytchinova & Flower, 2005), which accounts for immunological (Sidney et al., 1996) and laboratory (Sandberg et al., 1998) concepts. To this end, we performed six tests to identify positively selected sites (PSSs), namely putative amino acid sites involved in peptide binding, and conservatively retained only PSSs detected by at least two tests among the Bayes empirical Bayes (BEB; Yang et al., 2005); the fast unconstrained Bayesian approximation (FUBAR; Murrell et al., 2013), the mixed-effects model of evolution (MEME; Murrell et al., 2012), the single-likelihood ancestor counting (SLAC), the fixed effects likelihood (FEL), and the randomeffects likelihood (REL) analyses (Kosakovsky Pond & Frost, 2005).
The analyses involved novel alleles plus 98 unique sequences retrieved from GenBank (see above) to improve detection accuracy and robustness of PSSs. Then, detected PSSs were extracted from sequences and each PSS was characterized according to five physicochemical descriptors as in Sandberg et al. (1998). Supertypes were defined using the K-means clustering algorithm implemented in the adegenet R-package (Jombart, 2008) and the most probable number of supertypes was that returning the lowest value of the Bayesian information criterion (BIC). Finally, all alleles within each cluster were collapsed into a single MHC supertype. For downstream analyses, we used MHC genotypes resulting from either alleles or supertypes (referred to as MHC-DAB and MHC-ST "loci," respectively).

| Population structure and diversity at MHC and STR
As a preliminary step, we performed the STR-based Bayesian clustering analysis in STRUCTURE v.2.3.4 (Pritchard et al., 2000) to ensure that studied populations were genetically unrelated. The STRUCTURE clustering confirmed our a priori assumption and, consequently, populations were treated separately for most of the downstream analyses. Methods and results concerning the STRUCTURE analysis are reported in Appendix S3.
In FreeNA (Chapuis & Estoup, 2007), we explored the frequency of null alleles for each locus and population according to the Dempster's EM algorithm. For both MHC-DAB and STR, departures from Hardy-Weinberg equilibrium and linkage disequilibrium were tested in GenePop v.4.7 (Rousset, 2008) within each population, adjusting p-values for multiple testing (Holm-Bonferroni correction).
To estimate population genetic variability, we calculated the number of alleles, the number of private alleles, and the observed and expected heterozygosity in GenAlEx 6.5 (Peakall & Smouse, 2012).
Additionally, to compare diversity among populations with different sample sizes (see Table 1), we calculated allelic richness by rarefaction (minimum sample size = 23) in the PopGenReport R-package (Adamack & Gruber, 2014). All indices were computed for STR,

MHC-DAB, and MHC-ST separately.
We fitted linear models (lm function in R; R Core Team, 2020) to evaluate the contribution of neutral variation (STR richness) and admixture (LDH-admixture) to determine MHC diversity (MHC-DAB richness and MHC-ST richness) across populations. The fitted models followed the general formula: Note that LDH-admixture indicates the frequency of the less common LDH allele (typically the domestic *90 allele) within a population. As LDH is a biallelic locus, a population is completely admixed between native and domestic lineages when LDH-admixture is 50% (i.e., 50% of *100 and 50% of *90). Conversely, low levels of LDHadmixture indicate a "pure" population, either native or domestic. To circumvent the collinearity among STR richness and LDH-admixture (Spearman r s = 0.77, p = .10), we used a sequential regression approach (Dormann et al., 2013;Graham, 2003). Here, LDH-admixture was regressed against STR richness and the residuals of this model were used to represent the portion of variation in LDH-admixture not accounted for by STR richness. Both LDH-admixture and STR richness were scaled to zero mean and unit variance before modeling to allow a full comparability of regression coefficients.
In GenAlEx, we calculated overall and population-pairwise genetic differentiation (i.e., the Nei's G ST for highly polymorphic markers with correction for a small number of populations) for STR loci, the MHC-DAB, and the MHC-ST. Also, since collapsing MHC alleles into supertypes may ultimately cause biased patterns of population structuring, we tested whether the observed supertype-based population G ST was significantly different than expected by chance.
Following Lighten et al. (2017), we randomly assigned 10,000 times MHC alleles to supertypes, calculating each time per-population G ST . Observed G ST values significantly higher or lower than the null distribution of G ST may reflect a biological phenomenon (i.e., strong diversifying or balancing selection, respectively), while nonsignificant differences may reflect a bioinformatic artifact (Herdegen-Radwan et al., 2020;Lighten et al., 2017), thus suggesting caution in interpreting the observed degree of functional population structuring.

| Relationship between body condition and MHC
We used mixed-effects linear models (performed in the R-package lme4; Bates et al., 2015) to look for signatures of current balancing selection acting on MHC in the form of heterozygote advantage, divergent allele advantage, and rare allele advantage. Specifically, we assessed the effect of MHC zygosity and frequency of MHC variants (model 1) and MHC dissimilarity (model 2) on body condition, while controlling for the genetic background (STR heterozygosity and LDH type) and origin (population) of individuals-the latter was included in the models as a random intercept to account for possible ecological differences (e.g., food availability) across populations that may affect body condition. Model formulas were as follows: MHC richness ∼ STR richness + LDH admixture where: • Body condition refers to the Fulton's condition factor (computed as 100 · weight/standard length 3 ), a commonly used index of general condition which assumes that heavier fish of a given length are in better condition.
• MHC frequency index is simply the sum of the within-population in homozygotes-we adopted the previously described sequential regression approach (Dormann et al., 2013;Graham, 2003).
MHC zygosity was regressed against MHC frequency index using a generalized linear model with binomial error distribution, and the residuals from this model were used to represent the portion of variation in MHC zygosity not accounted for by MHC frequency index.
• MHC dissimilarity is the pairwise distance between MHC-DAB sequences carried by an individual, measured in MEGA7 as the nucleotide Kimura two-parameter distance (correcting for multiple hits) and the Poisson-corrected amino acid distance (adjusting for multiple substitutions). Note that MHC dissimilarity is set to 0 for MHC-DAB homozygous individuals.
• STR heterozygosity is the proportion of heterozygous STR loci, as performed in GenAlEx, approximating the genome-wide neutral diversity.
• LDH type is a three-level factor indicating the genotype at the LDH locus, namely "native" = 100/100, "hybrid" = 90/100, and "domestic" = 90/90. Finally, we examined the influence of specific MHC variants on body condition following the approach described in Sepil et al. (2012). In brief, we compared the baseline model (i.e., the previously defined model 1, which accounts for the effects of frequency of MHC variants carried by an individual, MHC zygosity status, genetic background, and population origin) to models including the baseline and the effect of single specific MHC variants (coded as presence/absence). Then, we ranked models by increasing AICc (AIC corrected for small sample size), namely from the best to the worst, and inspected the significance of beta coefficients (effect sizes) of models fitting the data substantially better than the baseline model, that is, those with an AICc at least 2 points lower than the baseline (Burnham & Anderson, 2002).

| MHC characterization
We

| Genetic variation in studied populations
The frequency of null alleles was negligible for all loci and populations, being barely >10% in MEL and SCR at the ONEPHI2 locus only (Appendix S7). A significant heterozygote deficiency was found in RAP at the STR locus SSOSL417 and in SIM at the MHC-DAB (Appendix S8). Linkage disequilibrium was statistically significant (p < .05) in only 10 out of 396 pairwise comparisons, nine of which occurred in RAP, hence indicating that there was no physical linkage between any loci. Consequently, we retained all STR loci as the reference for neutral variation for subsequent analyses.
Model 1: Body condition ∼ MHC frequency index + MHC zygosity + STR heterozygosity + LDH type + population Model 2: Body condition ∼ MHC dissimilarity + STR heterozygosity + LDH type + population A summary of population genetic variation is shown in Table 1. MHC-DAB allelic richness ranged between 6.3 and 19.0. Private alleles per population were on average 40.4% (± 11.7% SE). In each population, one to four alleles had a higher frequency (>10%), while the others occurred at lower frequencies (Figure 1a; Appendix S5).
We found 6-10 MHC supertypes per population; supertype-12 occurred only in RAP, but others were shared among at least three populations (Figure 1b; Appendix S5). The distribution of alleles and supertypes did not reveal any evident geographic pattern (Figure 1 (Table 1). Particularly, the *90 allele prevailed in MEL, likely due to the impact of massive stocking, while the LDH composition in SCR and CAR revealed the nearly pure native origin of these populations.

| Comparison between STR and MHC population diversity and structuring
Linear models revealed a strong, positive relationship between neutral diversity (STR richness) and MHC diversity (both MHC-DAB richness and MHC-ST richness) ( Table 2) (Table 3). Variation within populations accounted for the largest fraction of variance (>74%), as expected for highly polymorphic markers.

| Relationship between body condition and MHC
The relationship between multiple components of MHC diversity and body condition was tested through mixed-effect models (Table 4), which showed that (a) the within-population frequency of MHC alleles and supertypes carried by an individual was negatively associated with its body condition (model 1, p < .01), while zygosity was unrelated with body condition for either MHC-DAB and MHC-ST (model 1, p = .13 and 0.18, respectively); and (b) nucleotide and amino acid dissimilarity between MHC sequences carried by an individual was not a significant predictor of body condition (model 2, p = .10 and 0.11, respectively). The effect of other genetic predictors included in the models (i.e., LDH type and STR heterozygosity) was statistically nonsignificant (Table 4).
Mixed-effect models were also used to explore the relationship between specific MHC supertypes and body condition-this was not doable for MHC alleles, as only three alleles occurred with a reasonable sample size (>15 individuals). Among the five models with AICc lower than the baseline model, only those including supertype-7, supertype-5, and supertype-11 were substantially better supported (ΔAIC < −2; Table 5). Beta coefficients of the three top models were significant, with a positive effect (i.e., better body condition) for individuals carrying supertype-7, and a negative effect (i.e., worse body condition) for those carrying supertype-5 and supertype-11.

| D ISCUSS I ON
We characterized adaptive (MHC) and neutral (STR) variation in six brown trout populations from central Italy, which are/were exposed to introgression from domestic lineages. Our ultimate goal was to assess the contribution of selection versus neutral forces to drive MHC variation and testing for the occurrence of balancing selection in its forms. We also evaluated the impact of anthropogenic introgression on MHC diversity.

| Characterization of MHC polymorphism and population variability
The analysis of MHC-DAB polymorphism confirmed typical features of expressed MHC genes, most of which were previously reported in studies targeting brown trout (Campos et al., 2006;Monzón-Argüello et al., 2014;Schenekar & Weiss, 2017): high levels of sequence variability (e.g., 48.2% of polymorphic amino acid sites), lack of stop codons, and the strong signature of historical positive selection for amino acid replacement, especially in some putative antigen-binding sites. Identified PSSs were consistent with PSS detected in previous studies focusing on brown trout (Monzón-Argüello et al., 2014;Schenekar & Weiss, 2017) and significantly matched with corresponding human MHC antigen-binding sites (p < .001, Fisher exact test) described in Brown et al. (1993), hence supporting their robustness. Multiple circular patterns (splits) connecting nodes in the network of MHC-DAB alleles (Figure 2) indicate, at least in part, reticulate evolution by micro-recombination (Klitz et al., 2012)-we indeed detected recombination in our sequences-even if they may also imply conflicting phylogenetic signals (Makarenkov et al., 2006). Note: The estimates of intercept and coefficients for predictors, the coefficient of determination (adjusted R 2 ), and the number of observations (no. of observations) are given for the two models. Note that, since we applied a sequential regression approach, the coefficient estimates for LDH-admixture indicate the relative effect that admixture has once accounted for the positive effect of STR richness-thus negative values indicate a relative effect, while the overall effect (i.e., without accounting for STR richness) is positive. Standard errors are reported in parentheses; statistical significance (*p < .05; **p < .01; and ***p < .001) is referred to t tests and F test for intercept/predictors and R 2 , respectively.

| Neutral evolutionary forces mostly shaped MHC population structure and variation
We found that the pattern of MHC-DAB variability agreed with neutral expectations: STR population diversity strongly correlated with MHC diversity as expressed by both alleles and functional supertypes, despite the limited number of observations (Table 2); MHC-DAB population structure substantially matched that of STR loci. Both these pieces of evidence indicate the importance of neutral processes in driving spatial variation of MHC (Garrigan & Hedrick, 2003;Spurgin & Richardson, 2010). Among nonselective forces, we may argue that genetic drift and population demographic history have played a major role to determine the observed differentiation and diversity at the MHC, explaining, at least in part, the among-population differences in MHC allele frequencies and the occurrence of private alleles. Conversely, migration provided a minor contribution since little gene flow was observed among the examined populations at the investigated spatial scale (results of the STRUCTURE analysis).
The comparative analyses of population structure at MHC-DAB and neutral loci (outlier analyses, Mantel test) failed to detect strong signals of (pathogen-mediated) balancing or diversifying selection on the MHC-DAB, possibly indicating their weak effect in investigated brown trout populations in the recent past (Garrigan & Hedrick, 2003)-but see the discussion section below for balancing selection on the current timescale. Particularly, we did not find signatures of MHC local adaptation, namely a stronger genetic structure for the MHC compared to neutral expectations (STR), which is quite unexpected since pathogens are likely to differ among ecologically heterogeneous biotopes as our study locations (see Table 1, e.g., variation in elevation, but note that pathogen composition information is not available for examined populations).
The analysis of functional structuring did not provide evidence for selection either: randomizations showed that the observed supertype-based population structure was not substantially different from random expectations; therefore, the lack of correlation between population differentiation at MHC-ST and STR should be interpreted as an artifact rather than a biological phenomenon.
These results also imply that supertypes and MHC alleles "behave" similarly, hence not supporting selection acting differently for

| Current balancing selection on MHC: The signature of RAA rather than DAA and HA
We also investigated the occurrence of balancing selection in its various forms on a recent timescale, by testing the association between multiple measures of MHC diversity and fish body condition. This latter is a measure of the general condition of a fish (fitness-related trait), and it was also found to be negatively related to parasite load in freshwater fish Seifertová et al., 2016). Thus, it was assumed here as a simple, noninvasive proxy of the parasite load of trout.
After controlling for genetic background (i.e., STR heterozygosity and LDH type), we found the individuals carrying locally rare MHC alleles/supertypes being in better body condition, irrespective of the MHC zygosity or the dissimilarity between MHC variants carried ( Note: MHC diversity is measured as: within-population frequency of alleles (MHC-DAB frequency index) and supertypes (MHC-ST frequency index) carried by an individual and its zygosity (based on alleles, MHC-DAB zygosity, and supertypes, MHC-ST zygosity) in model 1; MHC dissimilarity based on nucleotide (Kimura two-parameter) and amino acid (Poisson-corrected) distances in model 2. Estimates for each model version are provided along with standard errors in parentheses (SE), degree of freedom (df), t-statistics (t), associated p-values (p), marginal R-squared (R 2 m, i.e., the variance explained by only fixed effects), and conditional R-squared (R 2 c, i.e., the variance explained by both fixed and random effects). Reference level for factor LDH type is the "domestic" genotype 90/90. Variables with a significant effect (p < .05) are in bold.
allowed disentangling the relative contribution of RAA and HA, which is notoriously challenging because sufficiently rare variants favored under RAA mostly occur in heterozygotes Spurgin & Richardson, 2010).
In brown trout, the analysis of the relationship between specific MHC functional variants and body condition revealed that individuals carrying the supertype-7 were generally in better body condition, suggesting that possessing this functional variant may confer a fitness advantage, possibly an improved resistance against pathogens.
The opposite could be hypothesized for individuals carrying the supertype-5 or supertype-11, although the latter likely reflects local rather than general effects as the vast majority of individuals carrying supertype-11 occurred in a single population (MEL). Positive and/or negative correlations between fitness-related traits and certain MHC alleles (reviewed in Sommer, 2005), as well as supertypes (Schwensow et al., 2007;Sepil et al., 2012Sepil et al., , 2013, are usually believed evidence of RAA, reflecting a snapshot of the MHC-pathogen cyclical interactions within a population (Spurgin & Richardson, 2010).
This could be the case of supertype-11 in our study. On the other hand, the occurrence of associations in multiple unrelated populations simultaneously (i.e., supertype-5 and supertype-7) is not entirely consistent with RAA, since MHC-pathogen dynamics under RAA should be spatially unsynchronized. However, they may be consistent even with other mechanisms of pathogen-mediated selection, such as fluctuating selection Spurgin & Richardson, 2010). In addition, RAA does not necessarily exclude interactions between specific MHC molecules and pathogens (Trachtenberg et al., 2003

| Effects of human-mediated introgression on MHC diversity and conservation remarks
Introgression from domestic lineages was detected in all examined populations, as indicated by the various frequencies of hybrids at the LDH gene. Noticeably, MHC-DAB richness (as well as STR richness) increased with increasing LDH-based admixture, thus supporting introgression as a non-negligible source of adaptive diversity in wild populations (Hedrick, 2013). Because adaptive genes undergoing balancing selection (RAA) are predicted to permeate into the gene pool of sink populations more readily than neutral ones (Schierup et al., 2000), introgression, even if limited, should be more pronounced for MHC than neutral markers (Dudek et al., 2019).
Such a marker-dependent pattern qualitatively emerged from the breakdown of AMOVAs (Table 3): the difference between groups of relatively pure and massively stocked populations was substantial although marginally significant for STR (9%, p = .063), while reduced for MHC-DAB (2%, p = .069) and negligible for MHC-ST (< 1%, p = .332). Therefore, we may speculate that the relatively recent hybridization together with RAA could have promoted the rapid introgression and the maintenance of hatchery-introduced, initially rare MHC (functional) variants into native populations, including those only marginally affected by stocking (i.e., SCR and CAR). This may have eventually resulted in an apparently reduced differentiation among roughly pure native and massively stocked populations for MHC, but not for STR. Such hypothesis, however, deserves further investigations.
In a conservation perspective, genetic introgression from hatchery-reared, non-autochthonous fish may entail a decreased viability of wild native populations (Bert et al., 2007;Pinter et al., 2019). In the case of adaptive genes, such as MHC, variation may have undergone intentional or unintentional domestication selection during the propagation in captivity, possibly implying nonadaptive evolution due to relaxed pathogen-mediated selection or, conversely, adaptation induced by a different parasite pressure.
According to the latter hypothesis, stronger signals of selection on MHC were observed in hatchery versus wild populations of brown trout in Austria (Schenekar & Weiss, 2017). On the other hand, the maintenance of a broader MHC repertoire within populations may contribute to coping more effectively with novel pathogens, possibly including those introduced by released domestic trout-anyway, if experimentally supported, this would require a careful cost-benefit evaluation. Whatever the actual consequences of anthropogenicintroduced MHC diversity on brown trout population viability, our study emphasizes the general value of MHC for conservation genetics of wild-ranging animals. MHC variation, particularly that expressed by functional supertypes, can be highly informative to assess the evolutionary potential and viability of populations because of its relationship with fitness-related traits (Sommer, 2005;Ujvari & Belov, 2011). Also, MHC may provide a valuable supportive tool for defining conservation and management units (e.g., Zhu et al., 2013).
Yet, an integrated conservation approach should take into account both neutral (classical markers) and adaptive (MHC) variation (Sommer, 2005;Ujvari & Belov, 2011), but this is still rarely adopted, unfortunately.

| CON CLUS ION
Our study contributes to understanding the evolution of MHC genes in wild populations: we provide evidence for the interplay between neutral (genetic drift, demography) and selective (balancing selection) processes to shape the pattern of adaptive (MHC) diversity in brown trout, disentangling between balancing selection mechanisms actually involved in maintaining MHC (functional) variation. We also demonstrated how human-mediated introgression with domestic trout coupled with balancing selection may affect patterns of adaptive diversity. Finally, our findings support the power and the biological relevance of the supertyping approach to investigate selection and accurately analyze the functional variation of highly polymorphic genes of MHC (Naugler & Liwski, 2008;Schwensow et al., 2007;Sepil et al., 2012Sepil et al., , 2013Sette & Sidney, 1998;Trachtenberg et al., 2003).

ACK N OWLED G M ENTS
We thank P.T. Colombari and G. Moccia for their kind help in fieldwork. We are grateful to the Associate Editor O. Manlik, J. Radwan, and three other anonymous reviewers for constructive comments and useful suggestions that contributed to improving the quality of this paper.

CO N FLI C T O F I NTE R E S T
None declared. Investigation (supporting). Lorenzo Tancioni: Conceptualization (lead); funding acquisition (lead); investigation (equal); resources (equal); writing-review and editing (supporting).

DATA AVA I L A B I L I T Y S TAT E M E N T
The following data have been deposited in Dryad (https://doi. org/10.5061/dryad.xsj3t x9fg): paired-end Fastq files and individual tag for demultiplexing; individual genetic and phenotypic data; the