Mapping of the gene network that regulates glycan clock of ageing

Results

To enable meta-analysis of GWAS of galactosylation phenotypes measured by different analytical platforms and thereby increase the total number of samples, we first developed a protocol for data harmonisation for glycan data generated using different analytical platforms (UPLC and LC-MS). The protocol was established based on the estimated correlation of the glycan traits derived in one cohort that had glycans measured using both UPLC and LC-MS platforms. The aim was to combine IgG subclass-specific glycan information obtained from LC-MS (i.e. galactosylation measured separately on IgG1-4 subclasses) in an appropriate manner to obtain information corresponding to the total IgG galactosylation measured by UPLC (analysis of glycans from both Fab and Fc fragments from all IgG subclasses). As described in Supplementary Methods, we compared different normalization and subclass information summary approaches. In a subset of individuals that were analysed using both methods, the Pearson correlation coefficients ranged from 0.95 and 0.97 for G0, 0.73 to 0.80 for G1, and 0.90 and 0.91 for G2 trait (Appendix Table 1). Different pre-processing approaches performed similarly, so we selected median quotient normalization due to previous recommendation [36], followed by the calculation of the derived traits without IgG subclass weighting in LC-MS data.

Discovery and replication GWAS

We performed a discovery GWAS of IgG galactosylation in seven cohorts of European descent (N=13,705). The association between HRC-imputed SNPs and three phenotypes (G0, G1 and G2) was studied under the assumption of an additive linear model. The inverse-variance meta-analysis of GWAS summary statistics resulted in 16 genome-wide significant (p ≤ 2.5 × 10^-8) genomic loci associated with at least one galactosylation trait (Table 1). Quantile-quantile plots of associations for three traits are shown in Supplementary Figure 1.

Table 1. List of genome-wide significant loci in IgG galactosylation GWAS.

Locus

Top SNP

Trait

EA

OA

EAF

Beta

SE

P

Beta repl

SE repl

P repl

Closest gene

Prioritized genes

1:24699711-25493756

rs188468174

G0

T

C

0.015

0.693

0.048

3.20E-47

0.865

0.076

2.96E-30

RUNX3

RUNX3

2:26109539-26149988

rs10177977

G0

T

C

0.334

-0.062

0.01

3.80E-10

-0.049

0.015

1.03E-03

KIF3C

KIF3C

4:103390496-103567348

rs3774964

G2

A

G

0.37

-0.065

0.01

1.56E-11

-0.038

0.015

9.20E-03

NFKB1

NFKB1; MANBA

6:31107733-31164511

rs1265109

G0

T

G

0.295

-0.063

0.011

1.57E-09

0.013

0.019

4.97E-01

HLA region

HLA region

6:74168723-74285118

rs3822960

G2

T

C

0.306

0.061

0.01

5.84E-10

0.028

0.015

6.24E-02

EEF1A1

EEF1A1; MTO1

6:143088071-143206826

rs7758383

G2

A

G

0.486

0.092

0.009

6.09E-23

0.145

0.014

9.44E-25

HIVEP2

HIVEP2

7:150856165-150906453

rs113745074

G0

T

C

0.094

0.12

0.015

6.19E-16

0.134

0.022

8.47E-10

SMARCD3

ABCF2; CHPF2; SMARCD3

8:103542538-103550211

rs13250010

G0

T

G

0.356

0.063

0.01

7.71E-11

0.062

0.015

3.47E-05

KB-1980E6.3

UBR5; RRM2B; ODF1; KB-1980E6.3

9:32933492-33385427

rs13297246

G2

A

G

0.184

0.195

0.013

4.21E-55

0.202

0.019

2.36E-26

B4GALT1

B4GALT1

11:65555524-65555524

rs10896045

G1

A

G

0.301

0.079

0.013

2.61E-09

0.079

0.017

2.42E-06

OVOL1

OVOL1; AP5B1

17:16813994-16875636

rs34562254

G1

A

G

0.106

0.166

0.019

1.48E-18

0.092

0.024

9.95E-05

TNFRSF13B

TNFRSF13B

17:43856639-44863413

rs199516

G0

T

C

0.221

0.068

0.012

1.20E-08

0.058

0.018

1.25E-03

WNT3

ARHGAP27; CRHR1; SPPL2C; MAPT; KANLS1; ARL17B; LRRC37A; NSF; WNT3

17:45766846-45870129

rs1808192

G2

A

G

0.351

0.061

0.01

5.94E-10

0.051

0.015

6.78E-04

TBKBP1

TBKBP1; TBX21

17:56404349-56418136

rs2526377

G2

A

G

0.463

0.062

0.009

5.59E-11

0.057

0.015

8.05E-05

BZRAP1

BZRAP1; SUPT4H1; RAD5C1

17:79158040-79268562

rs2659005

G2

T

C

0.436

0.079

0.009

4.36E-17

0.097

0.016

6.99E-10

SLC38A10

AZI1; ENTHD2; SLC38A10; C17orf89

21:36564553-36665202

rs4817708

G0

T

C

0.229

0.073

0.011

2.39E-11

0.043

0.017

1.15E-02

RUNX1

RUNX1

Locus- chromosome:start-end position in GRCh37 (hg19) build; Top SNP- rsID identifier for the SNP with the strongest association with the glycan trait; Trait- trait with the lowest p-value for the top SNP in the genomic locus; EA- effect allele for which the effect is reported; OA- non-effect allele; EAF- effect allele frequency; Beta- effect estimate for the effect allele of the top SNP; SE- standard error of the effect estimate; P– p-value for the top SNP-trait association; Beta repl- effect estimate for the effect allele of the top SNP in replication analysis; SE repl- standard error of the effect estimate in the replication analysis; P repl- p-value for the top SNP in the replication analysis; Closest gene- gene found closest to the top SNP in the locus; Prioritized genes- symbol for genes prioritized in the locus.

For twelve loci, the same trait-SNP association reached the significance threshold p-value ≤ 0.0031 (P ≤ 0.05/16 loci) in replication meta-analysis (N=7,775). The association between the top SNP in the HLA region (chr6:31107733-31164511) was not replicated due to the unavailability of the SNP and SNPs in LD in all four replication cohorts. However, the association of the HLA region and IgG glycan patterns was observed in the previous GWAS of IgG glycome [30, 31, 33, 34]. For top SNPs in three regions (chr6:74168723-74285118, chr4:103390496-103567348 and chr21:36564553-36665202) the association was not significant after correction for multiple testing in the replication meta-analysis, but the direction of the effect was consistent with the discovery analysis (Table 1).

Due to lack of compatibility in glycan trait definition in the current study and previous GWAS of IgG glycome [30–34], we checked for the overlap of the associated genomic regions across any previously examined IgG glycan traits. Thirteen loci overlap with previously identified and replicated genomic regions (Appendix Table 2). The remaining three associations in chr4:103390496-103567348, chr17:56404349-56418136 and chr6:74168723-74285118 are considered as novel.

Among 16 loci, three loci were G2-specific (rs3822960, rs1808192 and rs2526377), two were G1-specific (rs10896045 and rs34562254) and two were G0-specific (rs13250010 and rs199516) (Figure 1A), which is intriguing since conversion of G0 to G1, and then G2 is a consequence of the activity of the single enzyme (product of B4GALT1 gene). Before functional evaluation of GWAS hits, careful inspection of the genomic regions was undertaken to prioritize and choose the genes with strong evidence for a potential role in a biological pathway. We applied a set of in silico methods to prioritise genes for functional analysis. First, we mapped associated variants to genes using mapping strategies based on position, overlap with eQTL associations and chromatin interaction as implemented in FUMA [37], and identified 107 candidate genes in total (Appendix Table 3). The number of genes prioritised by two or more approaches is shown in Figure 1B.

Figure 1. Genome-wide significant associations with galactosylation phenotypes. (A) Top SNP in identified genomic regions for each associated trait, (B) Venn diagram showing the number of genes mapped by positional mapping, chromatin interaction mapping, eQTL mapping and genome-wide gene-based association analysis (MAGMA), (C) Manhattan plot of genome-wide significant associations in IgG galactosylation GWAS with prioritized genes in each locus. Plot shows -log10(p-values) of association on y-axis and SNPs ordered by chromosomal location on x-axis. Red line indicates the genome-wide significance threshold (2.5 ×10^-8). Orange gene names indicate novel loci associated with IgG glycosylation.

The gene analysis in MAGMA is based on a multiple linear principal components regression model in which the statistic is computed by combining the p-values of the individual SNPs in a gene while taking into account their LD structure. In gene analysis, SNP data is aggregated to the whole gene level to test the joint association of all SNPs in the gene with the phenotype, thereby making it possible to detect effects consisting of multiple weaker SNP-phenotype associations that would otherwise be missed. The gene-based association test for G0, G1 and G2 identified 38 genes significant at FDR < 0.05 (Appendix Table 4).

Given that missense variants have a direct impact on protein function and therefore a more straightforward effect on the phenotype, we also assessed the functional consequences of associated variants using Variant Effect Predictor [38]. For two genes, SLC38A10 and MAPT, variants with potentially pathogenic amino acid (aa) change were detected by both SIFT and Polyphen2 algorithms. For additional four genes, MANBA, SPPL2C, CRHR1 and TNFRSF13B, the variants were predicted to result in a benign aa change. The details on position and aa changes are listed in Supplementary Table 1.

Finally, to assess the potential pleiotropic effects of the variants on IgG galactosylation and gene expression we performed colocalization analysis with the whole blood eQTLgen dataset [39]. By comparing the regional association patterns with the approximate Bayesian Factor approach in coloc [40] we identified five genes (KIF3C, NFKB1, MIR-142, EEF1A1 and MTO1) whose expression patterns colocalize with IgG galactosylation (PP4 > 75%) (Supplementary Table 2 and Supplementary Figure 2). Previous prioritization efforts in GWAS of IgG glycosylation, as well as implications for gene's potential function in glycosylation process were taken into account in the final gene prioritization, which resulted in 37 credible gene candidates in the identified regions (Figure 1C). Supplementary Table 1 contains details on prioritization evidence for each gene.

IgG galactosylation and diseases

Due to existing evidence of altered IgG glycosylation patterns in multiple conditions [11], we assessed the pleiotropy of the IgG galactosylation-associated loci with a range of autoimmune, inflammatory and infectious diseases. The threshold of 75% and higher was applied for the posterior probability for colocalization of traits (PP4). A shared association pattern between SLE and monogalactosylation (PP4 = 82%) was detected in a locus on chromosome 11 (chr11:65555524), while agalactosylation-associated locus chr17:43856639-44863413 colocalized with breast cancer (PP4 = 95%), COVID-19 hospitalization (PP4 = 93%), SLE (PP4 = 89%) and schizophrenia (PP4 = 84%) (Supplementary Table 3). Because of the extensive LD pattern, there are more than twelve candidate genes in the chromosome 17 locus. It is important to note that high PP4 indicates at least one shared causal variant and low PP4 does not mean that two phenotypes do not share one, especially if the probability of both phenotypes being associated with the region but having different causal variants (PP3) is high.

Polygenic score

A polygenic score (PGS) aggregates the effects of many genetic variants into a single number which predicts genetic predisposition for the phenotype. In the standard approach, a PGS is a linear combination of linear regression effect estimates and allele counts at single-nucleotide polymorphisms (SNPs). For each of the G0, G1 and G2 traits, we created two PGS models based on genome-wide association meta-analysis (GWAMA) results. Specifically, we derived the first model based on the SBayesR [41] method, and the second model using a clumping and thresholding (C+T) method. These models were tested within the CEDAR dataset, which has aggregated genotype data and galactosylation information on 162 participants of European ancestry and which was not used in the GWAMA. The estimated SNP-based heritability (h2SNP, SBayesR package [41]) was 36.5% for the G0, 27.2% for G1, and 36.8% for the G2 trait. Our best PGS models explained 7.4% (p-value: 4.79 x 10^-4), 2.2% (p-value: 0.059) and 7.4% (p-value: 4.77 x 10^-4) of G0, G1 and G2 variance in CEDAR samples. PRS-Phenotype scatter plots and plots showing the mean of trait per group can be found in Supplementary Figures 5–7. The complete results of the implementation of SBayesR and C+T models can be found in Appendix Table 5.

Functional validation of GWAS hits

The genes for functional follow-up were selected based on novelty in GWA studies for IgG galactosylation and by reviewing the evidence for prioritization, expression in B cells (Human Protein Atlas), as well as known roles of the encoded proteins. Newly discovered GWAS hits, EEF1A1 and MANBA-NFKB1, where NFKB1 and MANBA are found in the same locus, were selected and functional follow-up was performed to decipher which one is indeed involved in galactosylation. We also included HIVEP2 because it was previously assigned in the gene network as a transcription factor potentially involved in the regulation of B4GALT1 gene expression [33], and TNFRSF13B for its highly specific expression in B cells and function in humoral immunity. SLC38A10 was included as one of the genes prioritized in chr17:79158040-79268562 locus to test its potential function as a transporter protein in the IgG glycosylation pathway. We also included KIF3C as it was the prioritized gene in chr2:26109539-26149988 locus.

To functionally validate the potential role of the selected GWAS hits in IgG galactosylation, we used the newly developed transient expression system HEK293FS for IgG production with stably integrated dSaCas9-VPR or dSpCas9-KRAB expression cassettes for target gene upregulation or downregulation [35]. Therefore, the candidate genes were manipulated in dCas9-VPR or dCas9-KRAB containing monoclonal cell lines transiently transfected with a plasmid carrying genes for IgG heavy and light chains and specific gRNA. Subsequently, secreted IgG was analysed for glycan phenotype.

Utilizing dCas9-VPR cell line and three gRNAs for each locus, we successfully upregulated all selected loci except SLC38A10 and HIVEP2 (Supplementary Figure 3A and Supplementary Table 4). Successful downregulation in the dCas9-KRAB cell line was accomplished using three gRNAs for all selected loci except NFKB1 and HIVEP2 (Supplementary Figure 3B). A significant change of IgG glycome composition was observed after dCas9-VPR targeting of HIVEP2, MANBA, TNFRSF13B and EEF1A1 (Figure 2 and Appendix Table 6). Interestingly, manipulation of the HIVEP2 promotor did not upregulate HIVEP2 transcript level sufficiently to reach statistical significance at transcript level, but it resulted in a significant increase of digalactosylated structures and a concomitant decrease of agalactosylated structures on IgG (Figure 2A). Upregulation of the MANBA locus resulted in a significant decrease of monogalactosylated IgG glycans, but no change in levels of agalactosylated and digalactosylated glycans was observed (Figure 2B). Targeting TNFRSF13B in dCas9-VPR cell line resulted in almost 400-fold increase of its expression and a moderate increase in agalactosylated IgG glycans (Figure 2C). Upregulation of the EEF1A1 locus resulted in a significant decrease of monogalactosylated glycans and in a concomitant increase of agalactosylated IgG glycans (Figure 2D).

Figure 2. Changes in IgG glycan composition upon targeting of selected GWAS loci associated with IgG galactosylation (HIVEP2, MANBA, TNFRSF13B and EEF1A1) in dCas9-VPR monoclonal cell lines. Samples containing non-targeting gRNAs served as controls. Changes in transcript levels are given as fold change values and changes in IgG phenotype are given as a relative change compared to control samples. (A) Manipulation of the HIVEP2 gene did not result in a statistically significant change in HIVEP2 transcript level, however, did induce a significant increase of digalactosylated structures (G2) with a concomitant decrease of agalactosylated IgG glycan structures (G0). (B) Targeting of MANBA by dCas9-VPR elevated transcription level of this gene which resulted in decrease of monogalactosylated IgG glycan structures (G1) (C) Targeting TNFRSF13B by dCas9-VPR resulted in ~ 400-fold increase of transcript levels and significant change of IgG agalactosylated IgG glycans (G0). (D) Successful upregulation of the EEF1A1 locus was followed by an increase of agalactosylated IgG glycans (G0) with a concomitant decrease of monogalactosylated IgG glycans (G1). Nominal p-value: *<0.05; **<0.01; ***<0.001; ns, not significant.

Discussion

Using GWAS approach we identified 16 genomic regions associated with IgG galactosylation, which is the molecular basis of the glycan clock of ageing. Thirteen out of 16 loci were associated with traits related to IgG glycosylation in our previous GWAS studies, while three were novel [30–34]. Contrary to previous studies, in which we were looking for genes that were associated with any aspect of IgG glycosylation, in this study we focused on galactosylation, a functionally well-defined element of IgG glycosylation that is the basis of increased inflammation in people with accelerated glycan ageing. This was done by converting percentages of individual glycans in the glycome into sums of glycan structures containing either one, two or no galactose residues. This enabled us to combine GWA studies for which the phenotype analysis was performed using different analytical methods. We also estimated the SNP-based heritability to be 36.5%, 27.2% and 36.8%, for G0, G1 and G2, respectively, an estimate that very closely replicates recent heritability estimate for biological age estimate based on these glycan structures [42]. Upon prioritization efforts which resulted in 37 candidate genes, we evaluated functional relevance of these genes using HEK293FS transient expression system for IgG production with stably integrated dCas9-KRAB or dCas9-VPR fusion. Results of this study expanded the list of genes for which the functional relevance for in vitro expression and glycosylation of IgG was confirmed with MANBA, TNFRSF13B and EEF1A1 (Figure 2). Thus, out of 13 replicated GWAS hits, for six of them we have functionally confirmed that they act at IgG glycosylation in a way that is conserved between plasma cells and our artificial in vitro expression system in HEK293FS cells (Figure 3). It is likely that some additional GWAS hits may be plasma cell specific and will be confirmed once an efficient system for gene manipulation and IgG production in B cells is developed. However, it is important to note that IgG glycosylation can be affected at multiple other levels (e.g. B-cell selection, maturation, expansion, etc.), thus the absence of functional confirmation of other loci in this assay does not mean that they are not implicated in IgG glycosylation in some other, more complex, manner.

Figure 3. A graphical review of up-to-date functional validation of GWAS hits associated with IgG glycosylation using HEK293FS transient expression system with stably integrated dCas9-VPR or dCas9-KRAB fusions. When specific gRNA targeted dCas9-VPR to two estrogen associated genes, RUNX3 and SPINK4 (Mijakovac et al. 2021), as well as to the TNFRSF13B, MANBA and EEF1A1 loci, transcription was upregulated that resulted in decreased levels of galactosylated (red arrow within box), increased levels of agalactosylated IgG glycan structures (green arrow within box) or both. Downregulated transcription of SPPL3 resulted in an increase in galactosylated structures with a concomitant decrease in agalactosylated IgG glycan structures. Protein structures do not depict true protein structures in humans, generic protein shapes are chosen for easier visualization. Figure was created with BioRender.com. Accessed on 16 November 2021.

The total heritability of G0, G1, G2 is not known, however, the heritability of individual IgG glycan peaks ranges between 27% and 76% [29]. More than half of the studied glycan traits have additive genetic component estimate of 50% which would imply that half of the variance in their values could be explained by common variants, an observation similar to other complex traits [29]. The constructed PGS models explained 7.4% of the variance in G0 and G2 traits in an out-of-sample testing, thus capturing about 20% of the estimated SNP heritability. However, for G1 the PGS model explained a minor proportion of SNP heritability. The explanation why for G1, that had comparable SNP heritability, the proportion of variance explained was much less may lie in a flatter distribution and hence individually smaller genetic effects onto G1; and observation which may be described as “higher polygenicity”. Indeed, while for G0 and G2 lead associations explained about 1.5% of variance, for G1 the lead variant near TNFRSF13B explained only about 0.5% of its variation. With smaller individual effects, one needs larger discovery GWAS to construct powerful PGS.

Given that the majority of identified candidate SNPs are located in non-coding regions, we applied several in silico approaches to prioritize potential causal genes in associated genomic regions, including exploring functional consequences of associated variants, pleiotropy with gene expression, gene-based association analysis, as well as evidence from previous GWA studies and prior biological knowledge of gene function. Among 37 candidate genes from 16 associated loci, for the two genes, RUNX3 and B4GALT1, the effect on IgG glycosylation has recently been demonstrated using direct gene expression modulation in HEK293FS IgG transient expression system [35, 43].

Due to existing evidence of altered IgG glycosylation patterns in multiple conditions [11], we performed colocalization analysis with a range of diseases and observed a high probability for pleiotropic effects in the AP5B1/OVOL1 locus on chromosome 11 for SLE and monogalactosylation. Interestingly, for the same locus, we observed considerable probability for colocalization (PP4=0.56) with asthma and allergy. Although this is below our pre-defined posterior probability for colocalization (PP4) threshold of PP4 > 0.75, the finding is in line with previous results linking this locus with variation of specific digalactosylated structures (FA2FG2S1), expression of the OVOL1 gene and the risk of asthma [34]. OVOL1 gene encodes a putative zinc finger containing transcription factor with a role in the development and differentiation of epithelial and germ cells [44]. The locus on chromosome 17 was pleiotropic for agalactosylation and SLE, breast cancer, COVID-19 hospitalization and schizophrenia. The colocalized region on chromosome 17 (chr17:43856639-44863413) contains at least twelve genes and is known for numerous copy number variants. It also includes a megabase-long inversion polymorphisms (H1 and inverted H2 forms) which contain partial duplication of KANSL1 gene with inverted H2 isoform having a higher frequency than H1 isoform in the European population [45]. The complexity of the locus makes it challenging to prioritise a specific gene without further investigation.

The functional follow-up was based on the previously developed HEK293FS transient expression system for IgG production which was shown to be an excellent tool for studying the role of the genes identified in glycosylation GWAS [35]. Despite the fact that we could not demonstrate upregulation of HIVEP2 mRNA using dCas9-VPR, a significant increase in digalactosylated structures with a concomitant decrease in agalactosylated structures was measured from IgG secreted from the same cells. HIVEP2 (HIVEP zinc protein 2) encodes a transcription factor that interacts with numerous viral and cellular promoters, and has an established role in immunity [46, 47]. Even though we prioritised HIVEP2, the locus interacts with variants in six genes via chromatin interactions (Supplementary Figure 4), including FUCA2 (alpha-L-fucosidase) which is located 0.6 megabases away. It is plausible that dCas9-VPR targeting of the HIVEP2 promoter disrupted these chromatin interactions which in turn altered the expression of other genes mapping to this region. In the case of the locus on chromosome 4, the functional study was utilized as one of the steps to decipher which of the two genes might influence galactosylation process as both NFKB1 and MANBA were prioritised with in silico approaches. The MANBA/NFKB1 region was identified in previous studies as a disease susceptibility locus for primary biliary cholangitis (PBC) [48] and a recent study found a decreased level of IgG galactosylation in PBC patients [49]. The NFKB1 gene encodes a subunit of the nuclear factor of kappa light polypeptide gene enhancer in B-cells (NF-κB) transcription factor family, known for its critical role in cell survival and inflammation [50]. The upregulation of NFKB1 using dCas9-VPR resulted in a significant change in transcript level, but no change in IgG glycoprofile was observed. On the other hand, the upregulation of the MANBA gene expression resulted in a decrease of IgG monogalactosylated species. The MANBA gene encodes beta mannosidase which is an exoglycosidase enzyme cleaving beta-mannose residues from the non-reducing end of N-linked glycans [51]. As such, it has a known role in glycosylation process but currently we cannot speculate about the exact mechanism by which it affects IgG galactosylation. Interestingly, mutations in MANBA are known to affect kidney function [52], blood pressure [53] and cardiovascular disease [54], while the monogalactosylated IgG glycan (G1) (that MANBA affects in our in vitro expression system) was previously identified as the best predictor of future cardio vascular events in women [55, 56]. While at the moment we do not know putative mechanisms that could explain the link between MANBA, G0 and cardiovascular or renal diseases, these intriguing associations warrant future studies.

The EEF1A1 (Eukaryotic Translation Elongation Factor 1 Alpha 1) gene encodes for an isomer of an alpha subunit of the eukaryotic translation elongation factor-1 complex that has a role in binding of aminoacyl tRNA to the A site on ribosomes during protein synthesis [57]. Together with EEF1A2, it is the second most abundant protein (1%-3%) in the cell with other roles outside of protein synthesis, including protein degradation, apoptosis modulation, oncogenesis and viral pathogenesis [58, 59]. Previous study has shown that EEF1A1 also serves as a signal transducer during inflammation, where it can enhance interleukin 6 expression through STAT3 and PKCδ [60]. Interleukin 6 is implicated in different autoimmune and inflammatory diseases [61] and its high levels lead to decreased IgG galactosylation [62]. By targeting dCas9-VPR to EEF1A1, we increased its expression level, which resulted in a decrease of galactosylated IgG glycans. Therefore, our results are concordant with the previous findings, suggesting that EEF1A1 upregulation could lead to pro-inflammatory modulation of IgG by decreasing its galactosylation.

A locus harbouring TNF receptor superfamily member 13B (TNFRSF13B) gene was previously identified in the multivariate IgG glycosylation GWAS [34], while here it was associated with IgG monogalactosylation. TNFRSF13B encodes a transmembrane activator calcium modulator and cyclophilin ligand interactor (TACI) protein, a lymphocyte-specific member of the tumour necrosis factor (TNF) receptor superfamily, which has a role in signalling pathway leading to B cell differentiation and antibody production [63, 64]. Common genetic variation in TNFRSF13B locus was previously associated with different levels of immunoglobulins in serum [65], while rare deleterious variants were implicated in common variable immunodeficiencies and IgA deficiencies [66, 67]. Therefore, we speculate that TNFRSF13B might impact glycosylation profile of the sample by controlling secretion of certain IgG glycoforms. Despite extensive increase in TNFRSF13B expression upon activation with dCas9-VPR, we observed only a moderate increase in IgG agalactosylation and no changes in mono- or digalactosylation. Our experimental setting might be limited by the B-cell specific function of TNFRSF13B and thus, we might not be able to observe changes comparable to changes in the gene expression level.

The role of four other prioritized genes, NFKB1, KIF3C, HIVEP2 and SLC38A10, could not be validated in our HEK293FS transient expression system. Even though this cell system represents a good model for studies of IgG glycosylation regulation, it has some limitations. HEK293FS cells do not inherently secrete immunoglobulins and their epigenetic context is not equivalent to plasma and B cells. In addition, alternative galactosylation might depend on gene expression in other cell types and their interactions, hence we might not observe the natural effects after manipulation of gene expression in this model system, due to the lack of the right biological context.

In conclusion, in the present study we mapped 16 genetic loci regulating glycan structures that are not only the basis of the glycan clock of ageing, but also functional effectors of inflammation at multiple levels. Furthermore, we were able to confirm the functional role of three genes (MANBA, TNFRSF13B and EEF1A1) in the IgG galactosylation pathway by use of the HEK293FS transient expression system with stably integrated dCas9 fusion designed for this purpose. Figure 3 summarizes all GWAS hits with previously unknown link to IgG glycosylation, whose functional relevance has been shown using this cell system. Further research is needed to fully elucidate functional mechanism behind their role in ageing and to reveal the complete network of gene interactions regulating the complex process of IgG glycosylation. However, this study is an important step in this direction since it clearly indicates that regulation of IgG galactosylation extends far beyond the expression of the galactosyltransferase gene that adds galactose.

Author Contributions

AFH- data analysis and interpretation, visualisation, writing- original draft preparation, writing- review and editing. KM- Functional study experiments, data analysis and interpretation, visualisation, writing- original draft preparation, writing- review and editing. AM- functional study experiments, data analysis and interpretation, visualisation, writing- original draft preparation, writing- review and editing. AN- PGS calculation, visualisation, writing- original draft preparation, writing- review and editing. SZS- PGS calculation, data analysis. AL- data analysis, TH- performing GWAS in EGCUT cohort; EA- performing GWAS in LLS cohort. SS- performing GWAS in KORA F4 cohort. RRCC- performing GWAS in EPIC cohort. MM- performing GWAS in EPIC cohort. YL- performing GWAS in GCKD cohort. TK- IgG glycan quantification. NR- IgG glycan quantification. TŠ- IgG glycan quantification, writing- review and editing. MPB- IgG glycan quantification, writing- review and editing. ITA- IgG glycan quantification, writing- review and editing. IG-IgG glycan quantification. JŠ- IgG glycan quantification. TP- IgG glycan quantification. BR- IgG glycan quantification. PT- IgG glycan quantification. KF- genomic and demographic data provider for EGCUT cohort. MB- genomic and demographic data provider for LLS cohort. MW- genomic and demographic data provider for LLS cohort, writing- review and editing. CG- genomic and demographic data provider for KORA F4 cohort. MBS- genomic and demographic data provider for EPIC cohort, writing- review and editing. CW- Genomic and demographic data provider for EPIC cohort, writing- review and editing. OP- genomic and demographic data provider for CROATIA cohorts. CH- genomic and demographic data provider for ORCADES cohort. JFW- genomic and demographic data provider for VIKING cohort, writing- review and editing. TDS- genomic and demographic data provider for TwinsUK cohort. AK- genomic and demographic data provider for GCKD cohort, writing- review and editing. FV- IgG glycan data quality control, supervision, YSA- PGS calculation supervision, writing- original draft preparation, writing- review and editing. AV- data analysis, visualization, supervision of functional study. JK- IgG glycan quantification in functional study, IgG glycan quantification in cohorts, supervision, writing- original draft preparation, writing- review and editing. LK- Conceptualisation, supervision, data interpretation, writing—original draft preparation, writing- review and editing, VZ- supervision for the functional study, data interpretation, writing- original draft preparation, writing- review and editing., GL- conceptualisation, supervision, data interpretation, glycan data provider, writing-original draft preparation, writing- review and editing.

Conflicts of Interest

YSA is a cofounder and a co-owner of PolyOmica and PolyKnomics, private organizations providing research services in the field of quantitative, computational and statistical genomics. GL is the founder and CEO of Genos Ltd, a private research organization that specializes in high throughput glycomics analysis and has several patents in this field. AFH, TŠ, MPB, ITA, IG, JŠ, TP, BR, PT, FV and JK are employees of Genos Ltd. The remaining authors declare no Conflicts of interest.

Ethical Statement and Consent

The participating cohorts were approved by local research ethics committees for academic and commercial use of the study. The informed consent was obtained from the participants of the study.

Funding

This work was supported by European Structural and Investment Funds grant for the Croatian National Centre of Competence in Molecular Diagnostics grant KK.01.2.2.03.0006, Croatian National Centre of Research Excellence in Personalized Healthcare grant KK.01.1.1.01.0010, IRI “CardioMetabolic” grant KK.01.2.1.02.0321. The work was co-funded by the European Union (ERC, GlycanSwitch, 101071386). AFH was supported by H2020-MSCA-ITN IMforFUTURE grant 721815. The work of AN was supported by the Ministry of Education and Science of the Russian Federation via the state assignment of the Novosibirsk State University (project “Graduates 2020”). The work of SZS and YSA was partially supported by the Research Program at the MSU Institute for Artificial Intelligence. The work of LK was supported by an RCUK Innovation Fellowship from the National Productivity Investment Fund (MR/R026408/1). The work of YL was supported by the German Research Foundation (grant KO_3598/4-2 to AK). The work of AK was supported by the German Research Foundation (grant KO_3598/5-1). RRCC, CW and MBS were supported by German Ministry of Education and Research (BMBF) and the State of Brandenburg DZD grants 82DZD00302 and 82DZD03D03. CW was also supported by SciLifeLab and Wallenberg Data Driven Life Science Program grant KAW 2020.0239. TwinsUK study was funded by Wellcome Trust grant 212904/Z/18/Z, Medical Research Council AIMHY grant MR/M016560/1, European Union H2020 grant 733100. TwinsUK and MM were also supported National Institute for Health Research (NIHR)-funded BioResource, Clinical Research Facility and Biomedical Research Centre based at Guy’s and St Thomas’ NHS Foundation Trust in partnership with King’s College London. Views and opinions expressed are, however, those of the author(s) only and do not necessarily reflect those of the European Union or the European Research Council. Neither the European Union nor the granting authority can be held responsible for them.

References

• 1. Moremen KW, Tiemeyer M, Nairn AV. Vertebrate protein glycosylation: diversity, synthesis and function. Nat Rev Mol Cell Biol. 2012; 13:448–62. https://doi.org/10.1038/nrm3383 [PubMed]
• 2. Wada R, Matsui M, Kawasaki N. Influence of N-glycosylation on effector functions and thermal stability of glycoengineered IgG1 monoclonal antibody with homogeneous glycoforms. MAbs. 2019; 11:350–72. https://doi.org/10.1080/19420862.2018.1551044 [PubMed]
• 3. Zheng K, Bantog C, Bayer R. The impact of glycosylation on monoclonal antibody conformation and stability. MAbs. 2011; 3:568–76. https://doi.org/10.4161/mabs.3.6.17922 [PubMed]
• 4. Masuda K, Kubota T, Kaneko E, Iida S, Wakitani M, Kobayashi-Natsume Y, Kubota A, Shitara K, Nakamura K. Enhanced binding affinity for FcgammaRIIIa of fucose-negative antibody is sufficient to induce maximal antibody-dependent cellular cytotoxicity. Mol Immunol. 2007; 44:3122–31. https://doi.org/10.1016/j.molimm.2007.02.005 [PubMed]
• 5. Shields RL, Lai J, Keck R, O’Connell LY, Hong K, Meng YG, Weikert SHA, Presta LG. Lack of fucose on human IgG1 N-linked oligosaccharide improves binding to human Fcgamma RIII and antibody-dependent cellular toxicity. J Biol Chem. 2002; 277:26733–40. https://doi.org/10.1074/jbc.M202069200 [PubMed]
• 6. Dekkers G, Treffers L, Plomp R, Bentlage AEH, de Boer M, Koeleman CAM, Lissenberg-Thunnissen SN, Visser R, Brouwer M, Mok JY, Matlung H, van den Berg TK, van Esch WJE, et al. Decoding the Human Immunoglobulin G-Glycan Repertoire Reveals a Spectrum of Fc-Receptor- and Complement-Mediated-Effector Activities. Front Immunol. 2017; 8:877. https://doi.org/10.3389/fimmu.2017.00877 [PubMed]
• 7. Dekkers G, Rispens T, Vidarsson G. Novel Concepts of Altered Immunoglobulin G Galactosylation in Autoimmune Diseases. Front Immunol. 2018; 9:553. https://doi.org/10.3389/fimmu.2018.00553 [PubMed]
• 8. van Osch TLJ, Nouta J, Derksen NIL, van Mierlo G, van der Schoot CE, Wuhrer M, Rispens T, Vidarsson G. Fc Galactosylation Promotes Hexamerization of Human IgG1, Leading to Enhanced Classical Complement Activation. J Immunol. 2021; 207:1545–54. https://doi.org/10.4049/jimmunol.2100399 [PubMed]
• 9. Peschke B, Keller CW, Weber P, Quast I, Lünemann JD. Fc-Galactosylation of Human Immunoglobulin Gamma Isotypes Improves C1q Binding and Enhances Complement-Dependent Cytotoxicity. Front Immunol. 2017; 8:646. https://doi.org/10.3389/fimmu.2017.00646 [PubMed]
• 10. Wei B, Gao X, Cadang L, Izadi S, Liu P, Zhang HM, Hecht E, Shim J, Magill G, Pabon JR, Dai L, Phung W, Lin E, et al. Fc galactosylation follows consecutive reaction kinetics and enhances immunoglobulin G hexamerization for complement activation. MAbs. 2021; 13:1893427. https://doi.org/10.1080/19420862.2021.1893427 [PubMed]
• 11. Gudelj I, Lauc G, Pezer M. Immunoglobulin G glycosylation in aging and diseases. Cell Immunol. 2018; 333:65–79. https://doi.org/10.1016/j.cellimm.2018.07.009 [PubMed]
• 12. Shkunnikova S, Mijakovac A, Sironic L, Hanic M, Lauc G, Kavur MM. IgG glycans in health and disease: Prediction, intervention, prognosis, and therapy. Biotechnol Adv. 2023; 67:108169. https://doi.org/10.1016/j.biotechadv.2023.108169 [PubMed]
• 13. Krištić J, Vučković F, Menni C, Klarić L, Keser T, Beceheli I, Pučić-Baković M, Novokmet M, Mangino M, Thaqi K, Rudan P, Novokmet N, Sarac J, et al. Glycans are a novel biomarker of chronological and biological ages. J Gerontol A Biol Sci Med Sci. 2014; 69:779–89. https://doi.org/10.1093/gerona/glt190 [PubMed]
• 14. Dall’Olio F, Vanhooren V, Chen CC, Slagboom PE, Wuhrer M, Franceschi C. N-glycomic biomarkers of biological aging and longevity: a link with inflammaging. Ageing Res Rev. 2013; 12:685–98. https://doi.org/10.1016/j.arr.2012.02.002 [PubMed]
• 15. Franceschi C, Capri M, Monti D, Giunta S, Olivieri F, Sevini F, Panourgia MP, Invidia L, Celani L, Scurti M, Cevenini E, Castellani GC, Salvioli S. Inflammaging and anti-inflammaging: a systemic perspective on aging and longevity emerged from studies in humans. Mech Ageing Dev. 2007; 128:92–105. https://doi.org/10.1016/j.mad.2006.11.016 [PubMed]
• 16. Huhn C, Selman MHJ, Ruhaak LR, Deelder AM, Wuhrer M. IgG glycosylation analysis. Proteomics. 2009; 9:882–913. https://doi.org/10.1002/pmic.200800715 [PubMed]
• 17. Gornik O, Pavić T, Lauc G. Alternative glycosylation modulates function of IgG and other proteins - implications on evolution and disease. Biochim Biophys Acta. 2012; 1820:1318–26. https://doi.org/10.1016/j.bbagen.2011.12.004 [PubMed]
• 18. Pucić M, Knezević A, Vidic J, Adamczyk B, Novokmet M, Polasek O, Gornik O, Supraha-Goreta S, Wormald MR, Redzić I, Campbell H, Wright A, Hastie ND, et al. High throughput isolation and glycosylation analysis of IgG-variability and heritability of the IgG glycome in three isolated human populations. Mol Cell Proteomics. 2011; 10:M111.010090. https://doi.org/10.1074/mcp.M111.010090 [PubMed]
• 19. Krištić J, Zaytseva OO, Ram R, Nguyen Q, Novokmet M, Vučković F, Vilaj M, Trbojević-Akmačić I, Pezer M, Davern KM, Morahan G, Lauc G. Profiling and genetic control of the murine immunoglobulin G glycome. Nat Chem Biol. 2018; 14:516–24. https://doi.org/10.1038/s41589-018-0034-3 [PubMed]
• 20. Parekh RB, Dwek RA, Sutton BJ, Fernandes DL, Leung A, Stanworth D, Rademacher TW, Mizuochi T, Taniguchi T, Matsuta K. Association of rheumatoid arthritis and primary osteoarthritis with changes in the glycosylation pattern of total serum IgG. Nature. 1985; 316:452–7. https://doi.org/10.1038/316452a0 [PubMed]
• 21. Tomana M, Schrohenloher RE, Koopman WJ, Alarcón GS, Paul WA. Abnormal glycosylation of serum IgG from patients with chronic inflammatory diseases. Arthritis Rheum. 1988; 31:333–8. https://doi.org/10.1002/art.1780310304 [PubMed]
• 22. Novak J, Tomana M, Shah GR, Brown R, Mestecky J. Heterogeneity of IgG glycosylation in adult periodontal disease. J Dent Res. 2005; 84:897–901. https://doi.org/10.1177/154405910508401005 [PubMed]
• 23. Leirisalo-Repo M, Hernandez-Munoz HE, Rook GAW. Agalactosyl IgG is elevated in patients with active spondyloarthropathy. Rheumatol Int. 1999; 18:171–6. https://doi.org/10.1007/s002960050080 [PubMed]
• 24. Holland M, Yagi H, Takahashi N, Kato K, Savage COS, Goodall DM, Jefferis R. Differential glycosylation of polyclonal IgG, IgG-Fc and IgG-Fab isolated from the sera of patients with ANCA-associated systemic vasculitis. Biochim Biophys Acta. 2006; 1760:669–77. https://doi.org/10.1016/j.bbagen.2005.11.021 [PubMed]
• 25. Parekh RB, Roitt IM, Isenberg DA, Dwek RA, Ansell BM, Rademacher TW. Galactosylation of IgG associated oligosaccharides: reduction in patients with adult and juvenile onset rheumatoid arthritis and relation to disease activity. Lancet. 1988; 1:966–9. https://doi.org/10.1016/s0140-6736(88)91781-3 [PubMed]
• 26. Siekman SL, Pongracz T, Wang W, Nouta J, Kremsner PG, da Silva-Neto PV, Esen M, Kreidenweiss A, Held J, Trapé ÁA, Fendel R, de Miranda Santos IKF, Wuhrer M, and ImmunoCovid Consortium. The IgG glycome of SARS-CoV-2 infected individuals reflects disease course and severity. Front Immunol. 2022; 13:993354. https://doi.org/10.3389/fimmu.2022.993354 [PubMed]
• 27. Bones J, Byrne JC, O’Donoghue N, McManus C, Scaife C, Boissin H, Nastase A, Rudd PM. Glycomic and glycoproteomic analysis of serum from patients with stomach cancer reveals potential markers arising from host defense response mechanisms. J Proteome Res. 2011; 10:1246–65. https://doi.org/10.1021/pr101036b [PubMed]
• 28. Arnold JN, Saldova R, Hamid UMA, Rudd PM. Evaluation of the serum N-linked glycome for the diagnosis of cancer and chronic inflammation. Proteomics. 2008; 8:3284–93. https://doi.org/10.1002/pmic.200800163 [PubMed]
• 29. Menni C, Keser T, Mangino M, Bell JT, Erte I, Akmačić I, Vučković F, Pučić Baković M, Gornik O, McCarthy MI, Zoldoš V, Spector TD, Lauc G, Valdes AM. Glycosylation of immunoglobulin g: role of genetic and epigenetic influences. PLoS One. 2013; 8:e82558. https://doi.org/10.1371/journal.pone.0082558 [PubMed]
• 30. Lauc G, Huffman JE, Pučić M, Zgaga L, Adamczyk B, Mužinić A, Novokmet M, Polašek O, Gornik O, Krištić J, Keser T, Vitart V, Scheijen B, et al. Loci associated with N-glycosylation of human immunoglobulin G show pleiotropy with autoimmune diseases and haematological cancers. PLoS Genet. 2013; 9:e1003225. https://doi.org/10.1371/journal.pgen.1003225 [PubMed]
• 31. Shen X, Klarić L, Sharapov S, Mangino M, Ning Z, Wu D, Trbojević-Akmačić I, Pučić-Baković M, Rudan I, Polašek O, Hayward C, Spector TD, Wilson JF, et al. Multivariate discovery and replication of five novel loci associated with Immunoglobulin G N-glycosylation. Nat Commun. 2017; 8:447. https://doi.org/10.1038/s41467-017-00453-3 [PubMed]
• 32. Wahl A, van den Akker E, Klaric L, Štambuk J, Benedetti E, Plomp R, Razdorov G, Trbojević-Akmačić I, Deelen J, van Heemst D, Slagboom PE, Vučković F, Grallert H, et al. Genome-Wide Association Study on Immunoglobulin G Glycosylation Patterns. Front Immunol. 2018; 9:277. https://doi.org/10.3389/fimmu.2018.00277 [PubMed]
• 33. Klarić L, Tsepilov YA, Stanton CM, Mangino M, Sikka TT, Esko T, Pakhomov E, Salo P, Deelen J, McGurnaghan SJ, Keser T, Vučković F, Ugrina I, et al. Glycosylation of immunoglobulin G is regulated by a large network of genes pleiotropic with inflammatory diseases. Sci Adv. 2020; 6:eaax0301. https://doi.org/10.1126/sciadv.aax0301 [PubMed]
• 34. Shadrina AS, Zlobin AS, Zaytseva OO, Klarić L, Sharapov SZ, D Pakhomov E, Perola M, Esko T, Hayward C, Wilson JF, Lauc G, Aulchenko YS, Tsepilov YA. Multivariate genome-wide analysis of immunoglobulin G N-glycosylation identifies new loci pleiotropic with immune function. Hum Mol Genet. 2021; 30:1259–70. https://doi.org/10.1093/hmg/ddab072 [PubMed]
• 35. Mijakovac A, Miškec K, Krištić J, Vičić Bočkor V, Tadić V, Bošković M, Lauc G, Zoldoš V, Vojta A. A Transient Expression System with Stably Integrated CRISPR-dCas9 Fusions for Regulation of Genes Involved in Immunoglobulin G Glycosylation. CRISPR J. 2022; 5:237–53. https://doi.org/10.1089/crispr.2021.0089 [PubMed]
• 36. Benedetti E, Pučić-Baković M, Keser T, Wahl A, Hassinen A, Yang JY, Liu L, Trbojević-Akmačić I, Razdorov G, Štambuk J, Klarić L, Ugrina I, Selman MHJ, et al. Network inference from glycoproteomics data reveals new reactions in the IgG glycosylation pathway. Nat Commun. 2017; 8:1483. https://doi.org/10.1038/s41467-017-01525-0 [PubMed]
• 37. Watanabe K, Taskesen E, van Bochoven A, Posthuma D. Functional mapping and annotation of genetic associations with FUMA. Nat Commun. 2017; 8:1826. https://doi.org/10.1038/s41467-017-01261-5 [PubMed]
• 38. McLaren W, Gil L, Hunt SE, Riat HS, Ritchie GRS, Thormann A, Flicek P, Cunningham F. The Ensembl Variant Effect Predictor. Genome Biol. 2016; 17:122. https://doi.org/10.1186/s13059-016-0974-4 [PubMed]
• 39. Võsa U, Claringbould A, Westra HJ, Bonder MJ, Deelen P, Zeng B, Kirsten H, Saha A, Kreuzhuber R, Yazar S, Brugge H, Oelen R, de Vries DH, et al, and BIOS Consortium, and i2QTL Consortium. Large-scale cis- and trans-eQTL analyses identify thousands of genetic loci and polygenic scores that regulate blood gene expression. Nat Genet. 2021; 53:1300–10. https://doi.org/10.1038/s41588-021-00913-z [PubMed]
• 40. Giambartolomei C, Vukcevic D, Schadt EE, Franke L, Hingorani AD, Wallace C, Plagnol V. Bayesian test for colocalisation between pairs of genetic association studies using summary statistics. PLoS Genet. 2014; 10:e1004383. https://doi.org/10.1371/journal.pgen.1004383 [PubMed]
• 41. Lloyd-Jones LR, Zeng J, Sidorenko J, Yengo L, Moser G, Kemper KE, Wang H, Zheng Z, Magi R, Esko T, Metspalu A, Wray NR, Goddard ME, et al. Improved polygenic prediction by Bayesian multiple regression on summary statistics. Nat Commun. 2019; 10:5086. https://doi.org/10.1038/s41467-019-12653-0 [PubMed]
• 42. Mijakovac A, Frkatović A, Hanić M, Ivok J, Martinić Kavur M, Pučić-Baković M, Spector T, Zoldoš V, Mangino M, Lauc G. Heritability of the glycan clock of biological age. Front Cell Dev Biol. 2022; 10:982609. https://doi.org/10.3389/fcell.2022.982609 [PubMed]
• 43. Mijakovac A, Jurić J, Kohrt WM, Krištić J, Kifer D, Gavin KM, Miškec K, Frkatović A, Vučković F, Pezer M, Vojta A, Nigrović PA, Zoldoš V, Lauc G. Effects of Estradiol on Immunoglobulin G Glycosylation: Mapping of the Downstream Signaling Mechanism. Front Immunol. 2021; 12:680227. https://doi.org/10.3389/fimmu.2021.680227 [PubMed]
• 44. Nair M, Teng A, Bilanchone V, Agrawal A, Li B, Dai X. Ovol1 regulates the growth arrest of embryonic epidermal progenitor cells and represses c-myc transcription. J Cell Biol. 2006; 173:253–64. https://doi.org/10.1083/jcb.200508196 [PubMed]
• 45. Boettger LM, Handsaker RE, Zody MC, McCarroll SA. Structural haplotypes and recent evolution of the human 17q21.31 region. Nat Genet. 2012; 44:881–5. https://doi.org/10.1038/ng.2334 [PubMed]
• 46. Yamashita J, Iwamura C, Mitsumori K, Hosokawa H, Sasaki T, Takahashi M, Tanaka H, Kaneko K, Hanazawa A, Watanabe Y, Shinoda K, Tumes D, Motohashi S, Nakayama T. Murine Schnurri-2 controls natural killer cell function and lymphoma development. Leuk Lymphoma. 2012; 53:479–86. https://doi.org/10.3109/10428194.2011.625099 [PubMed]
• 47. Schumann K, Raju SS, Lauber M, Kolb S, Shifrut E, Cortez JT, Skartsis N, Nguyen VQ, Woo JM, Roth TL, Yu R, Nguyen MLT, Simeonov DR, et al. Functional CRISPR dissection of gene networks controlling human regulatory T cell identity. Nat Immunol. 2020; 21:1456–66. https://doi.org/10.1038/s41590-020-0784-4 [PubMed]
• 48. Hitomi Y, Nakatani K, Kojima K, Nishida N, Kawai Y, Kawashima M, Aiba Y, Nagasaki M, Nakamura M, Tokunaga K. NFKB1 and MANBA Confer Disease Susceptibility to Primary Biliary Cholangitis via Independent Putative Primary Functional Variants. Cell Mol Gastroenterol Hepatol. 2019; 7:515–32. https://doi.org/10.1016/j.jcmgh.2018.11.006 [PubMed]
• 49. Zeng X, Li S, Tang S, Li X, Zhang G, Li M, Zeng X, Hu C. Changes of Serum IgG Glycosylation Patterns in Primary Biliary Cholangitis Patients. Front Immunol. 2021; 12:669137. https://doi.org/10.3389/fimmu.2021.669137 [PubMed]
• 50. Cartwright T, Perkins ND, L Wilson C. NFKB1: a suppressor of inflammation, ageing and cancer. FEBS J. 2016; 283:1812–22. https://doi.org/10.1111/febs.13627 [PubMed]
• 51. Zhu M, Lovell KL, Patterson JS, Saunders TL, Hughes ED, Friderici KH. Beta-mannosidosis mice: a model for the human lysosomal storage disease. Hum Mol Genet. 2006; 15:493–500. https://doi.org/10.1093/hmg/ddi465 [PubMed]
• 52. Gu X, Yang H, Sheng X, Ko YA, Qiu C, Park J, Huang S, Kember R, Judy RL, Park J, Damrauer SM, Nadkarni G, Loos RJF, et al. Kidney disease genetic risk variants alter lysosomal beta-mannosidase (MANBA) expression and disease severity. Sci Transl Med. 2021; 13:eaaz1458. https://doi.org/10.1126/scitranslmed.aaz1458 [PubMed]
• 53. Guo Y, Rist PM, Daghlas I, Giulianini F, Kurth T, Chasman DI, and International Headache Genetics Consortium, and 23andMe Research Team. A genome-wide cross-phenotype meta-analysis of the association of blood pressure with migraine. Nat Commun. 2020; 11:3368. https://doi.org/10.1038/s41467-020-17002-0 [PubMed]
• 54. Rouillard AD, Gundersen GW, Fernandez NF, Wang Z, Monteiro CD, McDermott MG, Ma’ayan A. The harmonizome: a collection of processed datasets gathered to serve and mine knowledge about genes and proteins. Database (Oxford). 2016; 2016:baw100. https://doi.org/10.1093/database/baw100 [PubMed]
• 55. Wittenbecher C, Štambuk T, Kuxhaus O, Rudman N, Vučković F, Štambuk J, Schiborn C, Rahelić D, Dietrich S, Gornik O, Perola M, Boeing H, Schulze MB, Lauc G. Plasma N-Glycans as Emerging Biomarkers of Cardiometabolic Risk: A Prospective Investigation in the EPIC-Potsdam Cohort Study. Diabetes Care. 2020; 43:661–8. https://doi.org/10.2337/dc19-1507 [PubMed]
• 56. Birukov A, Plavša B, Eichelmann F, Kuxhaus O, Hoshi RA, Rudman N, Štambuk T, Trbojević-Akmačić I, Schiborn C, Morze J, Mihelčić M, Cindrić A, Liu Y, et al. Immunoglobulin G N-Glycosylation Signatures in Incident Type 2 Diabetes and Cardiovascular Disease. Diabetes Care. 2022; 45:2729–36. https://doi.org/10.2337/dc22-0833 [PubMed]
• 57. Browne GJ, Proud CG. Regulation of peptide-chain elongation in mammalian cells. Eur J Biochem. 2002; 269:5360–8. https://doi.org/10.1046/j.1432-1033.2002.03290.x [PubMed]
• 58. Li D, Wei T, Abbott CM, Harrich D. The unexpected roles of eukaryotic translation elongation factors in RNA virus replication and pathogenesis. Microbiol Mol Biol Rev. 2013; 77:253–66. https://doi.org/10.1128/MMBR.00059-12 [PubMed]
• 59. Abbas W, Kumar A, Herbein G. The eEF1A Proteins: At the Crossroads of Oncogenesis, Apoptosis, and Viral Infections. Front Oncol. 2015; 5:75. https://doi.org/10.3389/fonc.2015.00075 [PubMed]
• 60. Schulz I, Engel C, Niestroj AJ, Kehlen A, Rahfeld JU, Kleinschmidt M, Lehmann K, Roßner S, Demuth HU. A non-canonical function of eukaryotic elongation factor 1A1: regulation of interleukin-6 expression. Biochim Biophys Acta. 2014; 1843:965–75. https://doi.org/10.1016/j.bbamcr.2014.01.022 [PubMed]
• 61. Hirano T. IL-6 in inflammation, autoimmunity and cancer. Int Immunol. 2021; 33:127–48. https://doi.org/10.1093/intimm/dxaa078 [PubMed]
• 62. Troelsen LN, Jacobsen S, Abrahams JL, Royle L, Rudd PM, Narvestad E, Heegaard NH, Garred P. IgG glycosylation changes and MBL2 polymorphisms: associations with markers of systemic inflammation and joint destruction in rheumatoid arthritis. J Rheumatol. 2012; 39:463–9. https://doi.org/10.3899/jrheum.110584 [PubMed]
• 63. Ozcan E, Garibyan L, Lee JJY, Bram RJ, Lam KP, Geha RS. Transmembrane activator, calcium modulator, and cyclophilin ligand interactor drives plasma cell differentiation in LPS-activated B cells. J Allergy Clin Immunol. 2009; 123:1277–86.e5. https://doi.org/10.1016/j.jaci.2009.03.019 [PubMed]
• 64. Castigli E, Wilson SA, Elkhal A, Ozcan E, Garibyan L, Geha RS. Transmembrane activator and calcium modulator and cyclophilin ligand interactor enhances CD40-driven plasma cell differentiation. J Allergy Clin Immunol. 2007; 120:885–91. https://doi.org/10.1016/j.jaci.2007.06.012 [PubMed]
• 65. Jonsson S, Sveinbjornsson G, de Lapuente Portilla AL, Swaminathan B, Plomp R, Dekkers G, Ajore R, Ali M, Bentlage AEH, Elmér E, Eyjolfsson GI, Gudjonsson SA, Gullberg U, et al. Identification of sequence variants influencing immunoglobulin levels. Nat Genet. 2017; 49:1182–91. https://doi.org/10.1038/ng.3897 [PubMed]
• 66. Martinez-Gallo M, Radigan L, Almejún MB, Martínez-Pomar N, Matamoros N, Cunningham-Rundles C. TACI mutations and impaired B-cell function in subjects with CVID and healthy heterozygotes. J Allergy Clin Immunol. 2013; 131:468–76. https://doi.org/10.1016/j.jaci.2012.10.029 [PubMed]
• 67. Poodt AEJ, Driessen GJA, de Klein A, van Dongen JJM, van der Burg M, de Vries E. TACI mutations and disease susceptibility in patients with common variable immunodeficiency. Clin Exp Immunol. 2009; 156:35–9. https://doi.org/10.1111/j.1365-2249.2008.03863.x [PubMed]
• 68. Johnson WE, Li C, Rabinovic A. Adjusting batch effects in microarray expression data using empirical Bayes methods. Biostatistics. 2007; 8:118–27. https://doi.org/10.1093/biostatistics/kxj037 [PubMed]
• 69. Leek JT. svaseq: removing batch effects and other unwanted noise from sequencing data. Nucleic Acids Res. 2014; 42:e161. https://doi.org/10.1093/nar/gku864 [PubMed]
• 70. R Core Team. R: A language and environment for statistical computing. [Internet]. Vienna, Austria: R Foundation for Statistical Computing. 2017.
• 71. Aulchenko YS, Ripke S, Isaacs A, van Duijn CM. GenABEL: an R library for genome-wide association analysis. Bioinformatics. 2007; 23:1294–6. https://doi.org/10.1093/bioinformatics/btm108 [PubMed]
• 72. Haller T, Kals M, Esko T, Mägi R, Fischer K. RegScan: a GWAS tool for quick estimation of allele effects on continuous traits and their combinations. Brief Bioinform. 2015; 16:39–44. https://doi.org/10.1093/bib/bbt066 [PubMed]
• 73. Zhou X, Stephens M. Genome-wide efficient mixed-model analysis for association studies. Nat Genet. 2012; 44:821–4. https://doi.org/10.1038/ng.2310 [PubMed]
• 74. McCarthy S, Das S, Kretzschmar W, Delaneau O, Wood AR, Teumer A, Kang HM, Fuchsberger C, Danecek P, Sharp K, Luo Y, Sidore C, Kwong A, et al, and Haplotype Reference Consortium. A reference panel of 64,976 haplotypes for genotype imputation. Nat Genet. 2016; 48:1279–83. https://doi.org/10.1038/ng.3643 [PubMed]
• 75. Winkler TW, Day FR, Croteau-Chonka DC, Wood AR, Locke AE, Mägi R, Ferreira T, Fall T, Graff M, Justice AE, Luan J, Gustafsson S, Randall JC, et al, and Genetic Investigation of Anthropometric Traits (GIANT) Consortium. Quality control and conduct of genome-wide association meta-analyses. Nat Protoc. 2014; 9:1192–212. https://doi.org/10.1038/nprot.2014.071 [PubMed]
• 76. Willer CJ, Li Y, Abecasis GR. METAL: fast and efficient meta-analysis of genomewide association scans. Bioinformatics. 2010; 26:2190–1. https://doi.org/10.1093/bioinformatics/btq340 [PubMed]
• 77. Bycroft C, Freeman C, Petkova D, Band G, Elliott LT, Sharp K, Motyer A, Vukcevic D, Delaneau O, O’Connell J, Cortes A, Welsh S, Young A, et al. The UK Biobank resource with deep phenotyping and genomic data. Nature. 2018; 562:203–9. https://doi.org/10.1038/s41586-018-0579-z [PubMed]
• 78. Wang K, Li M, Hakonarson H. ANNOVAR: functional annotation of genetic variants from high-throughput sequencing data. Nucleic Acids Res. 2010; 38:e164. https://doi.org/10.1093/nar/gkq603 [PubMed]
• 79. Schmiedel BJ, Singh D, Madrigal A, Valdovino-Gonzalez AG, White BM, Zapardiel-Gonzalo J, Ha B, Altay G, Greenbaum JA, McVicker G, Seumois G, Rao A, Kronenberg M, et al. Impact of Genetic Polymorphisms on Human Immune Cell Gene Expression. Cell. 2018; 175:1701–15.e16. https://doi.org/10.1016/j.cell.2018.10.022 [PubMed]
• 80. Fairfax BP, Humburg P, Makino S, Naranbhai V, Wong D, Lau E, Jostins L, Plant K, Andrews R, McGee C, Knight JC. Innate immune activity conditions the effect of regulatory variants upon monocyte gene expression. Science. 2014; 343:1246949. https://doi.org/10.1126/science.1246949 [PubMed]
• 81. Momozawa Y, Dmitrieva J, Théâtre E, Deffontaine V, Rahmouni S, Charloteaux B, Crins F, Docampo E, Elansary M, Gori AS, Lecut C, Mariman R, Mni M, et al, and International IBD Genetics Consortium. IBD risk loci are enriched in multigenic regulatory modules encompassing putative causative genes. Nat Commun. 2018; 9:2427. https://doi.org/10.1038/s41467-018-04365-8 [PubMed]
• 82. Rao SSP, Huntley MH, Durand NC, Stamenova EK, Bochkov ID, Robinson JT, Sanborn AL, Machol I, Omer AD, Lander ES, Aiden EL. A 3D map of the human genome at kilobase resolution reveals principles of chromatin looping. Cell. 2014; 159:1665–80. https://doi.org/10.1016/j.cell.2014.11.021 [PubMed]
• 83. Wakefield J. Bayes factors for genome-wide association studies: comparison with P-values. Genet Epidemiol. 2009; 33:79–86. https://doi.org/10.1002/gepi.20359 [PubMed]
• 84. de Leeuw CA, Mooij JM, Heskes T, Posthuma D. MAGMA: generalized gene-set analysis of GWAS data. PLoS Comput Biol. 2015; 11:e1004219. https://doi.org/10.1371/journal.pcbi.1004219 [PubMed]
• 85. Chang CC, Chow CC, Tellier LCA, Vattikuti S, Purcell SM, Lee JJ. Second-generation PLINK: rising to the challenge of larger and richer datasets. Gigascience. 2015; 4:7. https://doi.org/10.1186/s13742-015-0047-8 [PubMed]
• 86. Vink T, Oudshoorn-Dickmann M, Roza M, Reitsma JJ, de Jong RN. A simple, robust and highly efficient transient expression system for producing antibodies. Methods. 2014; 65:5–10. https://doi.org/10.1016/j.ymeth.2013.07.018 [PubMed]
• 87. Livak KJ, Schmittgen TD. Analysis of relative gene expression data using real-time quantitative PCR and the 2(-Delta Delta C(T)) Method. Methods. 2001; 25:402–8. https://doi.org/10.1006/meth.2001.1262 [PubMed]
• 88. Viechtbauer W. Conducting meta-analyses in {R} with the {metafor} package. J Stat Softw. 2010; 36:1–48. https://doi.org/10.18637/jss.v036.i03