Integrative genomic deconvolution of rheumatoid arthritis GWAS loci into gene and cell type associations

Abstract

Background: Although genome-wide association studies (GWAS) have identified over 100 genetic loci associated with rheumatoid arthritis (RA), our ability to translate these results into disease understanding and novel therapeutics is limited. Most RA GWAS loci reside outside of protein-coding regions and likely affect distal transcriptional enhancers. Furthermore, GWAS do not identify the cell types where the associated causal gene functions. Thus, mapping the transcriptional regulatory roles of GWAS hits and the relevant cell types will lead to better understanding of RA pathogenesis.

Results: We combine the whole-genome sequences and blood transcription profiles of 377 RA patients and identify over 6000 unique genes with expression quantitative trait loci (eQTLs). We demonstrate the quality of the identified eQTLs through comparison to non-RA individuals. We integrate the eQTLs with immune cell epigenome maps, RA GWAS risk loci, and adjustment for linkage disequilibrium to propose target genes of immune cell enhancers that overlap RA risk loci. We examine 20 immune cell epigenomes and perform a focused analysis on primary monocytes, B cells, and T cells.

Conclusions: We highlight cell-specific gene associations with relevance to RA pathogenesis including the identification of FCGR2B in B cells as possessing both intragenic and enhancer regulatory GWAS hits. We show that our RA patient cohort derived eQTL network is more informative for studying RA than that from a healthy cohort. While not experimentally validated here, the reported eQTLs and cell type-specific RA risk associations can prioritize future experiments with the goal of elucidating the regulatory mechanisms behind genetic risk associations.

Keywords: Epigenomics; Expression quantitative trait loci (eQTLs); Genome-wide association studies (GWAS); Rheumatoid arthritis.

Fig. 1

Overview of eQTLs from the blood of RA subjects and canonical pathway enrichment of genes with eQTLs. a Circos plot for eQTLs mapped in whole blood from RA subjects. Outermost track of gray dots shows a Manhattan plot of cis-eQTLs (FDR < 5 %) with egenes with lowest p values labeled. The inner links connect trans-eQTLs to the associated egene (only trans-eQTLs with nominal p values <10⁻¹⁴ are shown for clarity). Blue links represent associations within the same chromosome whereas gray links represent associations that span across chromosomes. b egenes associated with eQTLs in RA whole blood samples were tested for enrichment with the MSigDB term database separately for egenes associated in cis (left) or in trans (right)

Fig. 2

Comparison of RA eQTLs to published GWAS catalog and fold-enrichment of RA eQTLs in chromatin states from diverse tissues. a Comparison of RA cis-eQTLs to reported GWAS SNPs in the NHGRI catalog. The top 17 diseases or traits are shown, ranked by enrichment p value from one-sided Fisher’s exact test. The number of GWAS reported SNPs overlapping with RA eQTLs in our dataset are shown for each disease or trait. b Comparison of enrichment of GWAS SNPs in several eQTL studies. The top ten GWAS diseases or traits are shown. The eQTL studies are labeled by their first author and cell type or tissue where gene expression was measured. Only cis-eQTLs were considered. c Enrichment of RA eQTLs in chromatin states. Both cis- and trans-eQTLs were compared to positions of the following chromatin states in diverse tissues: transcription start sites (TSS), actively transcribed (Tx), enhancers (Enh), active enhancers (EnhA), heterochromatin (Het), bivalent enhancers (EnhBiv), bivalent/poised TSS (TssBiv), repressed Polycomb (ReprPC), and quiescent/low (Quies) regions. Only a subset of the cell or tissue types with the highest or lowest fold-enrichment in enhancer states are shown to ease representation. The full results are shown in Additional file 5: Figure S2

Fig. 3

Integration of epigenomics datasets allows cell-type annotation for RA GWAS and RA eQTL connections. a Schematic demonstrating the approach to connect reported GWAS SNPs of unknown function to putative causal genes using epigenomics and eQTL datasets. b Heatmap of genes identified to be associated with RA GWAS SNPs in primary monocyte, T cell, and B cell datasets

Fig. 4

Canonical pathway and cell-type enrichment of identified RA GWAS-associated genes and enrichment with RA gene expression disease profiles. Genes identified as associated with RA GWAS SNPs in monocytes, B cells, or T cells were tested for enrichment with the MSigDB term database a or a database of cell type-specific gene expression b. c The genes identified were compared to genes differentially expressed in three publicly available RA synovial tissue datasets with the GEO IDs indicated. The enrichment p value is the unadjusted p value from one-sided Fisher’s exact test. d The enrichment with the increased genes in the same RA synovial tissue datasets were compared with gene lists generated by choosing the closest genes to the RA GWAS SNPs (“closest genes”; the same number of genes are selected per SNP as with our method) or choosing any gene(s) within 5000 bp of the GWAS SNP (“within 5000 bp”)

Fig. 5

Example RA risk enhancers and their target genes. a A genome browser view showing a 250 kp window around the RA risk loci rs2317230. The RA GWAS hit and SNPs in LD are shown in red and overlap two immune enhancer regions. Yellow and orange are used to highlight enhancer locations and point to their target genes whose active TSSs are highlighted with gray. Histone modification tracks show the fold change between ChIP and control read counts. Active histone modifications have red bars next to their labels while repressive marks have blue. Genes on the positive and negative strands are shown in blue and cyan, respectively. b A genome browser view is shown as in a, but around the RA risk loci rs72717009. c Gene expression (normalized log2-transformed values) is plotted versus the genotype at representative enhancer-localized SNPs for the same genes shown in a and b