Next-Generation Sequencing Profiles of the Methylome and Transcriptome in Peripheral Blood Mononuclear Cells of Rheumatoid Arthritis

Abstract

Using next-generation sequencing to decipher methylome and transcriptome and underlying molecular mechanisms contributing to rheumatoid arthritis (RA) for improving future therapies, we performed methyl-seq and RNA-seq on peripheral blood mononuclear cells (PBMCs) from RA subjects and normal donors. Principal component analysis and hierarchical clustering revealed distinct methylation signatures in RA with methylation aberrations noted across chromosomes. Methylation alterations varied with CpG features and genic characteristics. Typically, CpG islands and CpG shores were hypermethylated and displayed the greatest methylation variance. Promoters were hypermethylated and enhancers/gene bodies were hypomethylated, with methylation variance associated with expression variance. RA genetically associated genes preferentially displayed differential methylation and differential expression or interacted with differentially methylated and differentially expressed genes. These differentially methylated and differentially expressed genes were enriched with several signaling pathways and disease categories. 10 genes (CD86, RAB20, XAF1, FOLR3, LTBR, KCNH8, DOK7, PDGFA, PITPNM2, CELSR1) with concomitantly differential methylation in enhancers/promoters/gene bodies and differential expression in B cells were validated. This integrated analysis of methylome and transcriptome identified novel epigenetic signatures associated with RA and highlighted the interaction between genetics and epigenetics in RA. These findings help our understanding of the pathogenesis of RA and advance epigenetic studies in regards to the disease.

Keywords: methylation; next-generation sequencing; rheumatoid arthritis.

Figure 1

Schematic representation of the next-generation sequencing data analytical workflow. After adjusting for cellular composition and batch effects, methyl-seq data first underwent Principal component analysis (PCA) and hierarchical clustering (HC) (Step 1), OmicCircos visualization (Step 2), CpG features mapping (Step 3), and genic characteristics annotation (Step 4). Methylation and expression profiles were then integrated for methylation-expression correlation (Step 5). Differentially methylated genes (FDR < 0.05) and differentially expressed genes (FDR < 0.05) were identified (Step 6a–6b) and intersected to yield genes with concomitant expression and methylation changes in enhancer/promoter/gene body (Step 7). These differentially methylated and differentially expressed genes underwent genetic–epigenetic interaction investigation (Step 8), IPA (Step 9), and upstream regulator deduction (Step 10). GEO dataset validation (Step 11) confirmed concomitant differential methylation and expression of 10 genes.

Figure 2

Methylation differences according to CpG features and genic characteristics. The bar charts showed the methylation difference (rheumatoid arthritis (R) minus healthy donor (H)) in CpG island, CpG shore, CpG shelf, open sea (a) and variance of methylation according to respective CpG features (b). Methylation difference in intergenic region, enhancer, promoter, gene body (c) and variance of methylation in respective genic characteristics (d) were also presented. * p < 0.001 for methylation difference and variance of methylation between different CpG features and genic characteristics.

Figure 3

Transcription factors identified through iRegulon analysis. The bubble chart showed the transcription factors associated with differentially methylated and differentially expressed genes identified by iRegulon. Y-axis label represented normalized enrichment score. The sizes of the bubbles were proportional to the number of regulated genes with concomitant differential methylation and differential expression for each transcription factor.

Figure 4

Validation of genes with differential methylation in enhancer and differential expression. (a) The results of methylation and expression obtained from next-generation sequencing (NGS meth, NGS exp), the cell subsets of validation dataset (Cell), the dataset of validation (Meth dataset, Exp dataset), and the probes of validation dataset (CpG probe, Exp probe). (b) Visualization of the methylation levels obtained from NGS in rheumatoid arthritis (RA) and healthy donors (HD) and location of validated CpG probe and enhancers. (c) Volcano plot of the −log₁₀(false discovery rate) on the Y-axis versus expression change (log₂ratio) on the X-axis. Of validated genes, (d) Methylation and (e) Expression levels of corresponding probes in the validation dataset.

Figure 5

Validation of genes with differential methylation in promoter and differential expression. (a) The results of methylation and expression obtained from next-generation sequencing (NGS meth, NGS exp), the cell subsets of validation dataset (Cell), the dataset of validation (Meth dataset, Exp dataset), and the probes of validation dataset (CpG probe, Exp probe). (b) Visualization of the methylation levels obtained from NGS in RA and healthy donors (HD) and location of validated CpG probe superposed onto the genomic locations of genes. (c) Volcano plot of the -log₁₀(false discovery rate) on the Y-axis versus expression change (log₂ratio) on the X-axis. Of validated genes, (d) Methylation and (e) Expression levels of corresponding probes in the validation dataset.

Figure 6

Validation of genes with differential methylation in gene body and differential expression. (a) The results of methylation and expression obtained from next-generation sequencing (NGS meth, NGS exp), the cell subsets of validation dataset (Cell), the dataset of validation (Meth dataset, Exp dataset), and the probes of validation dataset (CpG probe, Exp probe). (b) Visualization of the methylation levels obtained from NGS in RA and healthy donors (HD) and location of validated CpG probe superposed onto the genomic locations of genes. (c) Volcano plot of the -log₁₀(false discovery rate) on the Y-axis versus expression change (log₂ratio) on the X-axis. Of validated genes, (d) Methylation and (e) Expression levels of corresponding probes in the validation dataset.