The Rare IL22RA2 Signal Peptide Coding Variant rs28385692 Decreases Secretion of IL-22BP Isoform-1, -2 and -3 and Is Associated with Risk for Multiple Sclerosis

The IL22RA2 locus is associated with risk for multiple sclerosis (MS) but causative variants are yet to be determined. In a single nucleotide polymorphism (SNP) screen of this locus in a Basque population, rs28385692, a rare coding variant substituting Leu for Pro at position 16 emerged significantly (p = 0.02). This variant is located in the signal peptide (SP) shared by the three secreted protein isoforms produced by IL22RA2 (IL-22 binding protein-1(IL-22BPi1), IL-22BPi2 and IL-22BPi3). Genotyping was extended to a Europe-wide case-control dataset and yielded high significance in the full dataset (p = 3.17 × 10−4). Importantly, logistic regression analyses conditioning on the main known MS-associated SNP at this locus, rs17066096, revealed that this association was independent from the primary association signal in the full case-control dataset. In silico analysis predicted both disruption of the alpha helix of the H-region of the SP and decreased hydrophobicity of this region, ultimately affecting the SP cleavage site. We tested the effect of the p.Leu16Pro variant on the secretion of IL-22BPi1, IL-22BPi2 and IL-22BPi3 and observed that the Pro16 risk allele significantly lowers secretion levels of each of the isoforms to around 50%–60% in comparison to the Leu16 reference allele. Thus, our study suggests that genetically coded decreased levels of IL-22BP isoforms are associated with augmented risk for MS.


Introduction
Multiple sclerosis (MS) is a chronic inflammatory CNS disorder triggered by environmental factors in individuals with a susceptible genetic background. Current research implies more than 230 autosomal risk variants, many of which are located within or close to genes exerting functions in the peripheral immune system or in CNS-resident microglia [1]. Single nucleotide polymorphisms (SNPs) located at the interleukin-22 receptor subunit alpha-2 (IL22RA2) gene locus were originally reported to be associated with risk for MS in Scandinavian [2] and Basque [3] populations. A series of genome-wide association studies (GWAS)-based approaches have subsequently firmly established the association of IL22RA2 with risk for MS [1,[4][5][6].
The main known function of IL22RA2 is to produce interleukin-22 binding protein (IL-22BP), a secreted inhibitor of IL-22. IL-22, a member of the IL-10 family, is produced by a wide range of immune cells and can exert both pro-and anti-inflammatory effects [7,8]. Various lines of evidence suggest that the IL-22/IL-22BP axis has an important function in MS and neuroinflammation. Il22ra2-deficient mice experience a more benign course of disease in the experimental autoimmune encephalomyelitis (EAE) model [9]. The IL-22 receptor is highly expressed in blood-brain barrier (BBB) endothelial cells from Cells 2020, 9,175 3 of 23 patients with MS but not in healthy controls, and IL-22 contributes to the permeabilization of the BBB and recruitment of CD4 + T lymphocytes to the CNS [10]. A decrease of IL-22 levels correlates with the recovery phase of EAE in rats [11], and serum levels of IL-22 were found to be elevated in MS patients compared to healthy controls [12]. Perriard et al. [13] demonstrated that IL-22 targets astrocytes in the human brain and confers increased survival to these cells. They also found higher expression of IL-22BP mRNA in monocytes and monocyte-derived dendritic cells of MS patients compared to healthy controls.
IL22RA2 is capable of expressing three partially distinct isoforms that share an identical signal peptide (SP) at their N-terminus and lack intracellular and transmembrane domains but differ in their binding capacity of IL-22. Isoform 2 (UniProt nomenclature) shows the highest capacity of binding and inhibiting IL-22 [14,15], with a 20-to 1000-fold higher affinity than a soluble variant of the signal-transducing cell surface receptor [16][17][18]. Isoform 3 has also been demonstrated to bind IL-22, although with a similar affinity to that of the cell surface receptor [16,19]. Recently, we showed that the longest isoform, i.e., isoform 1, is not capable of binding IL-22 and displays hallmarks of a poorly secreted, intracellularly retained protein with intrinsic capacity to trigger the unfolded protein response (UPR; [20]).
Although the association of IL22RA2 with MS is now established and accumulating evidence points to an influence of IL-22 and IL-22RA2 in EAE and MS, the mechanism underlying the genetic association remains elusive. Here, we performed a SNP screen of the IL22RA2 locus in a Basque population in order to localize the most important association signal(s) within this locus and confirmed association of an infrequent coding SNP in a European cohort. We used dedicated in silico and wet experimentation methods to discover potentially causal variants that may explain the association of this gene with MS.

Patients and Controls
All patients were diagnosed with definite MS [21,22]. Written informed consent was obtained from all subjects, and the study was approved by the local ethics committees. Table 1 shows the clinical and demographic data of the patients and controls enrolled in this study. The fine-mapping was completed in the Bilbao dataset, comprising patients registered at the Basurto hospital (Bilbao, Basque Country, Spain) and controls provided by the Basque BioBank for Research-OEHUN (www.biobancovasco.org). Additionally, genotyping data of three SNPs (rs276466, rs10484798 and rs6570136) in the Bilbao cohort were available from the aforementioned screening [3], and these were included in the haplotype and logistic regression analyses.
The non-synonymous SNP, rs28385692, was further genotyped in additional patients and age-, gender-and ethnicity-matched controls from Donostia (Instituto Biodonostia, Basque Country, Spain), Barcelona (Hospital Vall d'Hebron), Madrid (Hospital Clínico S. Carlos), Andalucía (Instituto de Parasitología y Biomedicina "Lopez-Neyra"), Germany (various centers) and France (various centers) ( Table 1). Considering the SNP with the lowest MAF (rs28385692, 1000G European frequency = 0.032) and combining all datasets, this study had over 80% statistical power to detect allelic odds ratios ≥ 1.18 with a type-1-error rate of alpha = 0.05 (calculated with the Genetic Power Calculator allelic test [23]). When considering only the Bilbao dataset, we had over 80% power to detect odds ratios ≥ 1.8, assuming a MAF = 0.032 and a type-1-error of alpha = 0.05.
The lead SNP from the 2011 GWAS [6], rs17066096, was also genotyped in all the above-mentioned validation cohorts, excluding France (Table 1). Additionally, data on the most significantly associated SNP in the fine mapping, rs202573, were available from Andalucía, Barcelona and Madrid from our previous study [3], and this SNP was newly genotyped in the Donostia and Germany cohorts.

Selection and Genotyping of SNPs
Ten haplotype-tagging SNPs were selected based on genotype data from the CEU + TSI population (HapMap release #27) using an r 2 cut-off = 0.9 and a MAF = 0.1 (with the Multimarker Tagger Algorithm available on the HapMap website, www.hapmap.org, now defunct). In addition, SNP rs13217897 was selected based on its association with MS as reported by an independent group [2], and rs28385692 was included because of its potential functional effect, given its location in the coding region. Haplotype-tagging SNPs in the Bilbao cohort and rs28385692 in the Donostia cohort were genotyped using the iPLEX Sequenom MassARRAY platform in the Spanish National Genotyping Center (CEGEN, Santiago de Compostela, Spain, www.cegen.es). Genotyping of rs17066096, rs202573 and rs28385692 in the validation cohorts was performed using Taqman ® Genotyping Assays, following the manufacturer's instructions. The genotyping efficiency was above 95% in all datasets, and all SNPs were in Hardy-Weinberg equilibrium in controls in each dataset.

Statistical Analysis
The statistical analysis was performed using PLINK v1.07 (http://zzz.bwh.harvard.edu/plink/) [24]. For the association analysis, a logistic regression with an additive model, adjusted for sex, was applied. The independence of association signals was assessed using a conditional logistic regression analysis. A chi-squared test was used to check for the Hardy-Weinberg equilibrium. A haplotype analysis was performed using Haploview v.4.2 (https://www.broadinstitute.org/haploview/haploview) [25]. Only the samples that had genotyping data available for all the SNPs were used for the haplotype calculation (for the Bilbao cohort, 375 cases and 441 controls; for the joint datasets, 6545 cases and 5713 controls). Haplotype blocks were defined by the confidence interval method [26]. Permutations were applied (n permutations = 1000) to correct for multiple comparisons in the haplotype analysis. Statistical power was calculated using the CaTS power calculator at www.sph.umich.edu/csg/abecasis/CaTS/ [27]. Secretion levels of Leu16 IL-22BP isoforms compared to those of Pro16 variants were compared used Student's unpaired t-test.

p.Leu16Pro Mutagenesis of the Three IL-22BP Isoforms
The expression plasmid for IL22RA2v1 was constructed as described in our previous work [20], IL22RA2v2 and IL22RA2v3 expression plasmids were purchased from OriGene Technologies (RC219095, Rockville, MD, USA) and GenScript (Ohu00490, Piscataway, NJ, USA), respectively. The p.Leu16Pro mutants of IL-22BPi1, 2 and 3 were generated using the GENEART ® site-directed mutagenesis system (A13282, Invitrogen, Waltham, MA, USA) from the IL22RA2v1, 2 and 3 expression plasmids following the manufacturer's instructions. The site-directed mutagenesis primer design was also done following the manufacturer's instructions. Briefly, both primers contained the desired mutation centrally located and were 100% complementary with no overhangs, and with lengths between 30 and 45 nucleotides. The designed primers, purchased from IDT, were purified by HPLC to increase mutagenesis efficiency. PCR was performed using a Verity thermocycler (Applied Biosystems, Waltham, MA, USA) with the following primers: IL22RA2_p.Leu16Pro_FW: 5 -TCATCAGTTTCTTCCCTACTGGTGTAGCAGG-3 and IL22RA2_p.Leu16Pro _RV: 5 -CCTGCTACACCAGTAGGGAAGAAACTGATGA-3 . The PCR conditions used were: 1 cycle at 37 • C for 20 min and 94 • C for 45 s, 18 cycles at 94 • C for 45 s, 57 • C for 45 s and 72 • C for 6 min, 1 cycle at 72 • C for 10 min and a final holding stage at 4 • C. The resulting constructs were sequenced to confirm point mutations of IL22RA2 variant constructs. The p.Leu16Pro mutant plasmids were transformed into DH5α-T1R E. coli competent cells following the manufacturer's instructions. Plasmids were purified with EndoFree Plasmid Kit (Qiagen, Hilden, Germany), quantified by a NanoDrop Spectrophotometer (ThermoFisher Scientific, Waltham, MA, USA) and visually assessed via agarose gel electrophoresis. and 5% CO 2 in a humidified incubator. Cells were seeded on a 24-well plate at 500 µL/well at a density of 3 × 10 4 cells/well. Cells were transfected with the indicated expression plasmids when they reached 60%-70% confluency using MACSfectin Reagent (130-098-412, Miltenyi, Bergisch Gladbach, Germany). Twenty-four hours after transfection, conditioned media were collected and cells washed three times with cold PBS prior cell lysis in RIPA buffer (25 mM Tris-HCl, 150 mM NaCl, pH 7.6; plus 1% N-40, 1% sodium deoxycholate, and 0.1% SDS) in the presence of protease inhibitors (11697498001, Roche, Basel, Switzerland). All cell lysates were quantified for total protein using BCA kit (23225, ThermoFisher Scientific); equal amounts of total protein were resuspended in reducing loading buffer and resolved on SDS-PAGE for further immunoblotting. For acetone precipitation of conditioned media, 4 h prior media collection, cells were carefully washed five times with pre-warmed serum-free medium (SFM; 12-764Q, Lonza, Basel, Switzerland) to remove serum proteins and SFM supplemented with L-glutamine (G5792, Sigma) was added for a further 4 h. Conditioned SFM were acetone precipitated with four volumes of ice-cold acetone and incubated on ice for 10 min followed by centrifugation at 21,000× g a 4 • C. Pellets were resuspended in reducing loading buffer and resolved on SDS-PAGE for further immunoblotting. The antibodies used in this study are the following: anti-FLAG (1:1000; 2043-1-AP, Proteintech, Rosemont, IL, USA); anti-IL-22BP (1:1000; AF1087 and BAF1087, R&D Systems, Minneapolis, MN, USA); anti-IL-22BP (1:1000; ab133965, Abcam, San Francisco, CA, USA); anti-tubulin (1:1000; A01490, GenScript, Piscataway, NJ, USA) and all HRP-conjugated secondary antibodies were purchased from Jackson ImmunoResearch. IL-22BP ELISA capable of detecting the three isoforms was performed as previously described [20].

Replacement of the SP of IL-22BPi2 with the IL17A SP
The IL17SP_IL22RA2v2 construct was generated from two overlapping fragments; the first one containing the terminal SgfI restriction site followed by the sequence of the SP of IL-17A and the beginning of IL-22BPi2 mature protein, and the second one, containing the end of IL-17A SP followed by IL-22BPi2 mature protein flanked by Mlu restriction site. The overlapping fragment was amplified and digested with SgfI and MluI restriction enzymes and cloned into pCMV6-Entry vector. Primers for amplification of the IL17A_Signal peptide fragment were: FW: 5 -GCCGCGATCGCCATGACTCCTGGGAAGACC-3 and RV:5 -AGGCTTCAGAGACTCATGCG TTGACTGAGTGATTGTGATTCCTGCCTTCACTATGGCCTCCAGGCTC-3 . Primers for amplification of the IL22RA2v2 mature fragment were: FW 5 -GAGCCTGGAGGCCATAGTGAAGGCAGGAATCA CAATCACTCAGTCAACGCATGAGTC TCTG-3 and RV: 5 -CGTACGCGTTGGAATTTCCACACA TCTCTC-3 .

Flow Cytometry Analysis
HEK293 cells were transfected with expression vectors for wild-type or p.Leu16Pro mutant IL-22BPi2, collected after 24 h and washed with flow cytometry buffer (FC Buffer; PBS, 0.5% BSA and 2mM EDTA, pH 7.2). Single cell suspensions were fixed with 4% paraformaldehyde for 10 min at room temperature followed by permeabilization with 90% of ice-cold methanol for 30 min at 4 • C. Cells were blocked with 1% BSA for 15 min at RT, incubated with anti-IL-22BP primary antibody (1:250; 66190, Proteintech) for 30 min at room temperature and followed by another 30-minute incubation period with anti-mouse-FITC conjugated secondary antibody (1:500; AMI4608, ThermoFisher) protected from light. Two washes were performed after each step with FC buffer. Immunostained cells were analyzed using a MACSQuant Analyzer (Miltenyi).

The p.Leu16Pro coding SNP rs28385692 is Associated with Risk for Multiple Sclerosis
A total of 15 SNPs in a~100-kb interval around IL22RA2 were analyzed in the Bilbao cohort. Of these, five showed nominally significant (p < 0.05) association with MS ( Figure 1, Table 2). The most strongly associated SNP was rs202573 (OR = 1.27, p = 0.007), which had already been the index SNP in the initial screening [3]. This intronic SNP is in weak LD in the Bilbao dataset (r 2 = 0.06, D' = 0.31) with the most significant GWAS-derived top SNP found at the IL22RA2 locus, i.e., rs17066096 [1,[4][5][6]. Rs17066096 showed no association in the Bilbao dataset. Logistic regression analysis of all SNPs conditioning on rs202573 did not reveal any other independent SNP signal reaching statistical significance including rs17066096, suggesting that the main effect observed in this population is driven SNPs, depicted with dots in different colors depending on r 2 values with respect to the index SNP rs202573, are plotted as a function of their log-converted p-value (left Y axis) and their position on chromosome 6 according to hg19 assembly of the human genome (X axis). The recombination rate across the locus is provided on the right Y axis. The red line represents the significance threshold (p = 0.05). SNPs genotyped in the primary screening [3] are shown in italics. rs202573, which was genotyped both in the previous and in the present study, is underlined. Interestingly, an infrequent non-synonymous SNP (rs28385692; changes Leu to Pro at position 16 of the SP shared by the three IL-22BP isoforms) was at the limit of statistical significance in the Bilbao dataset ( Figure 1 and Table 2; p = 0.05, OR = 1.972, 95% CI = 0.983-3.954). The risk allele of this SNP (C) is comparatively rare (MAF = 0.015). Notably, haplotype analyses of all genotyped SNPs revealed that among all the haplotypes that were present with a frequency higher than 0.5% in the Bilbao cohort, the C allele was present only in one haplotype, which contained the risk alleles of four additional SNPs that displayed significant associations in the Bilbao dataset ( Figure A1). The association of rs28385692 with MS strengthened when adding 250 controls and 572 cases from nearby Donostia to our analysis, also located in the Basque Country (p = 0.03, OR = 1.88, 95% CI = SNPs, depicted with dots in different colors depending on r 2 values with respect to the index SNP rs202573, are plotted as a function of their log-converted p-value (left Y axis) and their position on chromosome 6 according to hg19 assembly of the human genome (X axis). The recombination rate across the locus is provided on the right Y axis. The red line represents the significance threshold (p = 0.05). SNPs genotyped in the primary screening [3] are shown in italics. rs202573, which was genotyped both in the previous and in the present study, is underlined. Interestingly, an infrequent non-synonymous SNP (rs28385692; changes Leu to Pro at position 16 of the SP shared by the three IL-22BP isoforms) was at the limit of statistical significance in the Bilbao dataset ( Figure 1 and Table 2; p = 0.05, OR = 1.972, 95% CI = 0.983-3.954). The risk allele of this SNP (C) is comparatively rare (MAF = 0.015). Notably, haplotype analyses of all genotyped SNPs revealed that among all the haplotypes that were present with a frequency higher than 0.5% in the Bilbao cohort, the C allele was present only in one haplotype, which contained the risk alleles of four additional SNPs that displayed significant associations in the Bilbao dataset ( Figure A1). The association of rs28385692 with MS strengthened when adding 250 controls and 572 cases from nearby Donostia to our analysis, also located in the Basque Country (p = 0.03, OR = 1.88, 95% CI = 1.08-3.282). Next, we validated this finding in independent Spanish and European cohorts comprising a total of 8960 cases and 7613 controls. The joint analysis of the discovery and validation cohorts confirmed the association of rs28385692 with MS (p = 3.6 × 10 −4 , OR = 1.262, 95% CI = 1.11-1.434 (Figure 2a).  To further delineate the findings of this mapping effort, our original index SNP, rs202573, and the GWAS-derived top SNP, rs17066096 [4,6], were also genotyped in all validation datasets except France. Upon combining all available datasets, we found a significant association of rs17066096 (p = 9.9 × 10 −5 , OR = 1.11, 95% CI = 1.054-1.172), but not rs202573 (p = 0.12, OR = 1.04, 95% CI = 0.99-1.09) (Figure 2b,c). Logistic regression analysis on the combined datasets revealed that the associations of rs17066096 and rs28385692 are statistically independent, since conditioning on one SNP did not abolish association of the other ( Table 3). None of the haplotypes formed by these three SNPs increased the statistical evidence for association compared to the single SNPs in the combined datasets, but, as seen in the dataset from Bilbao, the C allele of the infrequent SNP rs28385692 appeared only in the same haplotype as the risk alleles of the other two SNPs ( Figure A2). We next attempted to define which of the associated SNPs were most likely to have functional effects based on in silico predictions. To this purpose, we analyzed the five significant variants from the original Bilbao SNP screen and SNPs that were significantly associated in the validation effort, i.e., rs17066096 and rs28385692, using VEP and RegulomeDB ( Table 4). As this study was based on a haplotype-tagging method, the observed association signals were considered to be representative of SNPs in linkage disequilibrium (LD) with the associated markers, and therefore, proxies for the associated SNPs were also included in the analysis. Neither rs202573, the most associated SNP in the Bilbao dataset, nor the only proxy of this SNP which met the r 2 threshold (r 2 = 0.8) to be included in the analysis, were predicted to overlap with any regulatory feature. rs17066096 did not show evidence to be functional either, but two of its perfect proxies (r 2 = 1), rs17066063 and rs62420820, displayed a moderate regulatory potential based on RegulomeDB. rs17066063 was the SNP with the highest regulatory potential based on concordance between RegulomeDB (score 3a) and VEP (consequence: TF binding site variant). rs28385692 was predicted to lie within a regulatory region, although the level of evidence supporting this is modest (score 5 according to RegulomeDB). To further delineate the findings of this mapping effort, our original index SNP, rs202573, and the GWAS-derived top SNP, rs17066096 [4,6], were also genotyped in all validation datasets except France. Upon combining all available datasets, we found a significant association of rs17066096 (p = 9.9 × 10 −5 , OR = 1.11, 95% CI = 1.054-1.172), but not rs202573 (p = 0.12, OR = 1.04, 95% CI = 0.99-1.09) (Figure 2b,c). Logistic regression analysis on the combined datasets revealed that the associations of rs17066096 and rs28385692 are statistically independent, since conditioning on one SNP did not abolish association of the other ( Table 3). None of the haplotypes formed by these three SNPs increased the statistical evidence for association compared to the single SNPs in the combined datasets, but, as seen in the dataset from Bilbao, the C allele of the infrequent SNP rs28385692 appeared only in the same haplotype as the risk alleles of the other two SNPs ( Figure A2). We next attempted to define which of the associated SNPs were most likely to have functional effects based on in silico predictions. To this purpose, we analyzed the five significant variants from the original Bilbao SNP screen and SNPs that were significantly associated in the validation effort, i.e., rs17066096 and rs28385692, using VEP and RegulomeDB (Table 4). As this study was based on a haplotype-tagging method, the observed association signals were considered to be representative of SNPs in linkage disequilibrium (LD) with the associated markers, and therefore, proxies for the associated SNPs were also included in the analysis. Neither rs202573, the most associated SNP in the Bilbao dataset, nor the only proxy of this SNP which met the r 2 threshold (r 2 = 0.8) to be included in Cells 2020, 9, 175 9 of 23 the analysis, were predicted to overlap with any regulatory feature. rs17066096 did not show evidence to be functional either, but two of its perfect proxies (r 2 = 1), rs17066063 and rs62420820, displayed a moderate regulatory potential based on RegulomeDB. rs17066063 was the SNP with the highest regulatory potential based on concordance between RegulomeDB (score 3a) and VEP (consequence: TF binding site variant). rs28385692 was predicted to lie within a regulatory region, although the level of evidence supporting this is modest (score 5 according to RegulomeDB). Table 4. Functional predictions of associated SNPs and proxies. SNPs with significant associations in the mapping exercise or in the discovery + validation cohorts and their proxies (r 2 > 0.8 in 1000 Genomes Phase III CEU population) were assessed using VEP and RegulomeDB. 1 Minor allele frequency is based on European populations in the 1000 Genomes Project Phase III.

In Silico Analysis of the p.Leu16Pro Variant
We evaluated in silico whether rs28385692 is functionally neutral or deleterious by using computational tools based on different approaches [42]; i.e., PredictSNP [31], Meta-SNP [32] and Ensembl VEP tool [29,33]. A summary of the features of each tool used in this study is represented in Table A3. PredictSNP indicated that the overall result of the p.Leu16Pro transition in the SP was deleterious for all three isoforms ( Figure A3a). However, Meta-SNP [32] anticipated a neutral effect for the p.Leu16Pro point mutation in all IL-22BP isoforms ( Figure A3b). Finally, the pathogenicity prediction tools included in Ensembl VEP tool predicted the substitution not to be deleterious except for a sub-analysis based on the Mutation Assessor tool available in Ensembl VEP, which attributed a medium level of functional impact to this variant ( Figure A3c). Thus, the overall results obtained with the computational tools were not robust enough to unequivocally assign a functional effect to the p.Leu16Pro variant. As this variant occurs in the SP of IL-22BP protein, its potential effect on structural aspects pertinent to SP biological function per se was assessed in more detail by in silico methods. Coding SNPs in signal peptides may alter translocation efficiency, cleavage sites, as well as post-cleavage events [43]. We analyzed the charge, the hydrophobicity, and the helix-breaker amino acid residues comprised in the SP of IL-22PB, as well as the modifications introduced by the p.Leu16Pro variant (Figure 3). The Leu16 residue is one of three leucines that contribute to the IL-22BP SP hydrophobic core H-region, a stretch of hydrophobic amino acids with propensity to form a single alpha-helix. However, its replacement, Pro16, is a helix-breaker neutral residue with restricted conformational flexibility that has no free hydrogen to contribute to helix stability. Results obtained with various signal peptide cleavage site and secondary structure prediction tools are represented in Figure 3a. While SignalP-3.0 [34] and PSIPRED [36] did not predict any change in the cleavage site, Phobius [35] predicted the p.Leu16Pro polymorphism to shift the signal peptide cleavage site from the 21st to 20th residue. The three software applications coincided in predicting a shortening of the hydrophobic core of the signal peptide. Other prediction tools were applied as well and these predicted similar structural effects ( Figure A4).

The p.Leu16Pro Variant Decreases Secretion of IL-22BPi1, IL-22BPi2 and IL-22BPi3
To experimentally verify the effect of this variant on secretion, we performed site-directed mutagenesis to change Leu16 to Pro in the signal peptide of the three wild-type IL-22BP isoforms cloned in expression vectors reported before [20]. HEK293 cells were individually transfected with these vectors, and both intracellular and secreted levels of IL-22BPi1, IL-22BPi2 and IL-22BPi3 were measured by ELISA. Compared to Leu16, Pro16 strongly decreased the secreted levels of each isoform (Figure 4a). This was mirrored by nonsignificant trends towards decreased intracellular levels of IL-22BPi1 and IL-22BPi2. Cell lysates as well as acetone precipitates of conditioned medium of HEK293 cells transfected to produce Leu16 or Pro16 forms of IL-22BPi1 and IL-22BPi2 were analyzed by immunoblot and this revealed decreased intracellular (51 kDa for IL-22BPi1 and 48 kDa for IL-22BPi2) and secreted (56 kDa for IL-22BPi2) levels of the mutant forms (Figure 4b), while transfection efficiencies as measured by RT-qPCR through quantification of mRNA levels were similar (Figure 4c). These observations are compatible with the effects predicted by the above in silico analysis, suggesting that the Pro16 variant renders import of the precursor IL-22BP isoforms into the endoplasmic reticulum less efficient. Flow cytometry analysis also revealed a decrease of 41% to 33% in IL-22BP+ cells following transfection of HEK293 with Leu16 or Pro16 vectors encoding IL-22BPi2, respectively ( Figure A5). We also assessed the efficiency of the native signal peptide of IL-22BPi2 in facilitating secretion by replacing it with that of IL-17, an efficiently secreted protein [20]. This modification resulted in much higher levels of both intracellular and secreted IL-22BPi2 ( Figure 5), suggesting that the wild-type IL-22BP signal peptide is less efficiently engaged by the co-translational targeting and translocation machinery across the ER membrane than that of IL-17, which may affect overall IL-22BP biogenesis [44]. obtained with various signal peptide cleavage site and secondary structure prediction tools are represented in Figure 3a. While SignalP-3.0 [34] and PSIPRED [36] did not predict any change in the cleavage site, Phobius [35] predicted the p.Leu16Pro polymorphism to shift the signal peptide cleavage site from the 21st to 20th residue. The three software applications coincided in predicting a shortening of the hydrophobic core of the signal peptide. Other prediction tools were applied as well and these predicted similar structural effects ( Figure A4).

Discussion
In this study, we identified an infrequent functional variant (rs28385692) in IL22RA2, which confers risk to MS independently of the major GWAS-derived signal, and we were able to show that this variant has immediate functional effects by reducing secreted IL-22BP levels in transfected HEK293 cells. . For AP, the immunoblot membrane was subjected to longer exposure times. (c) HEK293 cells were transiently transfected with the indicated expression vectors (EV denotes empty vector), 24 h later cells were lysed and RNA purified. Intracellular IL-22BP protein was immunoaffinity-purified with FLAG agaroses and detected by WB following FLAG purification and in cell lysates (CL) and pass through fraction (PT). GRP94 detection and Ponceau staining served as loading controls. Transfection efficiency was measured by IL22RA2 RT-qPCR relative to the housekeeping gene ACTB. Mean ± SEM of three technical replicates. Note that as previously observed [20], immunoreactive bands corresponding to intracellular IL-22BP isoforms appear as a series of 43 to 56 kDa bands due to differential N-glycosylation, with secreted IL-22BPi2 gaining ~8 kDa (56 vs. 48 kDa) due to complex N-glycosylation. (a) HEK293 cells were transfected with the indicated expression plasmids, 24 h later cells were lysed and the conditioned medium collected. Intracellular and secreted IL-22BP isoform protein levels were measured by ELISA (mean ± SEM; n = 3; p-values by unpaired t-test). (b) HEK293 cells were transfected with the indicated expression plasmids, 24 h later cells were lysed and the conditioned medium was subjected to acetone precipitation (AP). Both cell lysates (CL) and AP were resolved by SDS-PAGE under non-reducing conditions and immunoblotted against FLAG (Ponceau staining served as loading control). For AP, the immunoblot membrane was subjected to longer exposure times. (c) HEK293 cells were transiently transfected with the indicated expression vectors (EV denotes empty vector), 24 h later cells were lysed and RNA purified. Intracellular IL-22BP protein was immunoaffinity-purified with FLAG agaroses and detected by WB following FLAG purification and in cell lysates (CL) and pass through fraction (PT). GRP94 detection and Ponceau staining served as loading controls. Transfection efficiency was measured by IL22RA2 RT-qPCR relative to the housekeeping gene ACTB. Mean ± SEM of three technical replicates. Note that as previously observed [20], immunoreactive bands corresponding to intracellular IL-22BP isoforms appear as a series of 43 to 56 kDa bands due to differential N-glycosylation, with secreted IL-22BPi2 gaining~8 kDa (56 vs. 48 kDa) due to complex N-glycosylation. SNP rs28385692, located in exon 2 of IL22RA2, changes the amino acid in position 16 of the protein from leucine to proline. This change occurs in the SP shared by the three IL-22BP isoforms produced by IL22RA2. In silico prediction methods diverged in their capacity to assign a functionally relevant effect on protein function to this variant. However, various in silico structural assessment methods coincided in predicting a shortening of the hydrophobic H-region of the IL-22BP SP by the 16P residue that may influence the precise location of the cleavage site with the mature portion of the IL-22BP proteins. Importantly, in this study, we went beyond the in silico predictions by assessing the effect of rs28385692 in vitro and showed experimentally using transfected cells that the p.Leu16Pro mutation significantly reduced the secreted levels of the IL-22BPi1, IL-22BPi2 and IL-22BPi3 isoforms. Specifically, the risk allele of p.Leu16Pro lowers secretion levels of each isoform to around 50%-60% of those of the respective wild-type isoforms as quantified using ELISA. Both cell lysates (CL) and AP were resolved by SDS-PAGE under non-reducing conditions and immunoblotted against FLAG, IL-22BP and using tubulin as loading control. Intracellular IL-22BP reactive bands relative to tubulin ones were scanned and represented as fold change to IL-22BP wild-type (mean ± SD, n = 2).

Discussion
In this study, we identified an infrequent functional variant (rs28385692) in IL22RA2, which confers risk to MS independently of the major GWAS-derived signal, and we were able to show that this variant has immediate functional effects by reducing secreted IL-22BP levels in transfected HEK293 cells. SNP rs28385692, located in exon 2 of IL22RA2, changes the amino acid in position 16 of the protein from leucine to proline. This change occurs in the SP shared by the three IL-22BP isoforms produced by IL22RA2. In silico prediction methods diverged in their capacity to assign a functionally relevant effect on protein function to this variant. However, various in silico structural assessment methods coincided in predicting a shortening of the hydrophobic H-region of the IL-22BP SP by the 16P residue that may influence the precise location of the cleavage site with the mature portion of the IL-22BP proteins. Importantly, in this study, we went beyond the in silico Figure 5. Native signal peptide of IL-22BP does not efficiently mediate secretion of IL-22BPi2. HEK293 cells were transfected with IL-17, IL-22BPi2, and IL-17SP_IL-22BPi2 expression plasmids, and 24 h later, cells were lysed and the conditioned medium subjected to acetone precipitation (AP). Both cell lysates (CL) and AP were resolved by SDS-PAGE under non-reducing conditions and immunoblotted against FLAG, IL-22BP and using tubulin as loading control. Intracellular IL-22BP reactive bands relative to tubulin ones were scanned and represented as fold change to IL-22BP wild-type (mean ± SD, n = 2).
The literature documents several examples of similar changes in signal peptides that have been proven to affect secretion and are related to human diseases or traits: a Leu->Pro mutation in the signal peptide of COL5A1 (Alpha 1 Type V Collagen), a subunit of type V collagen, caused retention of this subunit in the endoplasmic reticulum and was associated with Ehlers-Danlos syndrome, a heritable connective tissue disease, in two unrelated families [45]. In addition, a common SNP causing a Leu->Pro mutation in the prohormone region of NPY (Neuropeptide Y) caused an increased secretion of this protein in chromaffin cells [46]. The inverse case, i.e., a change from proline to leucine, has also been reported to affect protein secretion: a Pro->Leu mutation in the signal peptide of DSPP (Dentin Sialophosphoprotein) resulted in a defective secretion of this protein and was associated with dentinogenesis imperfecta type III in a Korean family [47]. Finally, a Thr->Pro variant associated with coronary artery disease risk, located in the SP of the lysosomal acid lipase gene (LIPA), yields reduced LIPA protein levels and activity due to enhanced degradation [48].
The logistic regression analysis showed that associations of rs28385692 and the main IL22RA2 GWAS SNP, rs17066096, with MS may be statistically independent. This suggests that there are at least two different variants causing association with MS within the IL22RA2 locus: one presumably corresponding to rs28385692 and the other one to a SNP in LD with rs17066096. rs17066096 itself is unlikely to be a causative variant, based on its intergenic location and its missing functional effect in our extensive in silico annotations. Recently, Lill and colleagues made an attempt to discover possible causal variants in LD with rs17066096. They performed an in silico and in vitro analysis of SNPs in the 3 UTR of IL22RA2 that were at least in moderate LD with rs17066096 (r 2 > 0.3) predicted to affect the binding of micro-RNAs (miRNAs; [49]). Although they successfully identified one SNP in a miRNA binding area, they did not find allele-specific differences on this binding; therefore, association of rs17066096 with MS remains unexplained. Despite these non-confirming functional data of Lill et al., rs17066063 might be a reasonable candidate for future functional studies given its location in a regulatory region of the gene.
Our study has several limitations. Firstly, in certain sub-analyses (e.g., individual populations), our study may have been underpowered to detect modest effects. Secondly, we restricted our analyses to individuals of European descent, which was self-reported. Thus, residual confounding due to undetected population stratification may have impacted some of our results and may account for some of the heterogeneity of the observed effect estimates. However, the effect size for the GWAS-derived SNP rs17066096 estimated in our datasets (OR = 1.11) was very similar to the one described in the original GWAS (OR = 1.14) [4][5][6], suggesting that the impact of population substructure on our results is minor. Lastly, observing both genetic risk association and an effect of the genotypes on gene expression or quantitative protein production does not imply causality. Therefore, additional studies are necessary to characterize the exact molecular mechanisms of the IL22RA2 association signal in MS.
In summary, we successfully identified a comparatively rare genetic variant in IL22RA2 that is strongly-and likely independently of the primary GWAS signal-associated with MS. It causes an amino acid change in the signal peptide of IL-22BP and correspondingly, functional data generated for this study suggest an inhibitory effect on the trafficking of IL-22BP. Functional research on this variant will certainly provide insight about the role of this gene in MS and a better understanding of the still fairly unexplored function of IL-22BP isoforms in the immune system.  Table A2. Association values of the haplotypes in blocks. LD blocks were calculated with the confidence interval method using Haploview [25,26].  Figure A1. Linkage disequilibrium structure and haplotypes of IL22RA2. On the left, LD plot of the 15 SNPs analyzed in the Bilbao cohort. Blocks were calculated using the confidence interval algorithm implemented in Haploview [25,26]. Only samples genotyped for all SNPs were considered for the haplotype calculation. The numbers inside the squares indicate r 2 values, and darker shades of gray represent higher degrees of linkage disequilibrium. On the right, all the haplotypes that were present in the Bilbao cohorts considering all SNPs in a single LD block are represented. The only haplotype containing the risk (C) allele of the non-synonymous SNP rs28385692 is boxed. SNPs with significant associations and rs28385692 are indicated on top of the right panel.  SNPs analyzed in the Bilbao cohort. Blocks were calculated using the confidence interval algorithm implemented in Haploview [25,26]. Only samples genotyped for all SNPs were considered for the haplotype calculation. The numbers inside the squares indicate r 2 values, and darker shades of gray represent higher degrees of linkage disequilibrium. On the right, all the haplotypes that were present in the Bilbao cohorts considering all SNPs in a single LD block are represented. The only haplotype containing the risk (C) allele of the non-synonymous SNP rs28385692 is boxed. SNPs with significant associations and rs28385692 are indicated on top of the right panel.  [25,26]. Only samples genotyped for all SNPs were considered for the haplotype calculation. The numbers inside the squares indicate r 2 values, and darker shades of gray represent higher degrees of linkage disequilibrium. On the right, all the haplotypes that were present in the Bilbao cohorts considering all SNPs in a single LD block are represented. The only haplotype containing the risk (C) allele of the non-synonymous SNP rs28385692 is boxed. SNPs with significant associations and rs28385692 are indicated on top of the right panel.