Genome-Wide Mapping of Plasma IgG N-Glycan Quantitative Trait Loci Identifies a Potentially Causal Association between IgG N-Glycans and Rheumatoid Arthritis

Rheumatoid arthritis (RA) is a systemic autoimmune disease characterized by chronic synovial joint inflammation, which can progress to a declined functional status, severe comorbidity, and shortened lifespan when left untreated or poorly controlled (1). RA is a complex multifactorial disorder, with contributions from both environmental triggers and genetic components (2–4). However, the exact etiology of RA is not fully understood, and there is pressing urgency to further explore the pathological mechanisms underlying RA to facilitate the design and implementation of efficient prevention strategies or novel treatments.

Glycosylation, which is an important and common posttranscriptional modification and is regulated by glycosyltransferases and glycosylase, has been shown to play a key role in many biological processes and diseases susceptibility (5, 6). Glycosylation of IgG influences its effector function by modulating binding to Fc receptors, acting as a switch between proinflammatory and anti-inflammatory IgG functionality (7–10). Serum IgG glycosylation was first reported to be implicated in RA (11). Subsequently, more and more studies have confirmed the role of IgG N-glycosylation in the development of RA (12–16). However, much like the observed associations in epidemiological studies, IgG N-glycosylation is prone to confounding factors and reverse causation. As a result, more studies are needed to explore causality between IgG N-glycosylation and RA.

Mendelian randomization (MR) is a popular approach to explore potentially causal associations between an exposure (e.g., IgG N-glycosylation) and an outcome (e.g., RA susceptibility) by using genetic variants as the instrumental variables (IVs) to reduce confounding and reverse causation. The increasing MR studies that integrate summarized genome-wide association study (GWAS) data for diseases and quantitative trait loci (QTL) data have been successful in identifying DNA methylation or gene expression or protein loci that are potentially causally associated with various phenotypes, such as cardiovascular diseases, inflammatory bowel disease, and arthritis (17–21).

However, the potentially causal relationship between IgG N-glycans and RA by linking IgG N-glycan-QTL variants with GWAS results has not been explored. Therefore, we adopted the MR method integrating summarized GWAS data for RA and IgG N-glycan-QTL data to prioritize IgG N-glycans that are potentially causally associated with the risk of RA in this study.

In our present study, we performed a standard two-sample MR analysis (Fig. 1A) to assess the association between exposure (each IgG N-glycan) and outcome (RA), in which the IVs-exposure (e.g., IVs–IgG N-glycan) and IVs-outcome (e.g., IVs-RA) associations were assessed in different samples. The previous studies reported that there were strong internal associations among IgG N-glycans, and that IgG N-glycosylations affected each other (8, 22, 23), which could be caused by the common regulatory mechanism (24, 25). As shown in (Fig. 1B, we combined all IgG N-glycans-QTL associated with all IgG N-glycans as IVs for each IgG N-glycan to further explore the association between exposure (each IgG N-glycan) and outcome (RA) by a multivariable MR analysis (26).

MR relies on some assumptions (27, 28). First, the genetic variant is robustly associated with the IgG N-glycosylation. Second, the genetic variant is independent of confounders of the IgG N-glycosylation-RA association. Third, the genetic variant is independent of the RA except via the IgG N-glycosylation. Fourth, all of the associations are linear and are unaffected by statistical interaction. For the first assumption, we filtered genetic variants with conservative multiple testing and combined all IgG N-glycans-QTL associated with all IgG N-glycans as IVs for each IgG N-glycan, alleviating the concern of weak IVs. For the second assumption, no association between genotype and confounders is assumed based on the common biological belief that the genotypes are not associated with socioeconomic and behavioral characteristics, which usually confound the effects of exposure on the outcome. With regard to the third assumption, we preformed several analyses to control for potential pleiotropic bias. We acknowledged that it is often difficult to validate the linearity assumption. However, the violation of the assumption was not essential when the aim was to test the null hypothesis of no effect of the exposure on the outcome.

In the MR analysis, IgG N-glycan-QTL genetic variants were used as the IVs for IgG N-glycosylation. Genome-wide analysis of IgG N-glycan was conducted in 536 blood donors of Chinese ancestry randomly selected from the Xuanwu Hospital in Beijing, China (29, 30). The plasma IgG N-glycosylation concentrations were quantified by hydrophilic interaction chromatography–ultra-performance liquid chromatography (31). Genotyping was conducted with Illumina OmniZhongHua chips (Illumina, San Diego, CA). Genotypes were imputed from the 1000 Genomes Project panel phase 3 based on the East Asian population using the Michigan Imputation Server (32). Genetic variants were excluded when they had call rate <99%, had minor allele frequency <0.01, and had an imputation quality ratio <0.7. Finally, 7,108,659 imputed single-nucleotide polymorphisms (SNPs) and 24 IgG N-glycans (glycan peaks [GPs]) were used for further IgG N-glycan-QTL mapping. The genetic associations were obtained in an additive genetic model, adjusted for age and sex. Based on not encountering the problem of population stratification, we did not correct the principal component. Conditionally uncorrelated variants (with the lowest p value having linkage disequilibrium r² <0.001) associated with IgG N-glycosylation (p < 1 × 10⁻⁵) were selected as IVs. We estimated the R² [2 × EAF × (1 − EAF) × β²] (33) of each SNP and summed them up to calculate the overall R² for each glycan. When the R² was greater than 1, we took the value of >0.999. Then, we calculated F statistics and the power of the study via https://shiny.cnsgenomics.com/mRnd/. A higher R², F statistic indicated a lower risk of weak IV bias.

The GWAS-summarized data for RA were reported in a Japanese population (34). The GWAS included a total of 212,453 individuals (4,199 cases and 208,254 controls) and 6,144,453 genetic variants. The GWAS incorporated age, sex, and the top five principal components as covariates. The GWAS summarized data can be downloaded at http://jenger.riken.jp/en/result/.

To explore the potentially causal associations of IgG N-glycans with RA, we used the inverse variance-weighted (IVW) method as the main analyses, in which the IVs-RA estimate is regressed on the IVs-IgG N-glycosylation estimate with the intercept term set to 0, weighted by the inverse variance of the IVs-RA estimate (35). The Cochran Q statistic based on the IVW method was conducted to test the heterogeneity of IVs. The random effects IVW method was used when heterogeneity existed; otherwise, the fixed-effects IVW method was used. We additionally conducted sensitivity analyses to account for potential violations of valid IV assumptions using the MR_Lasso, MR_Egger regression, and MR_weighted median. MR_Lasso extends the IVW model to include an intercept term for each genetic variant, in which these intercept terms represent associations of IVs-RA that bypass IgG N-glycosylation (36). The MR_Egger method may assess the robustness of estimates to potential violations of the standard IV assumptions because of directional pleiotropy, which was performed to valuate pleiotropy based on the intercept (37). The MR_weighted median method has greater robustness to individual genetic variants with strongly outlying causal estimates compared with the IVW and MR_Egger methods (38). The results of the IVW method were reported when pleiotropy bias (MR_Egger intercept p value >0.05) was also not observed. In addition, we performed a leave-one-out analysis in which we omitted one SNP in turn to investigate the influence of outlying or pleiotropic genetic variants on the result (35). Results are presented as β with SE or odds ratios (ORs) with their 95% confidence interval (CI) of outcomes per genetically predicted increase in each exposure factor. The data analyses were performed with the “MendelianRandomization” package.

Data cleaning and statistical analysis were performed using R version 4.0.2 (https://www.r-project.org/) and PLINK 1.9 (https://www.cog-genomics.org/plink/1.9/). A p value of <0.002 (0.05/24) was considered as suggestive of evidence for a potential association.

Written informed consent was obtained from each subject at the beginning of the study, and the study was approved by the Ethics Committee of the Capital Medical University (Beijing, China). The ethics approval was given in compliance with the Declaration of Helsinki. Participating studies of RA of GWAS summary statistic have received prior approval by relevant Institutional Review Boards, and informed consent was obtained from each study participant.

The average R² representing the variance of exposure explained by each SNP was 0.058, and the F statistics for each IgG N-glycan were >10, indicating that those included SNPs satisfied the strong relevance assumption and the weak IV bias would not likely jeopardize the causal inference (Supplemental Tables I, II). MR_Egger regression by the standard MR analysis showed no evidence of directional pleiotropy for the associations of IgG N-glycans with RA, except for GP18 (p = 0.012, Supplemental Table III). However, the threshold value of the pleiotropy test (p < 0.05) used was too conservative because a multiple test correction was not applied. No obvious heterogeneity was observed for any associations of IgG N-glycans with RA (all p values <0.05). Results of the IVW method demonstrated that there was one significant association of genetically instrumented IgG N-glycan (GP4) with RA after multiple testing (p < 0.002, Fig. 2). GP4 was negatively associated with RA, indicating that a higher GP4 was associated with the decreased risk of RA (OR_GP4 = 0.906, 95% CI = 0.857–0.958, p = 5.246 × 10⁻⁴; power = 0.97). As shown in (Fig. 3, results were broadly consistent in sensitivity analyses using MR_Lasso (OR_GP4 = 0.907, 95% CI = 0.862–0.955, p = 1.931 × 10⁻⁴), MR_weighted median (OR_GP4 = 0.901, 95% CI = 0.842–0.966, p = 3.095 × 10⁻³), and MR_Egger (OR_GP4 = 0.862, 95% CI = 0.754–0.985, p = 2.960 × 10⁻²; all p values <0.05). In addition, the results of leave-one-out sensitivity analyses showed that the associations of GP4 with RA were not substantially driven by any individual SNP (Fig. 4).

As shown in Supplemental Table IV, the observed associations between IgG N-glycans and RA by the multivariable MR analysis were limited evidence of horizontal pleiotropy based on the MR-Egger intercept test (all p values <0.05). Results of the IVW method demonstrated that there were seven significant associations of genetically instrumented IgG N-glycans (GP4, GP6, GP12, GP14, GP15, GP17, GP21) with RA after multiple testing (p < 0.002, Fig. 5). GP4 and GP6 were negatively associated with RA, whereas GP12, GP14, GP15, GP17, and GP21 were positively associated with RA. These observed associations could be confirmed in sensitivity analyses using MR_Lasso or MR_weighted median or MR_Egger (all p values <0.05).

Using two-sample MR analysis based on data from large-scale GWAS for RA and IgG N-glycan-QTL, our study demonstrated that a lower GP4 may lead to the risk of RA. In addition, the multivariable MR analysis showed that there were seven significant associations of genetically instrumented IgG N-glycans with RA. The findings were robust in sensitivity analyses with different statistical models and different instruments. To the best of our knowledge, this is the first study to examine a potentially causal association between IgG N-glycans and RA through an MR approach.

In our findings, we identified more associations between IgG N-glycans and RA in the multivariable MR analysis compared with those in the standard MR analysis, because the multivariable MR analysis can select more IgG N-glycan-QTL as IVs than the standard MR analysis to reduce the bias of weak IVs and improve the power of MR analysis. The associations of genetically instrumented IgG N-glycans (GP4, GP6, GP12, GP14, GP15, GP17, and GP21) with RA were identified by in our MR analysis, among which GP15, GP17, and GP21 were not reported by the observational studies (14, 16). GP4 is an agalactosylated and fucosylated glycan. The evidence showed that the decrease in galactosylation is generally inducing IgG Ab-driven inflammation, whereas IgG, deficient in the single fucose residue from the Fc glycan, can gain a 50-fold potency in terms of initiating Ab-dependent cellular cytotoxicity (8, 39). Higher GP4 was implicated as a risk factor for RA based on the European population (16). In addition, GP4 was positively associated with body mass index, and a higher GP4 was associated with the risk of obesity (39), which is a risk factor of RA (40). However, uncertainties remain over the nature of the association, perhaps complicated by confounders. To distinguish causation from correlation, we performed MR analysis to estimate a causal association between increased GP4 and RA risk independent of any confounding variables. We found that a genetically inherited higher GP4 level was associated with a lower risk of RA. This disparity suggested that the increased incidence of RA in higher GP4 is likely due to additional differences between the two groups. We hypothesize that higher GP4 or some environmental exposure associated with higher GP4, such as use of medications or a change in lifestyle, may reduce RA risk. More studies are needed to investigating what these differences might be.

A previous study in a large, well-powered cohort of RA patients showed that SNPs driving levels of N-glycosylation were not associated with RA susceptibility in the European population (41). Our finding, which was inconsistent with the previous study, might be explained by the fact that our study was based on an East Asian population, and thus may be due to race specificity and other differences. Second, we selected IgG N-glycan-QTL as the IVs for IgG N-glycan through genome-wide association, thereby increasing the power of MR analysis. A previous twin study found the heritability of RA to be ∼40–50% (42). Although GWASs have been successful in identifying genetic variants associated with RA (2, 3), biological interpretation of the findings remains largely undermined. A previous GWAS found that ∼88% of trait-associated genetic variants resided in noncoding regions of the genome and might perform regulatory functions on gene expression (43). The method incorporating QTL information into GWAS analyses has the potential to increase the power of GWASs in identifying loci associated with complex traits and improve the explanation of variances in traits (18, 44–46). Consistent with the previous studies that found that genetic variants play an important role in IgG N-glycosylation (24, 25, 47), our findings showed that IgG N-glycan-QTL variants and linking them to RA-associated genetic variants from GWASs might pinpoint molecular mechanisms underlying genetic susceptibility to RA that are due, at least in part, to altered IgG N-glycosylation. The IgG N-glycan-QTL resources provided by our study reveal a new richness of detail regarding genetic effects on IgG N-glycan patterns and the potentially causal relationship of IgG N-glycosylation with the RA phenotype. The above findings indicated that the therapeutic functions of Ab-based drugs should consider N-linked glycans on IgG and the related genetic variants.

However, there are several limitations in the current study. First, the MR analysis may be biased by potential violations of the IV assumptions. However, the pleiotropic effects were not observed in MR_Egger regression analysis, and sensitivity analyses with different statistical models and different instruments performed mostly similar results. It is difficult to completely exclude the potential influence of directional pleiotropy, which may lead to biased causal effect estimates. Second, we adopted correction for multiple testing to reduce the false-positive rate and retained independent IgG N-glycan-QTL to reduce linkage disequilibrium; however, we may have missed important genetic variants or IgG N-glycans. In addition, there are genetic differences in the genotypes for Chinese and Japanese people, which might lead to bias. The summarized GWAS analysis for RA did not control other confounding factors that might affect the outcome to some extent. We applied a relatively conservative Bonferroni correction to select IgG N-glycan-QTL as IVs for IgG N-glycan (i.e., p < 1 × 10⁻⁵), and we combined all IgG N-glycans-QTL associated with all IgG N-glycans as IVs for each IgG N-glycan to improve the power of MR analysis. The statistical powers ranged from 0.05 to 0.97 among GPs. We acknowledged this as a main limitation, considering the small sample size (n = 536) in the current study. Our findings were only based on participants of East Asian ancestry, which may not be generalizable to other populations. More studies with the same method are needed to validate our findings in a larger sample size and with independent populations.

In conclusion, our MR analysis incorporating GWAS and IgG N-glycan-QTL data revealed that IgG N-glycan showed a potentially causal association with RA. Our findings shed light on the role of IgG N-glycosylation in the development RA. Future studies are needed to validate our findings and to explore the underlying physiological mechanisms in the etiology of RA. The findings provide new insights into the mechanism underlying the association between genetic variants and RA. The therapeutic functions of Ab-based drugs should consider N-linked glycans on IgG and the related genetic variants.

We thank all of the research subjects for their participation and the Riken group, who reported the summary data, and we acknowledge the skillful work of the entire medical staff at Xuanwu Hospital.

The authors have no financial conflicts of interest.