Natural selection on synonymous genetic variation in the major histocompatibility complex

Synonymous Variation, Codon Usage, Genetics, Major Histocompatibility Complex, MHC, Molecular Evolution, Natural Selection, Birds, Disease Biology , Host -Pathogen Interactions

Protein coding DNA sequences harbor synonymous nucleotide variation that does not change amino acid sequences but influences phenotypes via multiple effects on the pathway from gene to protein. Synonymous variation has recently been shown to coevolve between viruses and their natural hosts, but its potential role in host immune defenses has not been explored. Here, I present evidence of natural selection on synonymous variation in the major histocompatibility complex (MHC), a highly polymorphic multigene locus that plays a crucial role in pathogen recognition by the adaptive immune system of vertebrate species. Using data from a wild population of Great Reed Warblers, I show that codon usage in exon 3 of MHC class I (MHC -I) genes is under strong purifying selection in 56 out of 87 codon sites. Scanning the Great Reed Warbler genome for tRNA genes revealed that, for most amino acids, bias towards preferred codons was associated with abundances of tRNA isotypes, indicating that the purifying selection is likely driven by selection for increased translational efficiency. However, spikes of synonymous variation appeared in 31 of the 87 sites in the MHC -I exon 3, and in those sites, codon usage bias and correlations with tRNA abundances were reduced. The distribution of the spikes of synonymous variation showed no consistent association with structural domains of the MHC -I protein, nor with s ites under positive selection for amino acid change , which are considered important for antigen binding properties. Intriguingly, the amount of synonymous variation in genotypes showed a positive correlation with Darwinian fitness, indicating that important evolutionary forces are at play that neutralize purifying selection in the 31 sites. From an ultimate perspective, the release of purifying selection among certain sites in MHC genes may indicate an arms race with pathogens, and I propose that the spikes of synonymous variation may reveal a footprint of natural selection on MHC genes to escape inhibitory molecular interactions between intracellular pathogens and MHC mRNA. Unravelling the mechanisms of such interactions should be of great importance to our understanding of this extremely important locus and I hope that the results and methodological advancements presented here will spark future studies of synonymous variation in the MHC and its biological effects. (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted February 23, 2026. ; https://doi.org/10.64898/2026.02.23.707394doi: bioRxiv preprint 2

A much-overlooked type of natural selection is that which concerns synonymous nucleotide variation that does not change amino acid sequence s. Approximately one quarter to one third of single nucleotide mutations are expected to be synonymous and have been presumed to have little or no effects on phenotypes (1–3). However, it is generally acknowledged that synonymous nucleotide variation can cause changes in protein expression, conformation, and function, and may be subject to both purifying and positive selection (3–7). Codon usage bias (CUB) is the common term describing the unequal use of synonymous codons in coding genomic regions (3, 5, 7). CUB is prevalent in highly expressed genes especially in species with large population sizes and is thought to be driven by optimization of codon usage to the relative concentrations of tRNA isotypes . Preferred codons are favored by natural selection because they speed up mRNA translation by ribosomes and reduce the risk of mistranslation, where incorporation of a different amino acid could introduce potentially detrimental changes to the protein (3, 8, 9). However, it has been argued that this view is simplistic, and that optimization of codon usage for fast gene expression is balanced by other factors such as correct folding of the protein, especially in higher organisms . In particular , co-translational folding of the emerging protein is thought to be regulated by codon usage, with slow folding β-sheets and interdomain bridges harboring more rare codons than fast folding α-helices (3, 5, 6, 10, 11). Codon usage can facilitate co-translational protein folding by slowing the progression of the ribosome to give time for a nascent protein domain to fold unaffected by downstream elements. Such slowing of the ribosome can be caused by codons with rare tRNAs or formation of local stable structures of the mRNA (3, 5, 6, 11). Furthermore, synonymous variation in specific sites along coding sequences may be preserved by molecular interactions e.g. , of DNA with transcription factors or of mRNA with microRNAs, splic ing enhancers, deaminases, ribosomes, and other molecules that regulate transcription and translation (3, 6, 7, 12, 13). Following the nearly neutral theory, it has been predicted that synonymous SNPs are unlikely to be affected by natural selection in species that have small population sizes, such as many vertebrate species (13, 14). However, that prediction rests on the assumption that synonymous nucleotide variation is mostly under weak selection , i.e. that the selection coefficient is less than half the effective population size: |s| < 1/(2Ne). While that appears to be the general trend observed across genomes, e.g. as shown in a comparison of 5,639 orthologous genes between Human, Chimpanzee, Macaque, Mouse, and Rat (4), some loci evolve under abnormal regimes of natural selection. Numerous disease associations involving synonymous SNPs have been detected and a recent study estimated that 25.9% of synonymous SNPs are under weak and 3.6% under strong negative selection in humans (6, 15, 16). A particularly interesting candidate locus for effects of synonymous nucleotide variation in vertebrate species is the major histocompatibility complex (MHC). The genes of the MHC encode molecules that present peptide antigens to T -cells, which is a decisive step in the induction of adaptive immune responses (17, 18). MHC genes are engaged in a coevolutionary arms race with pathogens, that has generated both a great diversity in these (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted February 23, 2026. ; https://doi.org/10.64898/2026.02.23.707394doi: bioRxiv preprint 3 genes and preservation of polymorphisms over evolutionary time (17, 19–23). Meanwhile, MHC genes have repeatedly been shown to also be subject to sexual selection, which can drive sexually antagonistic selection over MHC diversity (24–28). Among MHC studies, a broad consensus has prevailed to regard synonymous nucleotide variation as silent and to regard variation in protein sequences and gene expression as the only factors determining phenotypic variation. However, biological effects of synonymous nucleotide variation in the MHC have never been thoroughly investigated, leaving a major gap in our knowledge about a substantial body of genetic variation in these genes that serve a central function in vertebrate immune systems. Protein sequence variation in the MHC has been shown to affect the binding of antigen peptides and to be under positive selection, and it is evident that this variation is of evolutionary importance (22, 23, 29, 30). However, the assumption that synonymous variation in the MHC serves no biological function involves a severe risk of accumulating wrong or biased insights. A recent study observed that CUB in viruses tended to be more similar to that of symptomatic hosts than that of natural hosts (i.e., hosts that the virus coevolved with). Using an experiment in yeast, the authors showed that expression of exogenous genes with CUB similar to the host severely impeded translation of host genes, supporting a general deleterious effect of CUB similarity between viruses and hosts (31). This previously unrecognized complexity in coevolution between viruses and their hosts is likely a driving selective force behind the observed dissimilarity in CUB between viruses and their natural hosts , and it invites the question of how antagonistic coevolution of codon usage between viruses and hosts may have affected MHC genes, the expression of which is crucial to adaptive immune responses of vertebrate hosts . If hosts become more symptomatic by regulation of MHC expression by infecting viruses, then codon usage may be subject to an arms race similar to that between pathogen antigen epitopes and binding repertoires of host MHCs, but the aim of the host would be to escape regulatory mechanisms of the viruses to ensure gene expression. Here I characterized for the first time synonymous nucleotide variation in the MHC of a wild vertebrate species and investigated whether it is subject to natural selection. I used data from a 20- year study of a wild breeding population of Great Reed Warblers Acrocephalus arundinaceus, a socially polygynous, migratory songbird. I employed 559 genotypes of the highly polymorphic MHC class I (MHC-I) exon 3, which encodes a major section of the antigen binding groove in MHC-I molecules (32), along with the Great Reed Warbler genome assembly and field observations of individual life histories and lifetime reproductive success available from previous publications (25, 29, 33, 34). This detailed ecological data set provides a unique opportunity to combine analyses of synonymous nucleotide variation in the MHC with investigations of direct effects on individual Darwinian fitness. Specifically, I tested whether s ynonymous codon usage among 390 Great Reed Warbler MHC-I exon 3 sequences showed footprints of natural selection and analyzed whether synonymous nucleotide variation differs between sites along the exon. I test ed site specific selection on synonymous variation by a novel approach that compared the 390 empirical MHC-I exon 3 (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted February 23, 2026. ; https://doi.org/10.64898/2026.02.23.707394doi: bioRxiv preprint 4 sequences with 1,000 synonymous data sets simulated under the assumption of no selection on codon usage. Furthermore, I compared site-specific synonymous variation with sites predicted to be under positive selection for amino acid changes and with predicted structural domains of the folded MHC-I protein. Furthermore, I quantified CUB and analyzed whether codon usage in the MHC-I is correlated with frequencies of cognate tRNA isotypes in the Great Reed Warbler genome, a proxy for tRNA availability (10). Finally, I calculated the mean number of synonymous nucleotide changes among the MHC-I sequences in each individual and analyzed whether synonymous nucleotide variation in the Great Reed Warbler MHC-I was associated with life span or Darwinian fitness.

Selection analyses I employed mutation-selection models in codeml (35, 36) to test the null hypothesis that variation in codon usage in exon 3 of MHC -I genes (Table 1) is shaped by mutation bias alone and not by selection on synonymous nucleotide variation, cf. (4). I have previously demonstrated that the Great Reed Warbler MHC-I exon 3 sequences show evidence of positive selection for amino acid changes (29), and I therefore nested mutation-selection and null-models within the site model M8, which allows a category of sites evolving under positive selection (Ω > 1) (4, 37). The mutation- selection model (FMutSel), reflecting the alternative hypothesis that codon usage is affected by selection, fit the data significantly better than the null-model (FMutSel0), which assumes that codon usage is only affected by mutation bias (likelihood ratio test: p = 2.43e-15; Table 2). Notably, the null-model overestimated both the proportion of sites with Ω > 1 and the Ω estimate for that category compared to the mutation-selection model (Table 2). This indicates that synonymous nucleotide variation is not only subject to natural selection itself, but taking selection on synonymous variation into account a lso affects analyses of selection on amino acid variation in codeml. Repeating the analysis using the site model M2a confirmed the results from the M8 model (Table S1). Sites with Ω > 1 inferred by Bayes Empirical Bayes analysis in the M8 mutation - selection model are specified in Table S2 (38). Patterns of synonymous variation To further investigate the patterns of synonymous nucleotide variation, I quantified the synonymous changes per base and per codon site in pairwise comparisons of the 390 MHC-I exon 3 sequences. The frequency of synonymous codon variation differed greatly along the sequences, with most sites showing little or no synonymous variation, but with prominent spikes appearing in roughly one third of the sites (N=31) along the sequence (Fig. 1 a). The synonymous codon variation was mostly associated with variation in the third nucleotide position, except in site 61, where more than half of the observed synonymous variation involved a change in the first nucleotide position (Fig. 1a; Table S3). That was also the case for site 33, but that site harbored almost no synonymous variation (freq. ~0.005) . There was no consistent association between frequencies of synonymous codon variation and the location of encoded amino acids in bridge (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted February 23, 2026. ; https://doi.org/10.64898/2026.02.23.707394doi: bioRxiv preprint 5 domains, β-sheets, or α-helices of the folded protein structure predicted with AlphaFold 3 (39) (Fig. 1a; Fig. S1; Table S4). Sites that were inferred to evolve under positive selection for amino acid change (i.e., sites with Ω > 1 marked by asterisks in Fig. 1a) showed no obvious trend towards harboring high or low synonymous variation. Table 1 Codon usage counts among the 390 Great Reed Warbler MHC-I exon 3 sequences. Phe F TTT 350 Ser S TCT 263 Tyr Y TAT 242 Cys C TGT 394 TTC 986 TCC 1210 TAC 1601 TGC 392 Leu L TTA 2 TCA 1 Stop X TAA 0 Stop X TGA 0 TTG 480 TCG 0 TAG 0 Trp W TGG 1071 Leu L CTT 96 Pro P CCT 120 His H CAT 467 Arg R CGT 434 CTC 756 CCC 3 CAC 635 CGC 722 CTA 9 CCA 27 Gln Q CAA 14 CGA 256 CTG 1548 CCG 373 CAG 913 CGG 901 Ile I ATT 15 Thr T ACT 52 Asn N AAT 501 Ser S AGT 176 ATC 1170 ACC 723 AAC 33 AGC 674 ATA 35 ACA 173 Lys K AAA 259 Arg R AGA 853 Met M ATG 116 ACG 225 AAG 913 AGG 945 Val V GTT 448 Ala A GCT 1126 Asp D GAT 611 Gly G GGT 104 GTC 746 GCC 362 GAC 1428 GGC 775 GTA 1 GCA 27 Glu E GAA 1352 GGA 660 GTG 635 GCG 273 GAG 2379 GGG 1874 To understand how selection has shaped synonymous nucleotide variation along the MHC-I exon 3 sequences, I simulated 1,000 data sets of 390 nucleotide sequences synonymous to the empirical Great Reed Warbler MHC-I exon 3 sequences (i.e., with identical amino acid sequences), assuming the null hypothesis that variation in codon usage is due to mutation bias alone and not affected by selection acting at synonymous codons. The frequencies of synonymous variation a mong the sequences in the simulated data sets highlight the n on-random distribution of synonymous variation among the empirical MHC-I sequences (Fig. 1b ). 86 out of the 87 sites harbored significantly less synonymous differences than predicted from the sequences simulated under the Table 2 Summary of the M8 mutation-selection vs. null models in codeml. ln(L) specifies the log likelihood of the model s, Ω > 1 specifies the proportion of sites that were estimated to be under positive selection for amino acid change, and Ω est indicates the value of Ω estimated for that category of sites. P-value was calculated by likelihood ratio test of the nested FMutSel vs. FMutSel0 models. Site model Codon model Hypothesis ln(L) Ω > 1 Ω est P(H0) M8 FMutSel0 H0: No selection on synonymous variation -7521.6 21.2% 4.22 2.43e-15 FMutSel H1: Selection on synonymous variation -7443.6 15.5% 3.23 (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted February 23, 2026. ; https://doi.org/10.64898/2026.02.23.707394doi: bioRxiv preprint 6 Fig. 1 Frequency of synonymous changes for each amino acid site among the 390 Great Reed Warbler MHC-I exon 3 sequences. (a) Gray bars indicate the total frequency of synonymous variation for each site, with red bars indicating synonymous changes in first nucleotide position and blue bars synonymous changes in third nucleotide position. (b) Gray bars indicate the frequency of synonymous changes for each site in the empirical data set (same as the gray bars in a). Black diamonds indicate the mean frequency of synonymous changes for each site among 1,000 data sets simulated under the assumption of no selection, with error bars indicating 95% C.I. The consensus amino acid sequence indicates the most frequent amino acid for each site. Asterisks mark the sites predicted to evolve under positive selection for amino acid changes by B ayes Empirical Bayes analysis in the M8 FMutSel model. Shaded areas in (a) mark bridges between β- sheet and α-helix domains in the predicted protein structure. Frequency of synonymous changes in pairwise comparisons of MHC−I sequences Site Frequency 0.0 0.1 0.2 0.3 0.4 0.5 α α αβ β β β 0.0 0.1 0.2 0.3 0.4 0.50.0 0.1 0.2 0.3 0.4 0.5 WLRVYGCEL LSDGSVRGSYRFGYDGRDF I SFDLESGRFVAADSAAE I TRRRWEHEGTVAERWTNYLKHECPEWLQRHVRYGQKELER 1 3 5 7 9 11 13 15 17 19 21 23 25 27 29 31 33 35 37 39 41 43 45 47 49 51 53 55 57 59 61 63 65 67 69 71 73 75 77 79 81 83 85 871 3 5 7 9 11 13 15 17 19 21 23 25 27 29 31 33 35 37 39 41 43 45 47 49 51 53 55 57 59 61 63 65 67 69 71 73 75 77 79 81 83 85 0 0 Site total 1st codon position 3rd codon position Omega > 1 Frequency of synonymous changes in pairwise comparisons of MHC−I sequences Site Frequency 0.0 0.2 0.4 0.6 0.8 WLRVYGCE L L SDGSVRGSYRFGYDGRDF I SFDL ESGRFVAADSAAE I TRRRWEHEGTVAERWTNY L KHECPEWLQRHVRYGQKE L ER 1 3 5 7 9 11 13 15 17 19 21 23 25 27 29 31 33 35 37 39 41 43 45 47 49 51 53 55 57 59 61 63 65 67 69 71 73 75 77 79 81 83 85 87 0 0 Empirical frequencies Simulated mean, 95% C.I. Omega > 1 a b (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted February 23, 2026. ; https://doi.org/10.64898/2026.02.23.707394doi: bioRxiv preprint 7 assumption of no selection , providing strong e vidence for purifying selection on codon usage in the Great Reed Warbler MHC-I (Table S5). However, the spikes of synonymous variation along the sequences indicate that purifying selection is relaxed in about one third of the sites (1, 3, 5, 11, 12, 14, 17, 18, 23, 24, 35, 38, 40, 41, 53, 54, 57, 61, 63, 66, 68, 69, 71, 72, 74, 77, 78, 80, and 85). Site 31 even shows evidence of selection for increased synonymous variation (P < 0.001; Fig. 1b; Table S5; Table S6). Codon usage bias For most degenerate amino acids, CUB was greater in the 390 Great Reed Warbler MHC-I exon 3 sequences than in the 1,000 data sets simulated under the assumption of no selection on synonymous codon usage (Fig. 2; Table S7). The exceptions are histidine (His), glutamic acid (Glu), cysteine (Cys), and arginine (Arg), where the observed CUBs were smaller than expected, and lysine (Lys) and aspartic acid (Asp), where the observed CUB s were similar to those in the simulated data sets. CUB was on average greater among the sites where synonymous variation is under purifying selection compared to the sites where selection is relaxed (paired t-test, t = 2.22, d.f. = 16, p = 0.021). However, for Leu, Val, Asp, Arg, and glycine (Gly) it was smaller (Fig. 3; Table S8). Fig. 2 Codon usage bias (CUB) measured as the variance in relative synonymous codon usage for degenerate amino acids . Gray bars show the CUB s observed among the 390 empirical MHC-I exon 3 sequences . Black diamonds show the mean CUB s observed among 1,000 data sets simulated under the assumption of no selection on synonymous codon usage , with error bars indicating 95% C.I. Phe Leu Ile Val Ser Pro Thr Ala Tyr His Gln Asn Lys Asp Glu Cys Arg Gly Amino acid Codon usage bias 0.0 0.5 1.0 1.5 2.0 2.5 0 0 Empirical CUB Simulated mean, 95% C.I. (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted February 23, 2026. ; https://doi.org/10.64898/2026.02.23.707394doi: bioRxiv preprint 8 Fig. 3 Codon usage bias (CUB) measured as the variance in relative synonymous codon usage for degenerate amino acids among the 390 empirical MHC-I exon 3 sequences. Black bars show the CUBs measured in the 56 sites that evolve under purifying selection for codon usage. White bars show the CUBs measured in the 31 sites where purifying selection for codon usage was relaxed. Codon usage and abundance of tRNA isotypes To estimate tRNA isotype availability, I employed tRNAscan -SE (40) to scan the Great Reed Warbler genome for tRNA genes. The analysis detected 325 functional genes for cytosolic standard amino acid tRNAs. The counts of tRNA isotypes among the 325 genes are shown in Table 3. Fourteen codons encoding 13 amino acids were not matched by a tRNA isotype in the genome. Twelve amino acids were matched by more than one tRNA isotype in the Great Reed Warbler genome. For 7 of those amino acids ( Leu, Val, serine (Ser), alanine (Ala), glutamine (Gln), Lys, and Glu), relative synonymous codon usage (RSCU) in the 390 MHC -I exon 3 sequences was positively correlated with normalized relative frequencies of synonymous tRNA isotypes ( Fig. 4; Table S9; S10). For Leu, Val, Ser, and Ala the correlation coefficients were significantly more positive than expected from comparisons with sequences simulated under the assumption of no selection, while for Gln, Lys, and Glu the correlation coefficients were 1 with both the empirical and simulated data sets. When comparing the 56 sites that evolve under purifying selection for synonymous codon usage with the 31 sites where purifying selection for codon usage was relaxed , a different image appeared. Overall, the correlation coefficients were significantly greater among the 56 sites that evolved under purifying selection than among the 31 sites where selection was relaxed (paired t- Phe Leu Ile Val Ser Pro Thr Ala Tyr His Gln Asn Lys Asp Glu Arg Gly Amino acid Codon usage bias 0.0 0.5 1.0 1.5 2.0 2.5 3.0 3.5 0 0 Purifying selection Relaxed selection (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted February 23, 2026. ; https://doi.org/10.64898/2026.02.23.707394doi: bioRxiv preprint 9 Table 3 Isotype counts among the 325 tRNA genes in the Great Reed Warbler genome. Anticodons are specified in italics. Phe F TTT AAA 0 Ser S TCT AGA 12 Tyr Y TAT ATA 0 Cys C TGT ACA 0 TTC GAA 7 TCC GGA 0 TAC GTA 9 TGC GCA 28 Leu L TTA TAA 1 TCA TGA 3 Stop X TAA TTA 0 Stop X TGA TCA 0 TTG CAA 2 TCG CGA 6 TAG CTA 0 Trp W TGG CCA 5 Leu L CTT AAG 6 Pro P CCT AGG 7 His H CAT ATG 0 Arg R CGT ACG 6 CTC GAG 0 CCC GGG 0 CAC GTG 7 CGC GCG 0 CTA TAG 2 CCA TGG 3 Gln Q CAA TTG 3 CGA TCG 2 CTG CAG 5 CCG CGG 3 CAG CTG 12 CGG CCG 1 Ile I ATT AAT 46 Thr T ACT AGT 6 Asn N AAT ATT 0 Ser S AGT ACT 0 ATC GAT 12 ACC GGT 0 AAC GTT 7 AGC GCT 5 ATA TAT 2 ACA TGT 3 Lys K AAA TTT 6 Arg R AGA TCT 3 Met M ATG CAT 15 ACG CGT 2 AAG CTT 10 AGG CCT 3 Val V GTT AAC 4 Ala A GCT AGC 15 Asp D GAT ATC 0 Gly G GGT ACC 0 GTC GAC 1 GCC GGC 0 GAC GTC 7 GGC GCC 12 GTA TAC 1 GCA TGC 8 Glu E GAA TTC 4 GGA TCC 7 GTG CAC 4 GCG CGC 3 GAG CTC 5 GGG CCC 4 (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted February 23, 2026. ; https://doi.org/10.64898/2026.02.23.707394doi: bioRxiv preprint 10 test; t = 2.05, d .f. = 11, p = 0.032 , mean diff. = 0.63; Fig. 5; Table S11). The amino acids Val, threonine (Thr), Ala, Lys, and Glu showed strong positive correlations among the 56 sites thatevolve under purifying selection for codon usage, but negative correlations among the 31 sites where purifying selection for codon usage was relaxed. A contrasting pattern was also seen for Ser, where the correlation was negative among the 56 sites and positive among the 31 sites. Leu, isoleucine (Ile), proline (Pro), Arg, and Gly showed correlations in similar directions between the 56 sites and the 31 sites, although with different strength. For Gln the correlation coefficients were 1 among both the 56 sites and the 31 sites. Together, these observations suggest that synonymous codon usage in the Great Reed Warbler MHC-I is affected by availability of cognate tRNAs during translation , and that the associations differ both between amino acids and between sites. Fig. 4 The relationship between normalized relative frequencies of synonymous tRNA isotypes in the Great Reed Warbler genome and relative synonymous codon usage (RSCU). Gray bars show Pearson correlation coefficients for RSCU in the 390 empirical MHC -I exon 3 sequences. Black diamonds show the mean Pearson correlation coefficients for RSCU among 1,000 data sets simulated under the assumption of no selection on synonymou s codon usage, with error bars indicating 95% C.I. The numbers under the bars specify the number of codons matched by tRNA isotypes for each amino acid. Only amino acids with more than one tRNA isotype in t he genome were included in the analyses. tRNA−codon usage correlations Amino acid Pearson correlation coefficient −1.0 −0.5 0.0 0.5 1.0 Leu Ile Val Ser Pro Thr Ala Gln Lys Glu Arg Gly 5 3 4 4 3 3 3 2 2 2 5 3 0 0 Empirical data Simulated mean, 95% C.I. (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted February 23, 2026. ; https://doi.org/10.64898/2026.02.23.707394doi: bioRxiv preprint 11 Fig. 5 The relationship between normalized relative frequencies of synonymous tRNA isotypes in the Great Reed Warbler genome and relative synonymous codon usage (RSCU). Black bars show Pearson correlation coefficients for RSCU calculated among the 56 sites of the empirical 390 MHC-I exon 3 sequences that evolve under purifying selection for codon usage. White bars show Pearson correlation coefficients for RSCU calculated among the 31 sites where selection for codon usage was relaxed. The numbers under the bars specify the number of codons matched by tRNA isotypes for each amino acid. Only amino acids with more than one tRNA isotype in the genome were included in the analyses. Fitness effects To investigate whether synonymous nucleotide variation in the MHC-I has a measurable effect on the survival or fitness of individuals, I analyzed previously published life history data from a long- term study of Great Reed Warblers at lake Kvismare in Sweden (25). For statistical analyses of fitness effects, I calculated the mean number of synonymous nucleotide changes in pairwise comparisons among the MHC-I exon 3 sequences in each individual. The mean number of synonymous nucleotide changes in each individual ranged from 4.90 to 8.48 and followed a normal distribution with mean = 6.69 and s.d. = 0.64 (Fig. S2). Among male Great Reed Warblers, the mean number of synonymous nucleotide changes in each individual was positively associated with life span (glm, b = 0.18, t = 2.21, p = 0.030, d.f. = 76; Fig. 6a; Table S12) and lifetime number of fledged offspring (glm, b = 0.36, t = 2.33, p = 0.022, d.f. = 76; Fig. 6b; Table S13). However, no association was observed in females (life span: glm, b tRNA−codon usage correlations Amino acid Pearson correlation coefficient −1.0 −0.5 0.0 0.5 1.0 Leu Ile Val Ser Pro Thr Ala Gln Lys Glu Arg Gly 5 3 4 4 3 3 3 2 2 2 5 3 0 0 Purifying selection Relaxed selection (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted February 23, 2026. ; https://doi.org/10.64898/2026.02.23.707394doi: bioRxiv preprint 12 = -0.058, t = -0.61, p = 0.54, d.f. = 95; Fig. 6a; Table S12; lifetime number of fledg ed offspring: glm, b = -0.015, t = -0.13, p = 0.90, d.f. = 95; Fig. 6b; Table S13). As common in species with discrete breeding events, lifetime fitness was strongly associated with life span in the Great Reed Warbler. I therefore also analyzed the effects of the mean number of synonymous nucleotide changes on lifetime number of fledged offspring in a model that adjusted for individual life span. In this model, the mean number of synonymous nucleotide changes in each individual showed a trend towards a positive effect in both males and females (glm, b = 0.10, t = 1.77, p = 0.079, d.f. = 172; Fig. S3; Table S14). No effects were observed of the number of different MHC-I sequences or mean amino acid p- distance between the MHC -I sequenc es in each individual. These variables were included as covariates in the full models but removed by the model simplification step in each analysis. Results from models of the effects of number of different MHC-I sequenc es and mean amino acid p- distance on life span , lifetime number of fledg ed offspring, and lifetime number of fledg ed offspring adjusted for life span are presented in Tables S15–S17.

This study shows for the first time that MHC genes are subject to natural selection on synonymous nucleotide variation and that the synonymous variation has direct measurable effects on Darwinian fitness of individuals. Intriguingly, 56 sites in the Great Reed Warbler MHC-I exon 3 showed evidence of strong purifying selection on codon usage, while selection was relaxed in the remaining 31 sites (Fig. 1) . Both CUB and the correlation s between CUB and abundances of cognate tRNA isotypes in the genome were on average greater in the sites under strong purifying selection than in the sites where selection was relaxed , indicating that selection for translation speed and fidelity may generally restrict codon usage in the MHC (Fig. 3; Fig. 5) . The driving mechanisms behind such an effect may be optimization of codons to the availability of amino acyl tRNAs for translation and the advantage of recycling ribosomes (3, 7, 9, 41). Some other forces must be at play that allows selection to be relaxed in 31 sites of the MHC-I exon 3. A possible explanation is that unpreferred codons may temporarily slow translation, allowing time for the nascent protein to fold without interference from downstream elements (11), but the distribution of the sites with higher frequencies of synonymous variation showed no consistent association with predicted domains of the folded MHC-I protein (Fig. 1a). However, the observed patterns warrant further investigations, as three of the four interdomain bridges following β-sheets in the first half of the exon show ed elevated frequencies of synonymous variation. Temporary slowing of translation may not just be caused by use of rare tRNAs but may also be enforced by local stable mRNA structures that physically slow the progression of the ribosome. Synonymous nucleotide variation is known to affect mRNA structure and stability (5), but it is unclear how selection for local mRNA structures would affect its distribution in relation to the encoded protein domains. This is an interesting avenue for future investigations, but it is beyond the scope of the (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted February 23, 2026. ; https://doi.org/10.64898/2026.02.23.707394doi: bioRxiv preprint 13 Fig. 6 The relationship between synonymous nucleotide variation and (a) life span and (b) lifetime number of fledged offspring in adult Great Reed Warbler s. The lines show predictions from generalized linear model s that included sex and the interaction ‘mean number of synonymous nucleotide changes ´ sex’. Open circles and dashed lines indicate females. Black triangles and solid lines indicate males. Jitter was added to life span and lifetime number of fledgling s to distinguish discrete data points. Mean number of synonymous substitutions Life span 1 2 3 4 5 6 7 8 9 4.5 5 5.5 6 6.5 7 7.5 8 8.5 9 Mean number of synonymous substitutions Lifetime number of fledglings 1 11 21 31 41 51 4.5 5 5.5 6 6.5 7 7.5 8 8.5 9 a b (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted February 23, 2026. ; https://doi.org/10.64898/2026.02.23.707394doi: bioRxiv preprint 14 present study as proper modelling of the mRNA structure would require long range transcripts of complete MHC genes, which are challenging to obtain on a large scale with current sequencing technologies (42, 43). There was also no apparent association between the distribution of synonymous nucleotide variation and sites that were inferred to evolve under positive selection for amino acid change (i.e., sites with Ω > 1) . However, a recent study found that the Great Reed Warbler MHC-I binds antigen in a flat conformation due to a tweezer-like interaction between two Arg (R) that form a conserved restriction point within the peptide binding gr oove (44). Interestingly, both sites forming this restriction point show elevated synonymous variation (sites 3 and 61; Fig. 1a). The observed correlations with individual survival and fitness in the study population indicate that synonymous nucleotide variation in the MHC -I is of significant biological importance (Fig. 6). The fitness effects of synonymous variation were not accompanied by effects of number of different MHC-I sequences in each individual or amino acid distance between sequences, and the effects observed on life span and lifetime number of fledglings represent a stronger influence on individual fitness than the previously reported effects of number of different MHC-I sequences on offspring recruitment success (25). The positive correlation with life span suggests that the phenotypic advantage of high synonymous nucleotide variation in males is largely associated with survival, however the trend towards a positive effect on lifetime number of fledged offspring in a model that adjusted for life span suggested that both males and females with high synonymous variation may on average rear more offspring per breeding year (Fig. S3). From an ultimate perspective, the release of purifying selection among certain sites in MHC genes may be driven by an arms race with pathogens, that MHC genes are known to coevolve with (23, 45). Recent investigations uncovered a previously unrecognized antagonistic coevolution of codon usage between viruses and hosts (31). The authors showed that similarity of codon usage between viruses and hosts can be deleterious to the hosts as the expression of virus genes may impede translation of host genes. However, the patterns of synonymous variation observed in the Great Reed Warbler MHC-I exon 3 do not seem consistent with such antagonistic coevolution . Two thirds of the sites in the Great Reed Warbler MHC-I exon 3 show evidence of strong purifying selection on codon usage, and this does not suggest a coevolutionary response to avoid similarity with codon usage of viruses. If MHC expression would be affected by viruses matching host codon usage, then reducing the translation of two thirds of the sites would be efficient. It seems unlikely that synonymous variation in the remaining one third of the sites would do much to escape the depletion of tRNAs used by virus genes. Nonetheless, expression of MHC genes is crucial to vertebrate immunity, and it is obvious that any mechanism that obstructs MHC expression would be adaptive to pathogens . It is possible that intracellular parasites such as viruses evolved to evade recognition and inhibit host immunity by targeting mRNA of MHC genes. This could be mediated via molecular interactions, where virus RNA or proteins hybridize or bind to MHC mRNA to inhibit translation. Regulation at the mRNA level could prevent MHC molecules from being expressed, thereby reducing presentation of (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted February 23, 2026. ; https://doi.org/10.64898/2026.02.23.707394doi: bioRxiv preprint 15 antigen to T-cells. Such regulation by pathogens targeting MHC expression at the mRNA level would necessarily elicit an evolutionary response in the hosts to evade molecular interactions with virus molecules. It is well described that interactions between endogenous molecules and mRNA leave an evolutionary footprint in synonymous nucleotide variation (3, 5), and I propose the hypothesis that synonymous nucleotide variation in the MHC is engaged in a coevolutionary arms race to escape molecular interactions between MHC mRNA and pathog en molecules. Due to secondary and tertiary structures of the mRNA, i t is likely that not all sites would be equally exposed to molecular interactions with pathogenic RNA or proteins, and this may explain why the synonymous variation appears in spikes along the exon 3 sequences (Fig. 1). Only a single study has previously investigated synonymous variation in the MHC. Matsushita & Kano-Sueoka (2023) investigated sequence variation in the HLA -A (a human MHC -I locus) . While that study did not test for natural selection or address phenotypic or evolutionary effects, they did observe conspicuous non -random spikes of synonymous nucleotide variation along the coding sequences, that resemble the patterns in the Great Reed Warbler MHC-I exon 3 (46). Synonymous variation was more abundant in the Great Reed Warbler MHC-I exon 3 sequences compared to the HLA -A, but that may be due to the fact that (46) only studied variation in the HLA-A locus (disregarding the other class I loci HLA-B and -C), whereas the sequences in my study represent multiple MHC-I loci (29). The correlations between synonymous nucleotide variation in the MHC -I and individual survival and fitness in the study population differed between males and females (Fig. 6a, b). This difference in the effects of synonymous variation in the MHC-I is likely associated with regulatory effects of sex hormones on immune responses, that cause males to generally have weaker immune responses and increased risk of pathogen infection compared to females (e.g. (47–51)). I have previously proposed that differences between the immune response phenotypes of males and females can drive sexually antagonistic (SA) selection on genes that exert quantitative effects on immune responses, including the MHC , and found evidence for an unresolved sexual conflict over the number of different MHC-I alleles in Great Reed Warblers (24, 25). Additional evidence for sexual conflict on MHC genes recently emerged from a study on a social mammal, the Banded Mongoose Mungos mungo (28). If the observed sex difference in the effects of synonymous variation in the Great Reed Warbler MHC-I is to be explained by the SA selection hypothesis, then synonymous variation in the MHC-I should fulfill the assumption of exerting quantitative effects on immune responses, cf. (24). This is consistent with the hypothesis that synonymous variation serves an adaptive function in evading regulation of MHC expression by pathogens . It has long been recognized that polymorphism in MHC genes contributes to quantitative differences in overall antibody production (52). Hence, if pathogens regulate MHC expression by targeting MHC mRNA, host mechanisms that evade such regulation would indeed be predicted to exert quantitative effects on immune responses. Synonymous variation has been shown to exert multiple effects on the pathway from gene to protein, and the evidence that MHC genes are subject to natural selection on synonymous variation (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted February 23, 2026. ; https://doi.org/10.64898/2026.02.23.707394doi: bioRxiv preprint 16 invites reflection over some current practices in MHC research. If pathogens regulate MHC expression by targeting mRNA, then expression of MHC alleles is probably biased in infected cells. The practice of estimating MHC expression by abundance of RNA-derived cDNA needs to be evaluated, as mRNA abundance measured at the organism level may not reflect MHC expression in infected cells . Furthermore, it is common to analyze natural selection on protein sequences by comparing rates of non-synonymous to synonymous variation (dN/dS, aka. Ω). This analysis assumes that rates of synonymous variation (dS) can be used as a representation of the

rate of evolution. However, if synonymous variation deviates from the background rate of evolution, such models may fail to provide a meaningful estimate of natural selection on protein sequences (3). Yang & Nielsen disputed concerns that natural selection on synonymous nucleotide variation should compromise dN/dS analyses (4), but my present results indicate that disregarding selection on synonymous variation leads to substantial overestimation of positive selection in MHC sequences by codeml (Table 2).

Here I show that MHC genes are subject to natural selection on synonymous nucleotide variation and that the extent of synonymous variation has direct measurable effects on the Darwinian fitness of individuals. These observations are remarkable and call for a paradigm shift in the way that we approach studies of the MHC. Unraveling the driving forces behind selection on synonymous variation in the MHC is of great importance to our understanding of evolution in this extremely important locus. I have taken a first few steps on the road and shown that most sites in the Great Reed Warbler MHC-I exon 3 are under strong purifying selection on synonymous codon usage, which is associated with the abundance of tRNA isotypes in the genome and likely driven by selection for increased translational efficiency. However, that is only half of the story. Significant spikes of synonymous variation appear in about one third of the sites in the Great Reed Warbler MHC-I exon 3, consistent with recent observations in the HLA-A, and I propose that these are best explained in the context of a coevolutionary arms race, where synonymous nucleotide variation serve s an adaptive function to escape inhibitory molecular interactions, e.g. between MHC mRNA and pathogen molecules. The prevailing lack of knowledge on synonymous variation in the MHC is potentially problematic and further uncovering the mechanisms by which synonymous nucleotide variation in the MHC exerts its biological effects is central to advance our understanding of the MHC and host-pathogen coevolution. This study introduces the SynDist function in MHCtools to quantify synonymous nucleotide differences (53) and demonstrates a novel approach to test site specific selection on synonymous variation by comparing empirical sequences with synonymous data sets simulated under the assumption of no selection on codon usage. I hope that these methodological advancements will pave the way for future studies on synonymous variation in the MHC. (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted February 23, 2026. ; https://doi.org/10.64898/2026.02.23.707394doi: bioRxiv preprint 17

Data set This study employed data from previous studies of a wild population of Great Reed Warblers at lake Kvismare in Sweden (59°10’N, 15°25’E) . The data set is based on exhaustive field observations and samples collected in the period 1984–2004 and details on the field methods can be found in (25, 29, 54–59). For genetic analyses, I employed the Great Reed Warbler genome assembly (GenBank: ASM2153481v1) (33, 34) in conjunction with a data set of 390 MHC -I exon 3 DNA sequences (25, 29). The MHC-I data set originated from 559 individuals (141 adult males, 131 adult females, and 287 chicks) and the whole genome was sequenced from one of the individuals in the MHC -I data set. The size of the genome assembly was 1.2 Gb in 3,012 scaffolds with scaffold N50 = 21.4 Mb and coverage 45x. The 390 sequences of the MHC-I exon 3 were aligned and trimmed to open reading frame according to conserved residues and sequence motifs (60, 61). The sequences were 87 codons long and contained no gaps. Details on DNA extraction, sequencing, and genotyping of the MHC-I exon 3 can be found in (25, 29). The Great Reed Warbler MHC-I exon 3 sequences are available at GenBank (accession numbers: MH468831–MH469159; MT193762–MT193822). Fitness analyses were conducted on a subset of samples from the genetic data set for which fitness data has previously been published in (25). The fitness data set included lifetime observations of 88 adult males and 100 adult females. These individuals harbored 329 of the 390 MHC-I exon 3 sequences between them, with 6 –24 different sequences per individual. Further details can be found in (25). Selection analyses I tested 12 substitution models on the Great Reed Warbler MHC-I exon 3 sequence s in PhyML version 3.1 (62, 63) using maximum likelihood estimation of nucleotide frequencies and tree topology optimization. Among the 12 models, t he generalized time-reversible (GTR) model had the lowest Akaike Information Criterion (AIC) value (Table S18). I used codeml from the PAML software package (35, 36) to test for natural selection on synonymous codon variation among the sequences, using a GTR tree as input. I set codeml to assume one Ω (i.e., dN/dS) ratio for all branches and specified site models M2 a (estimating three categories of Ω with one category of sites evolving under positive selection (Ω > 1)) and M8 (estimating a beta distribution for Ω 1)) (37). Within each site model, I employed nested FMutSel vs. FMutSel0 models to test the null hypothesis that variation in codon usage is due to mutation bias alone and not affected by selection acting at synonymous sites (4). The log likelihoods of the nested FMutSel and FMutSel0 models were compared using a likelihood ratio test with the formula: 2 × Δln(L) ~ χ2, with 41 degrees of freedom of the χ2 distribution (reflecting the difference in number of parameters between the models ). Sites with Ω > 1 were inferred by Bayes Empirical Bayes analysis for each model (38). (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted February 23, 2026. ; https://doi.org/10.64898/2026.02.23.707394doi: bioRxiv preprint 18 Consensus sequence and prediction of protein folding I analyzed the frequencies of amino acids in each site among the 390 MHC -I exon 3 sequences and derived a consensus sequence of the amino acids most frequently observed at each site. I employed AlphaFold 3 (39) at alphafoldserver.com to predict the folded protein structure from the consensus amino acid sequence (Fig. S1). Synonymous genetic variation All downstream analyses were conducted in R v. 4.5.1 (64), except where otherwise specified. Fasta files were handled in R using the package seqinr (65). I used the SynDist function in MHCtools v. 1.6 (53) to quantify synonymous nucleotide variation among the Great Reed Warbler MHC-I exon 3 sequences. I specified analysis="codon" to obtain counts of synonymous changes per base and per (codon) site among all pairwise sequence comparisons in the data set. I also used the setting analysis="dist" to obtain the number of synonymous nucleotide changes in each pairwise sequence comparison in the data set as well as the mean number of synonymous nucleotide changes in pairwise comparisons among the sequences in each individual sample. I calculated relative frequencies of codon usage as the proportion of counts for a codon out of the total counts for all codons within each synonymous block, i.e., codons encoding the same amino acid. For comparisons between amino acids, I normalized the relative codon frequencies to 1 by multiplying each frequency with the number of synonymous codons in each block, thereby generating the measure “relative synonymous codon usage” (RSCU), cf. (10, 66). I quantified codon usage bias (CUB) for each amino acid as the variance of RSCU within each synonymous block, following the rationale that a bias towards certain codons will skew the values of RSCU within each synonymous block and increase the variance. Data simulations If one assumes that an observed sequence alignment is a snapshot of a Markov process of synonymous codon substitutions that has run for infinitely long time, the probability of observing a certain codon can be regarded as a function of the mutation bias among nucleotides and selection on codon usage, following (4). Hence, under the null hypothesis of no selection on codon usage, the probability of observing a certain codon depends only on mutation bias, and thus regarding the process of synonymous codon substitution as a Markov process enables simulation of nucleotide sequences that reflect the null hypothesis using mutation bias parameters alone. I employed custom scripts to perform simulations of MHC -I exon 3 sequences reflecting the null hypothesis that variation in codon usage is not affected by selection acting at synonymous sites . I first used the base frequencies by nucleotide position (3x4 table) derived from the models in codeml (Table S19) to compute a table of relative probabilities of codon usage under neutrality , i.e., reflecting the scenario where the probabilities of observing synonymous codon variants depend only on the mutation bias, cf. (4) (Table S20). I then ran 1,000 simulations of the data set of 390 MHC-I exon 3 sequences, where for each site in each sequence, a codon was randomly sampled among the (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted February 23, 2026. ; https://doi.org/10.64898/2026.02.23.707394doi: bioRxiv preprint 19 synonymous codons encoding the observed amino acid. In the random sampling step, the relative probabilities of synonymous codons encoding the observed amino acid in the computed table of probabilities of codon usage under neutrality (Table S20) were used as weights. The simulated codons were concatenated into nucleotide sequences and collated in 1,000 fasta files, each containing 390 simulated DNA sequences synonymous to the empirical Great Reed Warbler MHC- I exon 3 sequences. I then used the SynDist function in MHCtools v. 1.6 to quantify the synonymous genetic variation among the 390 sequences in each of the 1,000 simulated data sets. I specified analysis="codon" to obtain counts of synonymous changes per base and per (codon) site among all pairwise sequence comparisons in the data set. From the results, I extracted the proportion of pairwise sequence comparisons that harbor synonymous substitutions for each site. In addition, I calculated relative frequencies of codon usage, RSCU, and codon usage bias as described above. Finally, I calculated means and 95% confidence intervals among the 1,000 simulated data sets and calculated p-values as the proportion of observations in the simulated data sets more extreme than the value observed in the empirical data set. Abundance of tRNA isotypes I used tRNAscan -SE v. 2.0.12 (40) to scan the Great Reed Warbler genome assembly for tRNA genes. I ran tRNAscan-SE with the options ‘--mt vert’ to identify potential mitochondrial origin of detected tRNAs and ‘ --detail’ to obtain isotype -specific model classification results . tRNAscan- SE predicted 815 tRNAs in the first pass, 571 of which were confirmed by Infernal (second pass). Among the 571 detected tRNAs, 325 decoded standard amino acids, 1 was a selenocysteine tRNA, 7 had undetermined isotypes, 30 had mismatch isotypes, and 208 were pseudogenes. I used only the 325 confirmed cytosolic standard amino acid tRNAs in further analyses. I calculated relative frequencies of tRNA isotypes as the proportion of counts for an isotype out of the total counts f or isotypes for the same amino acid . For comparisons between amino acids, I normalized the relative tRNA isotype frequencies by multiplying each frequency with the number of synonymous codons in each block. To test whether synonymous codon usage was affected by tRNA isotype availability, I calculated Pearson correlation coefficients between RSCU and normalized relative frequencies of tRNA isotypes within each synonymous block. The analysis included only amino acids that were matched by more than one tRNA isotype in the genome and excluded codons that were not matched by a cognate tRNA isotype in the genome . The Pearson correlations were calculated using RSCU measured among all 87 sites of the MHC -I sequences as well as RSCU measured among subsets of sites where synonymous variation was under strong purifying selection (n = 56) or where purifying selection was relaxed (n = 31). Pearson correlations were also calculated using RSCU for each of the 1,000 data sets simulated under the assumption of no selection on synonymous codon usage, allowing derivation of 95% confidence intervals and p -values for the correlat ions between normalized relative frequencies of tRNA isotypes and RSCU in the empirical data set. (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted February 23, 2026. ; https://doi.org/10.64898/2026.02.23.707394doi: bioRxiv preprint 20 Fitness effects As estimates of individual Darwinian fitness, I used lifetime number of fledged offspring, defined as the number of fledged offspring over an individual’s lifetime, and offspring fledging success, defined as the lifetime number of fledged offspring corrected for life span. Previous studies have shown that Great Reed Warbler s that failed their first breeding attempt were unlikely to be observed in the study area in following years (25, 55, 67). Estimates of life span and lifetime reproductive success are therefore not reliable for unsuccessful first-time breeders and I excluded 17 individuals that failed to rear any offspring from analyses of fitness effects, resulting in a data set of 77 males and 96 females. I used generalized linear models to test the effects of synonymous genetic variation in the Great Reed Warbler MHC-I exon 3 on life span , lifetime number of fledged offspring , and offspring fledging success following the suggested model design in (68). As both life span and lifetime number of fledged offspring follow negative binomial distributions, the generalized linear models were run with negative binomial errors using the package MASS in R (69). The aggregation parameters of the negative binomial distributions were estimated by maximum likelihood and specified in the model formulae (70). Effects on offspring fledging success were modelled by including life span as covariate in a generalized linear model with lifetime number of fledg ed offspring as dependent variable. The purpose of the models was to test the effects of the mean number of synonymous nucleotide changes in pairwise comparisons among the sequences in each individual sample. The full models included the total number of different MHC -I sequences per individual as covariate, sex as fixed factor, and the interactions ‘total number of different MHC-I sequences ´ sex’ and ‘mean number of synonymous nucleotide changes ´ sex’. The step function in R was employed for model simplification and simulateResiduals and testDispersion from the DHARMa package (71) were used for model diagnostics. The final models on life span and lifetime number of fledged offspring showed significant or nearly significant interactions between the mean number of synonymous nucleotide changes and sex, and the models were therefore also run independently for each sex. To ascertain that the observed effects of synonymous variation were not influenced by natural selection on amino acid sequences, I calculated the mean proportion of amino acid changes between pairs of sequences in each individual (i.e., amino acid p -distance) using the DistCalc function in MHCtools v. 1.6 (53). I repeated the generalized linear models as above including amino acid p-distance as covariate and the interaction ‘ amino acid p-distance ´ sex’ in the full models. In addition, I repeated the models as described above, but replacing the mean number of synonymous nucleotide changes with amino acid p-distance as independent variable. (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted February 23, 2026. ; https://doi.org/10.64898/2026.02.23.707394doi: bioRxiv preprint 21

1. D. Agashe, “Evolutionary forces that generate SNPs: the evolutionary impacts of synonymous mutations” in Single Nucleotide Polymorphisms , Z. E. Sauna, C. Kimchi - Sarfaty, Eds. (Springer, Cham, Switzerland, 2022), pp. 15–36. 2. M. Kimura, Evolutionary rate at the molecular level. Nature 217, 624–626 (1968). 3. J. Zhang, W. Qian, Functional synonymous mutations and their evolutionary consequences. Nat. Rev. Genet. 26, 789–804 (2025). 4. Z. Yang, R. Nielsen, Mutation -selection models of codon substitution and their use to estimate selective strengths on codon usage. Mol. Biol. Evol. 25, 568–579 (2008). 5. S. A. Shabalina, N. A. Spiridonov, A. Kashina, Sounds of silence: Synonymous nucleotides as a key to biological regulation and complexity. Nucleic Acids Res. 41, 2073–2094 (2013). 6. Z. E. Sauna, C. Kimchi-Sarfaty, Understanding the contribution of synonymous mutations to human disease. Nat. Rev. Genet. 12, 683–691 (2011). 7. Y . Liu, Q. Yang, F. Zhao, Synonymous but not silent: the codon usage code for gene expression and protein folding. Annu. Rev. Biochem. 90, 375–401 (2021). 8. W. Qian, J. R. Yang, N. M. Pearson, C. Maclean, J. Zhang, Balanced codon usage optimizes eukaryotic translational efficiency. PLoS Genet. 8, 1–18 (2012). 9. M. Sun, J. Zhang, Preferred synonymous codons are translated more accurately: Proteomic evidence, among-species variation, and mechanistic basis. Sci. Adv. 8, 1–10 (2022). 10. P. S. Spencer, J. M. Barral, Genetic code redundancy and its influence on the encoded polypeptides. Comput. Struct. Biotechnol. J. 1, 1–8 (2012). 11. M. J. Moss, L. M. Chamness, P. L. Clark, The effects of codon usage on protein structure and folding. Annu. Rev. Biophys. 53, 87–108 (2024). 12. A. B. Stergachis, E. Haugen, A. Shafer, F. Wenqing, B. Vernot, A. Reynolds, A. Raubitschek, S. Ziegler, E. LeProust, J. M. Akey, J. A. Stamatoyannopoulos, Exonic transcription factor binding directs codon choice and affects protein evolution. Science 342, 1367–1372 (2013). 13. J. V . Chamary, J. L. Parmley, L. D. Hurst, Hearing silence: Non -neutral evolution at synonymous sites in mammals. Nat. Rev. Genet. 7, 98–108 (2006). 14. T. Ohta, J. H. Gillespie, Development of neutral and nearly neutral theories. Theor. Popul. Biol. 49, 128–142 (1996). 15. V . Bali, Z. Bebok, Decoding mechanisms by which silent codon changes influence protein biogenesis and function. International Journal of Biochemistry and Cell Biology 64, 58–74 (2015). 16. Y . F. Huang, A. Siepel, Estimation of allele-specific fitness effects across human protein - coding sequences and implications for disease. Genome Res. 29, 1310–1321 (2019). 17. J. Kaufman, Unfinished business: Evolution of the MHC and the adaptive immune system of jawed vertebrates. Annu. Rev. Immunol 36, 383–409 (2018). 18. J. Klein, A. Sato, The HLA system - Second of two parts. N. Engl. J. Med. 343, 782–786 (2000). (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted February 23, 2026. ; https://doi.org/10.64898/2026.02.23.707394doi: bioRxiv preprint 22 19. J. Klein, A. Sato, N. Nikolaidis, MHC, TSP, and the origin of species: From immunogenetics to evolutionary genetics. Annu. Rev. Genet. 41, 281–304 (2007). 20. M. Nei, X. Gu, T. Sitnikova, Evolution by the birth-and-death process in multigene families of the vertebrate immune system. Proceedings of the National Academy of Sciences 94, 7799–7806 (1997). 21. S. B. Piertney, M. K. Oliver, The evolutionary ecology of the major histocompatibility complex. Heredity (Edinb). 96, 7–21 (2006). 22. J. Lighten, A. S. T. Papadopulos, R. S. Mohammed, B. J. Ward, I. G. Paterson, L. Baillie, I. R. Bradbury, A. P. Hendry, P. Bentzen, C. Van Oosterhout, Evolutionary genetics of immunological supertypes reveals two faces of the Red Queen. Nat. Commun. 8, 1 –10 (2017). 23. M. J. Ejsmond, J. Radwan, Red queen processes drive positive selection on major histocompatibility complex (MHC) genes. PLoS Comput. Biol. 11, 1–14 (2015). 24. J. Roved, H. Westerdahl, D. Hasselquist, Sex differences in immune responses: Hormonal effects, antagonistic selection, and evolutionary consequences. Horm. Behav. 88, 95–105 (2017). 25. J. Roved, B. Hansson, M. Tarka, D. Hasselquist, H. Westerdahl, Evidence for sexual conflict over MHC diversity in a wild songbird. Proceedings of the Royal Society B 285, 1–9 (2018). 26. T. Kamiya, K. O’Dwyer, H. Westerdahl, A. Senior, S. Nakagawa, A quantitative review of MHC-based mating preference: The role of diversity and dissimilarity. Mol. Ecol. 23, 5151– 5163 (2014). 27. M. Milinski, The major histocompatibility complex, sexual selection, and mate choice. Annu. Rev. Ecol. Evol. Syst. 37, 159–186 (2006). 28. N. Schubert, H. J. Nichols, F. Mwanguhya, R. Businge, S. Kyambulima, K. Mwesige, J. I. Hoffman, M. A. Cant, J. C. Winternitz, Sex-dependent influence of major histocompatibility complex diversity on fitness in a social mammal. Mol. Ecol. 34, 1–15 (2025). 29. J. Roved, B. Hansson, M. Stervander, D. Hasselquist, H. Westerdahl, MHCtools – an R package for MHC high‐throughput sequencing data: genotyping, haplotype and supertype inference, and downstream genetic analyses in non‐model organisms. Mol. Ecol. Resour. 22, 2775–2792 (2022). 30. F. Pierini, T. L. Lenz, Divergent allele advantage at human MHC genes: Signatures of past and ongoing selection. Mol. Biol. Evol. 35, 2145–2158 (2018). 31. F. Chen, P. Wu, S. Deng, H. Zhang, Y . Hou, Z. Hu, J. Zhang, X. Chen, J. R. Yang, Dissimilation of synonymous codon usage bias in virus –host coevolution due to translational selection. Nat. Ecol. Evol. 4, 589–600 (2020). 32. S. Eltschkner, S. Mellinger, S. Buus, M. Nielsen, K. M. Paulsson, K. Lindkvist -Petersson, H. Westerdahl, The structure of songbird MHC class I reveals antigen binding that is flexible at the N-terminus and static at the C-terminus. Front. Immunol. 14, 1–14 (2023). 33. H. Sigeman, M. Strandh, E. Proux -Wéra, V . E. Kutschera, S. Ponnikas, H. Zhang, M. Lundberg, L. Soler, I. Bunikis , M. Tarka, D. Hasselquist, B. Nystedt, H. Westerdahl, B. (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted February 23, 2026. ; https://doi.org/10.64898/2026.02.23.707394doi: bioRxiv preprint 23 Hansson, Avian Neo-Sex Chromosomes Reveal Dynamics of Recombination Suppression and W Degeneration. Mol. Biol. Evol. 38, 5275–5291 (2021). 34. S. Ponnikas, H. Sigeman , M. Lundberg, B. Hansson, Extreme variation in recombination rate and genetic diversity along the Sylvioidea neo‐sex chromosome. Mol. Ecol. 31, 3566– 3583 (2022). 35. Z. Yang, PAML: a program package for phylogenetic analysis by maximum likelihood. Bioinformatics 13, 555–556 (1997). 36. Z. Yang, PAML 4: Phylogenetic analysis by maximum likelihood. Mol. Biol. Evol. 24, 1586–1591 (2007). 37. Z. Yang, R. Nielsen, N. Goldman, A. -M. K. Pedersen, Codon -substitution models for heterogeneous selection pressure at amino acid sites. Genetics 155, 431–449 (2000). 38. Z. Yang, W. S. W. Wong, R. Nielsen, Bayes empirical Bayes inference of amino acid sites under positive selection. Mol. Biol. Evol. 22, 1107–1118 (2005). 39. J. Abramson, J. Adler, J. Dunger, R. Evans, T. Green, A. Pritzel, O. Ronneberger, L. Willmore, A. J. Ballard, J. Bambrick, S. W. Bodenstein, D. A. Evans, C. C. Hung, M. O’Neill, D. Reiman, K. Tunyasuvunakool, Z. Wu, A. Žemgulytė, E. Arvaniti, C. Beattie, O. Bertolli, A. Bridgland, A. Cherepanov, M. Congreve, A. I. Cowen -Rivers, A. Cowie, M. Figurnov, F. B. Fuchs, H. Gladman, R. Jain, Y . A. Khan, C. M. R. Low, K. Perlin, A. Potapenko, P. Savy, S. Singh, A. Stecula, A. Thillaisundaram, C. Tong, S. Yaknee n, E. D. Zhong, M. Zielinski, A. Žídek, V . Bapst, P. Kohli, M. Jaderberg, D. Hassabis, J. M. Jumper, Accurate structure prediction of biomolecular interactions with AlphaFold 3. Nature 630, 493–500 (2024). 40. P. P. Chan, B. Y . Lin, A. J. Mak, T. M. Lowe, tRNAscan-SE 2.0: improved detection and functional classification of transfer RNA genes. Nucleic Acids Res. 49, 9077–9096 (2021). 41. X. Wu, M. Xu, J. R. Yang, J. Lu, Genome-wide impact of codon usage bias on translation optimization in Drosophila melanogaster. Nature Communications 15, 1–16 (2024). 42. V . Peona, M. H. Weissensteiner, A. Suh, How complete are “complete” genome assemblies?—an avian perspective. Mol. Ecol. Resour. 18, 1188–1195 (2018). 43. K. Näpflin, E. A. O’Connor, L. Becks, S. Bensch, V . A. Ellis, N. Hafer -Hahmann, K. C. Harding, M. T. Olsen, J. Roved, T. B. Sackton, A. J. Shultz, V . Venkatakrishnan, E. Videvall, H. Westerdahl, J. C. Winternitz, S. V . Edwards, Genomics of hosts-pathogen interactions: challenges and opportunities across ecological and spatiotemporal scales. PeerJ 7, 1–37 (2019). 44. R. Venskutonytė, S. Kjellström, E. A. O’Connor, H. Westerdahl, K. Lindkvist‐Petersson, MHC I of the great reed warbler promotes a flat peptide binding mode. Immunology 176, 508–519 (2025). 45. M. E. J. Woolhouse, J. P. Webster, E. Domingo, B. Charlesworth, B. R. Levin, Biological and biomedical implications of the co -evolution of pathogens and their hosts. Nat. Genet. 32, 569–577 (2002). (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted February 23, 2026. ; https://doi.org/10.64898/2026.02.23.707394doi: bioRxiv preprint 24 46. T. Matsushita, T. Kano -Sueoka, Non -random codon usage of synonymous and non - synonymous mutations in the human HLA-A gene. J. Mol. Evol. 91, 169–191 (2023). 47. S. L. Klein, Hormonal and immunological mechanisms mediating sex differences in parasite infection. Parasite Immunol. 26, 247–264 (2004). 48. S. L. Klein, K. L. Flanagan, Sex differences in immune responses. Nat. Rev. Immunol. 16, 626–638 (2016). 49. J. Alexander, W. H. Stimson, Sex -hormones and the course of parasitic infection. Parasitology Today 4, 189–193 (1988). 50. M. Zuk, K. A. McKean, Sex differences in parasite infections: Patterns and processes. Int. J. Parasitol. 26, 1009–1023 (1996). 51. Y . Z. Foo, S. Nakagawa, G. Rhodes, L. W. Simmons, The effects of sex hormones on immune function: a meta-analysis. Biological Reviews 92, 551–571 (2017). 52. A. Puel, D. Mouton, Genes responsible for quantitative regulation of antibody production. Crit Rev Immunol 16, 223–250 (1996). 53. J. Roved, MHCtools : Analysis of MHC data in non -model species. CRAN 1.6 (2026). https://cran.r-project.org/package=MHCtools. 54. S. Bensch, D. Hasselquist, Nest predation lowers the polygyny threshold - a new compesation model. American Naturalist 138, 1297–1306 (1991). 55. S. Bensch, D. Hasselquist, Territory infidelity in the polygynous great reed warbler Acrocephalus arundinaceus: The effect of variation in territory attractiveness. Journal of Animal Ecology 60, 857–871 (1991). 56. D. Hasselquist, S. Bensch, Trade -off between mate guarding and mate attraction in the polygynous great reed warbler. Behav. Ecol. Sociobiol. 28, 187–193 (1991). 57. S. Bensch, D. Hasselquist, B. Nielsen, B. Hansson, Higher fitness for philopatric than for immigrant males in a semi-isolated population of great reed warblers. Evolution (N Y). 52, 877–883 (1998). 58. D. Hasselquist, Polygyny in great reed warblers: A long-term study of factors contributing to male fitness. Ecology 79, 2376–2390 (1998). 59. M. Tarka, M. Akesson, D. Hasselquist, B. Hansson, Intralocus sexual conflict over wing length in a wild migratory bird. American Naturalist 183, 62–73 (2014). 60. P. J. Bjorkman, M. a Saper, B. Samraoui, W. S. Bennett, J. L. Strominger, D. C. Wiley, The foreign antigen binding site and T cell recognition regions of class I histocompatibility antigens. Nature 329, 512–8 (1987). 61. A. L. Hughes, M. Nei, Evolution of the major histocompatibility complex: independent origin of nonclassical class I genes in different groups of mammals. Mol. Biol. Evol. 6, 559– 579 (1989). 62. S. Guindon, O. Gascuel, A simple, fast, and accurate algorithm to estimate large phylogenies by maximum likelihood. Syst. Biol. 52, 696–704 (2003). (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted February 23, 2026. ; https://doi.org/10.64898/2026.02.23.707394doi: bioRxiv preprint 25 63. S. Guindon, J. F. Dufayard, V . Lefort, M. Anisimova, W. Hordijk, O. Gascuel, New algorithms and methods to estimate maximum -likelihood phylogenies: Assessing the performance of PhyML 3.0. Syst. Biol. 59, 307–321 (2010). 64. R Core Team, R: A language and environment for statistical computing. R Foundation for Statistical Computing 4.5.1 (2025). https://www.R-project.org/. 65. D. Charif, J. R. Lobry, “SeqinR 1.0-3: A contributed package to the R project for statistical computing devoted to biological sequences retrieval and analysis” in Structural Approaches to Sequence Evolution, U. Bastolla, M. Porto, H. E. Roman, M. Vendruscolo, Eds. (Springer Verlag, New York, 2007), pp. 207–232. 66. P. M. Sharp, T. M. F. Tuohy, K. R. Mosurski, Codon usage in yeast: cluster analysis clearly differentiates highly and lowly expressed genes. Nucleic Acids Res. 14, 5125–5143 (1986). 67. B. Hansson, S. Bensch, D. Hasselquist, B. Nielsen, Restricted dispersal in a long -distance migrant bird with patchy distribution, the great reed warbler. Oecologia 130, 536 –542 (2002). 68. W. P. Gilks, J. K. Abbott, E. H. Morrow, Sex differences in disease genetics: Evidence, evolution, and detection. Trends in Genetics 30, 453–463 (2014). 69. W. N. Venables, B. D. Ripley, Modern Applied Statistics with S (Springer, New York, NY , Fourth., 2002). 70. M. J. Crawley, The R Book (John Wiley & Sons, Ltd., Chichester, Second., 2013). 71. F. Hartig, DHARMa: Residual diagnostics for hierarchical (multi-level / mixed) regression models. (2025).

I wish to thank my partner Tianhao Zhao for patience and support. In memory of Bengt Olle Bengtsson (1946-2025). Funding: I received no funding in support for this research. Author contributions: This study was conceived, designed and carried out by J.R. MHCtools is developed and maintained by J.R. The data simulation approach to test site specific selection on synonymous variation was conceived and developed by J.R. Competing interests: I declare no competing interests. Data and materials availability: MHCtools v. 1.6 including user manual and documentation is available at CRAN: https://cran.r-project.org/package=MHCtools. The data sets are available at the Dryad repository: https://datadryad.org/dataset/doi:10.5061/dryad.b321hf1 and the Zenodo repository: https://doi.org/10.5281/zenodo.3716048. The great reed warbler genome assembly and MHC-I exon 3 sequences are available at GenBank: https://ncbi.nlm.nih.gov (accession numbers: MH468831–MH469159; MT193762–MT193822; genome assembly: ASM2153481v1). (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted February 23, 2026. ; https://doi.org/10.64898/2026.02.23.707394doi: bioRxiv preprint