HLA-G Haplotypes Are Differentially Associated with Asthmatic Features

Human leukocyte antigen (HLA)-G, a HLA class Ib molecule, interacts with receptors on lymphocytes such as T cells, B cells, and natural killer cells to influence immune responses. Unlike classical HLA molecules, HLA-G expression is not found on all somatic cells, but restricted to tissue sites, including human bronchial epithelium cells (HBEC). Individual variation in HLA-G expression is linked to its genetic polymorphism and has been associated with many pathological situations such as asthma, which is characterized by epithelium abnormalities and inflammatory cell activation. Studies reported both higher and equivalent soluble HLA-G (sHLA-G) expression in different cohorts of asthmatic patients. In particular, we recently described impaired local expression of HLA-G and abnormal profiles for alternatively spliced isoforms in HBEC from asthmatic patients. sHLA-G dosage is challenging because of its many levels of polymorphism (dimerization, association with β2-microglobulin, and alternative splicing), thus many clinical studies focused on HLA-G single-nucleotide polymorphisms as predictive biomarkers, but few analyzed HLA-G haplotypes. Here, we aimed to characterize HLA-G haplotypes and describe their association with asthmatic clinical features and sHLA-G peripheral expression and to describe variations in transcription factor (TF) binding sites and alternative splicing sites. HLA-G haplotypes were differentially distributed in 330 healthy and 580 asthmatic individuals. Furthermore, HLA-G haplotypes were associated with asthmatic clinical features showed. However, we did not confirm an association between sHLA-G and genetic, biological, or clinical parameters. HLA-G haplotypes were phylogenetically split into distinct groups, with each group displaying particular variations in TF binding or RNA splicing sites that could reflect differential HLA-G qualitative or quantitative expression, with tissue-dependent specificities. Our results, based on a multicenter cohort, thus support the pertinence of HLA-G haplotypes as predictive genetic markers for asthma.

inTrODUcTiOn Human leukocyte antigen-G (HLA-G), a member of the HLA class Ib (non-classical) family, modulates natural killer (NK) cells and cytotoxic T-lymphocyte-mediated activity as well as B-lymphocyte proliferation. HLA-G is also involved in epithelial cell differentiation (1)(2)(3). Human leukocyte antigen-G is expressed by a restricted pattern of cells including human bronchial epithelium cells (HBEC) and displays interindividual variation expression in physiological conditions. The mechanisms responsible for tissue-specific expression and interindividual HLA-G expression variation are of utmost importance since HLA-G expression is associated with many diseases involving immune system surveillance and inflammatory disease, such as asthma (3,4).
Asthma is characterized by structural changes in the airways with abnormalities of the epithelium and concomitant recruitment, survival, and activation of inflammatory cells. In asthma, soluble HLA-G (sHLA-G) expression was reported to be higher in bronchoalveolar lavage (BAL) from asthmatic individuals and in both plasma and serum from atopic asthmatic children, whereas others reported no difference in non-allergic asthmatics (5)(6)(7)(8). In BAL from asthmatic patients, however, sHLA-G expression was reported to be inversely correlated with markers of asthma inflammation (9). Our group reported an impairment of full-length protein HLA-G expression in HBEC from asthmatic patients and abnormal expression of alternatively spliced membrane-bound (HLA-G1 to 4) and soluble (HLA-G5 to 7) forms (10). Both HLA-G1 and -G5 isoforms reduced NK cell cytotoxicity with an additive effect, and HLA-G2, -G4, and -G6 were described as functional, whereas conflicting results were published concerning HLA-G3 (11)(12)(13)(14).
Human leukocyte antigen-G tissue-specific expression is reported to be driven by a 12 kb upstream enhancer site, defined as cis-regulatory Enhancer L (15). Concerning interindividual variation expression, HLA-G*01:04 has been associated with high HLA-G secretion and G*01:03 and G*01:05N with low HLA-G secretion (3). In the 3′ untranslated region (UTR), conflicting results were reported concerning association of HLA-G expression and the ex8 ins/del polymorphism (16); the 3142G allele was identified as a binding site for the microRNA-152 family (17) and 3187A/G and 3196C/G were reported to influence mRNA stability. 3′UTR haplotypes differentially influenced HLA-G expression in gene reporter assays conducted in vitro (18). However, in silico analysis of 3′UTR haplotypes suggested that most microRNAs targeting HLA-G in a strong or specific manner are not influenced by 3′UTR single-nucleotide polymorphisms (SNPs) and that 3′UTR haplotypes display equivalent means of microRNA binding (19). Polymorphism in the 5′ upstream regulatory region (URR) within regulatory elements, such as -1305G>A, -964G>A 725C>G>T, -716G>T and -56C>T, were reported to impact HLA-G expression (3,20).
Human leukocyte antigen-G protein dosage has its drawbacks, and therefore, many teams have studied the association between HLA-G SNP and pathological situations (3). Few studies have focused on HLA-G genetic haplotypes as predictive biomarkers (21). HLA-G haplotypes are defined as SNPs being in Linkage Disequilibrium; Castelli et al. first described a restricted number of HLA-G haplotypes using SNPs in the 5′URR, 3′UTR, and coding regions in a Brazilian population (22). This low variability suggests a stabilizing selective effect on HLA-G haplotypes which could reflect a synergic effect of HLA-G SNPs. We previously characterized a deleterious association of two HLA-G haplotypes, HLA-G*01:04-UTR3 and HLA-G*01:06-UTR2, in lung transplantation outcome and in lung inflammatory diseases such as cystic fibrosis (23). Although the HLA-G gene was identified as an asthma susceptibility gene and SNPs -964 and 3142 were associated with asthma (17,(24)(25)(26)(27), no study so far has described the haplotypes effect, shown to be more relevant than SNPs to predict sHLA-G dosage (28,29), in asthma clinical studies.
We aimed to investigate the association of HLA-G haplotypes and sHLA-G peripheral expression with clinical features of asthma. We then aimed to analyze phylogenetic relationships, variations at potential transcription factor (TF) binding sites and alternative splicing sites between HLA-G haplotypes that could account for their differential association with asthmatic features.

MaTerials anD MeThODs samples
This study was carried out in accordance with the French Public Health Code (art L1221-1). All subjects gave written informed consent in accordance with the Declaration of Helsinki. Three hundred and thirty volunteer bone marrow donors from France were recruited, blood donations were collected, and DNA was extracted from a 200 µl total blood sample using the QIAmp Blood DNA Mini kit (Qiagen, France) in accordance with BSL-2 practices. A medical interview was carried out prior to blood donation to exclude donors with medical contraindications such as severe asthma, respiratory failure, or allergies. Five hundred and ninety asthmatic patients enrolled in a French asthmatic multicenter cohort (Cohorte Obstruction Bronchique et Asthme, INSERM Hôpital Bichat, CNIL 28/01/2008) were included in the study. DNA samples with known genotype were used as control. Biological samples included 590 samples of DNA, 582 of serum, and a further 528 of serum collected 1 year after the first sample. Patient data were categorized into biological and clinical data and asthma severity according to the Guidelines for the Diagnosis and Management of Asthma, NHLBI 2007. Healthy individuals (HI) and asthmatic patient characteristics are detailed in Table 1.

serum shla-g Protein Measurement
Soluble HLA-G (sHLA-G), i.e., peripheral HLA-G expression, was quantified in sera from 582 asthmatic patients and 528 sera collected 1 year after the first sample. All serum samples were stored at −80°C prior to analysis. Measurement of soluble isoforms sHLA-G1 and sHLA-G5 was performed in duplicate on serum samples using the ELISA test (Biovendor ® , Czech Republic; capture antibody: MEM-G/9) according to the manufacturer's protocol.
association of HLA-G alleles and haplotypes with asthmatic Features and Peripheral shla-g Dosage Missing data led to the exclusion of the concerned sample from further analyses. No multiple imputation was used. Genetic data frequencies between HI and asthmatic patients were compared with Chi-square tests. Associations between patients' asthmatic features, sHLA-G expression, and genetic polymorphism were tested with logistic regressions when the predicted variable was categorical and general linear regression when the predicted variable was quantitative. Univariate analyses were performed on each of these asthmatic features: eosinophil count (/mm 3 ), history of near-fatal asthma, asthma control, asthma control test score, asthmatic exacerbations within the last 12 months, emergency consultation for asthma within the last 12 months, hospitalization within the last 12 months because of asthma, intensive care unit stays, according to biological parameters [gender, body mass index (BMI), occupational exposure, tobacco use, allergies, maternal asthma, and paternal asthma], HLA-G predictors (HLA-G alleles, HLA-G regulatory regions typed at low resolution, and HLA-G regulatory regions typed at high resolution), and sHLA-G serology and differences between sHLA-G dosage results in sera collected at a 1-year interval. Univariate analyses were also performed on sHLA-G dosage and differences between sHLA-G dosage results in sera collected at a 1-year interval according to clinical and biological parameters such as gender, height, weight, BMI, occupational exposure, tobacco use, allergies, associated pathologies, maternal asthma, paternal asthma, treatment for asthma or for allergy (CO intake, anticholinergics, anti-IgE, anti-leukotrienes, theophylline, antihistamines, anti-RGO, antihypertensives, statins, and fibrates), month of serum collection, and genetic parameters (HLA-G alleles, HLA-G UTR, and HLA-G haplotypes). The False Discovery Rate correction for multiplicity was applied; however, no power analysis was performed. Variables with p values below 0.05 and corrected p values below 0.20 in univariate analysis were further investigated in multivariate analysis. Association results with quantitative variables are expressed as estimates and standard errors (Std error), whereas association results with qualitative variables are expressed as odd ratios (OR) and 95% confidence intervals.

HLA-G Full-length sequence Phylogenetic analysis
HLA-G full-length sequences from position −1983 to +3447 were aligned using the multiple sequence alignment tool MUSCLE in the Molecular Evolutionary Genetics Analysis (MEGA) software version 7.0 (32). Evolutionary analyses were conducted in MEGA7. Evolutionary relationships among HLA-G full-length sequences were inferred using the neighbor-joining (NJ) method with 1,000 replicates. Evolutionary distances were computed using the p-distance method, and units were the number of base differences per site. Analysis involved gorilla HLA-G genomic sequence as an outgroup (accession number CU104658.1).

TF binding site analysis and TF expression in the lungs and in blood
Each SNP in the 5′UTR was analyzed for potential TF binding site variation with the AliBaba 2.1 and ConSite programs. TFs found to have differential binding sites between HLA-G haplotypes were checked for their tissue-specific expression (blood and lungs) Table 2 | HLA-G regulatory region haplotypes identifying number (ID), nucleotide variations, untranslated region (UTR) and allelic association, and frequencies estimated by EM algorithm in healthy individuals (HI) and asthmatic patients [Cohorte Obstruction Bronchique et Asthme (COBRA)] with a total frequency higher than 0.3% estimated by the EM algorithm: ID, nucleotide variations, UTR and allelic association, and frequencies in HI and asthmatic patients (COBRA      T    in NCBI and Tissue-specific Gene Expression and Regulation (TiGER) databases (33). NCBI data included the "HPA RNA-seq normal tissues" project, i.e., RNA-seq performed on 27 different tissue samples from 95 human individuals and the "Illumina bodyMap2 transcriptome" project, i.e., transcription profiling by high-throughput sequencing of RNA from 16 different human tissues.

splicing site Variation analysis
Variation in splicing sites was studied for alleles sharing identical or similar UTR haplotypes using the Human Splicing Finder program (34). This program combines different algorithms to identify and predict the effect of mutations on acceptor and donor splice site motifs.

HLA-G haplotypes are Differentially Distributed in healthy and asthmatic individuals and correlate with asthmatic Features
Next generation sequencing typing was successful in 317 (96%) HI and 582 (98%) asthmatic patients. Thirty-seven SNPs, 1 base deletion, 1 base insertion, and a 14bp insertion were validated through Xegen pipeline analysis in the regulatory regions 5′URR and 3′UTR from position −1983 to +3447; all were previously described (30,35,36). HLA-G sequencing data have been submitted to the NCBI public data base ( Table 2). HLA-G alleles were typed at maximum resolution (eight digits) and regulatory regions (UTR1-8) were typed at low resolution according to eight SNPs in 5′URR and at high resolution based on all 28 polymorphisms in the 5′URR and 12 polymorphisms in the 3′UTR. Fifteen haplotypes with frequency above 0.3% displayed 94% of the cumulated frequencies. Their description and frequencies are detailed in Table 2.
Correlations between HLA-G alleles, UTR, or haplotypes and asthmatic patient data investigated in multivariate analysis are shown in Tables 3A,B. Variables included in the multivariate

Peripheral shla-g expression is not associated with genetics or asthmatic Features
Serum sHLA-G measurements from asthmatic patients were 17.1 (0-493) UI/ml and 17.5 (0-404) UI/ml 1 year later (mean and range). Multivariate analyses revealed no association between peripheral sHLA-G expression and patients' asthmatic features, biological and clinical data, or HLA-G genetics (data not shown).

HLA-G haplotype Variation Defines Four Main groups of sequences
Phylogenetic relationships for the 19 combinations of HLA-G alleles and haplotypes with a frequency above 0.3% were investigated using the NJ method (Figure 1).

TF binding site Variations between HLA-G haplotypes
Thirty-four TFs were identified to bind differentially to each 5′URR sequence by in silico analysis (Tables 5 and 6), their classification is provided in Table 7; their specific expression in blood and lungs according to NCBI and TiGER databases is given in Table 6.
Interestingly, haplotype grouping based on their TF binding site profiles mirrored that observed on the phylogenetic tree (Table 6; Figure 1). Hence, haplotypes displayed different TF binding site profiles, notably in haplotype H23 (UTR3-G*01:04:01/G*01:04:04): the transition G>A at position −1155 impaired progesterone receptor (PR) and SOX5 binding sites. Besides this unique feature, this haplotype displayed similarities with H16 (UTR7-G*01:01:03:03) and H10 (UTR2-G*01:06/G*01:05N/G*01:01:02:XX). Comparing the most common sequences, some SNPs created new binding sites: one to two sites for SRF, one site for CUX1, two to three additional sites for POU2F1, one site for E2F, one site for glucocorticoid receptor, additional sites for SP1 (up to 10), one site for XFD-1, RXR, NR2C1, VDR, CDX1, and an additional site for HNF4. Conversely, some SNPs in H23, H16, and H10 disrupted binding sites for RUNX1, RFX1, and MYB. H16 and H10 differed from H23 and other haplotypes as their SNPs created a fixation site for Lmo complex, FOXI1, and GATA1 and disrupted a binding site for NF-E and TATA binding proteins.
Haplotypes H21, H19, and H20 (associated with UTR5-G*01:03:01:02) distinguished themselves by an additional fixation site for ELK and IRF1, and the disruption of a fixation site for STAT; H19 and H20 displayed a disruption of a IKZF1 fixation site. H21 displayed an increase in REL binding sites. Genetic polymorphism in H21, H19, and H20 created a site for NF-1 and SREBF1 which were also present in haplotypes H23, H16, and H10. H21, H19, H20, H16, and H10 had a genetic polymorphism that created a site for TEF. Haplotypes H21, H19, H20, H23, H16, H10, H04, and H02 displayed a fixation site for C/EBP.

splicing site Variation analysis
Variation in splicing sites was investigated for alleles with identical coding sequences but with variations at a six-or eight-digit level and sharing identical or similar haplotypes. Evolutionary relationships among HLA-G whole sequences were inferred using the neighbor-joining method with 1,000 replicates. Evolutionary distances were computed using the p-distance method, and units were the number of base differences per site. The percentage of trees in which the associated taxa clustered together is shown next to the branches. The sum of branch length equals 0.01616800 base differences per site. The percentage of replicate trees in which the associated taxa clustered together in the bootstrap test is shown next to the branches. The tree was drawn to scale, with branch lengths in the same units (base differences per site) as the evolutionary distances used to construct the phylogenetic tree. Analysis involved gorilla HLA-G genomic sequence as an outgroup (accession number CU104658.1).  Table 8B).

DiscUssiOn
HLA-G genetic polymorphism is extensively described in different human populations worldwide. However, there are only limited data describing the synergic effect of these SNPs, i.e., haplotypes, in clinical studies.
We fully characterized HLA-G haplotypes by NGS from position −1983 to +3447 both in a cohort of 330 HI and in 580 patients from a French asthmatic multicenter cohort. sHLA-G was quantified using the ELISA test (MEM-G/9) in sera from 582 asthmatic patients and in 528 serum samples collected 1 year after the first sample. Associations between patients' asthmatic features, sHLA-G expression, and genetic polymorphism were tested with logistic regressions.   HLA-G haplotypes displayed statistically significant differential distribution between HI and asthmatic patients. Furthermore, multivariate analyses showed significant association of HLA-G haplotypes with asthma inflammation markers (eosinophil count) and asthma severity (history of near-fatal asthma and exacerbation). We did not find any association between peripheral sHLA-G and genetic data in our asthmatic patient cohort; such an association described by us and others in physiological conditions may be masked in an inflammatory context as asthma. The lack of association between peripheral sHLA-G expression and asthmatic features supports the hypothesis that sHLA-G is not overexpressed as a systemic immune response to control local inflammation. These results reinforce the involvement of HLA-G in asthma and suggest that the modulation of local sHLA-G associated with asthma phenotypes and asthma inflammation markers described in Ref. (9) might be partly driven by HLA-G genetic polymorphism.
We also performed descriptive analyses on full length HLA-G sequences: phylogenetic analysis, in silico prediction of TF binding sites, and in silico prediction of splicing sites.
Phylogenetic analysis suggested a specific evolution of HLA-G haplotypes. HLA-G haplotype grouping based on their TF binding site profiles mirrored that observed on the phylogenetic tree.
UTR2-H10 was associated with asthmatic exacerbation within the previous 12 months and eosinophil count, and UTR7 (HLA-G*01:01:03:03-H16) was associated with acute asthma recurrence. These sequences, the closest within the phylogenetic tree and according to TF binding site profiles, could lead to similar expression regulation. These results support the hypothesis of an impaired immunomodulation of the UTR2-HLA-G molecule, previously associated with a worse evolution of cystic fibrosis (23) and pregnancy complications (37,38) for individuals carrying UTR2-HLA-G*01:06.
Transcription factor binding site analysis supported a specific evolution of haplotype H23 that displayed an impaired PR site compared to other haplotypes. Progesterone is thought to play a major role in HLA-G expression because of concomitant requirement of both molecules for embryo implantation and pregnancy success. In vitro studies showed that progesterone enhanced HLA-G expression (39,40). This unique feature of haplotype H23 may lead to an impaired activation of local HLA-G expression and a loss of local inflammation control. Haplotype H23 (Fq = 7.4%) was associated in the asthmatic patient cohort with higher eosinophil count in univariate analysis (p = 0.04); however, multiple correction testing did not reach statistical significance (pc = 0.22). One may speculate, however, that a larger cohort may lead to statistical significance; thus this haplotype might be a good candidate marker for inflammation, as we previously associated UTR3-HLA-G*01:04 with impaired long-term survival following lung transplantation (23). Haplotype H23 is split into H23-UTR3-G*01:04:01 and H23-UTR3-G*01:04:04, respectively, higher and lower in HI than in asthmatic patients (p = 0.044 and p = 0.012). These sequences differ at position 1827G>A in the last positions of exon 4 which, according to in silico predictive analysis, creates an exonic cryptic acceptor site. We previously reported alternatively spliced isoform expression of HLA-G in HBEC from asthmatic patients as compared to HI forms (10). These results support the hypothesis that the mRNA transcript of these HLA-G*01:04 alleles might be processed differently leading to varied protein expression but could also explain discrepancies of clinical studies with HLA-G typing at a four-digit level.
H04-UTR4-G*01:01:01:05 and H02-UTR4-G*01:01:01:05 were both higher in HI than in asthmatic patients (p = 0.015 and p = 0.009). No specific TF binding site was predicted by in silico analysis that could account for a putative protective effect; however, it could be due to a specific splicing process: when compared to its closest sequences (clade 4, Figure 1 In conclusion, our results suggest that HLA-G phylogeny reflects HLA-G haplotype-specific association with different clinical conditions: clades 1, 2 and 3 gathered HLA-G sequences associated with immune impairments in diverse pathological conditions (Figure 1) (3, 10, 23). This study would benefit from the replication of pertinent results in additional cohorts to support the pertinence of defining HLA-G haplotype as a predictive genetic marker for asthma. Whereas the mechanisms with which HLA-G expression modulates asthma severity or phenotype are still poorly understood, it appears clearer that bronchial HLA-G Table 6 | Transcription factor (TF) binding site variation between human leukocyte antigen-G (HLA-G) haplotypes ("0" describes a loss of binding site; "1" a gain of binding site the sites variation; "−"a binding site diminution; and "+" a binding site augmentation).   expression is a more reliable marker than systemic and peripheral expression [(9) and our data]. However, given the discomfort of taking BAL samples and cross-reactivity in HLA-G antibody capture assays, HLA-G haplotype sequencing offers a reproducible qualitative marker that can be obtained from a sampling method which is easy to implement. Our results support a phylogenetic evolution of HLA-G sequences associated with differential qualitative and/or quantitative levels of protein expression, potentially driven by variations in TF binding sites and alternative splicing sites. Post-transcriptional mechanisms, not explored in this study, might also modify translation and shedding of sHLA-G and thus may further explain the lack of correlation observed between peripheral sHLA-G and genetic data in our asthmatic patient cohort. This descriptive study does not address the putative role of HLA-G in impaired control of bronchial epithelium immune response and/or pro-inflammatory and pro-angiogenic responses described at the maternal fetal interface (2,41,42), nor did we provide any experimental evidence for HLA-G haplotype actual involvement in HLA-G expression.