Revealing the Immune Heterogeneity between Systemic Lupus Erythematosus and Rheumatoid Arthritis Based on Multi-Omics Data Analysis

Abstract

The pathogenesis of systemic lupus erythematosus (SLE) and rheumatoid arthritis (RA) are greatly influenced by different immune cells. Nowadays both T-cell receptor (TCR) and B-cell receptor (BCR) sequencing technology have emerged with the maturity of NGS technology. However, both SLE and RA peripheral blood TCR or BCR repertoire sequencing remains lacking because repertoire sequencing is an expensive assay and consumes valuable tissue samples. This study used computational methods TRUST4 to construct TCR repertoire and BCR repertoire from bulk RNA-seq data of both SLE and RA patients' peripheral blood and analyzed the clonality and diversity of the immune repertoire between the two diseases. Although the functions of immune cells have been studied, the mechanism is still complicated. Differentially expressed genes in each immune cell type and cell-cell interactions between immune cell clusters have not been covered. In this work, we clustered eight immune cell subsets from original scRNA-seq data and disentangled the characteristic alterations of cell subset proportion under both SLE and RA conditions. The cell-cell communication analysis tool CellChat was also utilized to analyze the influence of MIF family and GALECTIN family cytokines, which were reported to regulate SLE and RA, respectively. Our findings correspond to previous findings that MIF increases in the serum of SLE patients. This work proved that the presence of LGALS9, PTPRC and CD44 in platelets could serve as a clinical indicator of rheumatoid arthritis. Our findings comprehensively illustrate dynamic alterations in immune cells during pathogenesis of SLE and RA. This work identified specific V genes and J genes in TCR and BCR that could be used to expand our understanding of SLE and RA. These findings provide a new insight inti the diagnosis and treatment of the two autoimmune diseases.

Keywords: cell–cell interaction; immune repertoire; rheumatoid arthritis; single-cell RNA sequencing technology; systemic lupus erythematosus.

Figure 1

Flowchart for this paper. The main framework includes immune repertoire construction and single-cell RNA-seq data analysis. During the process of immune repertoire construction, bulk RNA-seq data were used as inputs for TRUST4 algorithm. After building both T-cell receptor and B-cell receptor repertoires, we used the extracted information to perform further statistical analysis. During the process of single cell RNA-seq data analysis, we utilized the original data to perform cell clustering and compared the proportions of various of immune cell types under SLE, RA and HC conditions. Finally, DEGs identification and functional enrichment analysis were conducted for different immune cell types, respectively, under SLE, RA and HC conditions. Cell–cell interaction analysis of MIF and GALECTIN family signaling pathways was also implemented between several immune cell types.

Figure 2

T-cell receptor repertoire construction and analysis. (A) Box plot showing the top 10 high-frequency T-cell clonotypes for each sample under SLE, RA and HC conditions in this study. (B) Box plots of Shannon–Weiner index for each sample under SLE, RA and HC conditions are used to compare the TCR diversity between SLE, RA and HC groups. (C) Distributions of TCR CDR3 amino acid sequence length in SLE, RA and HC groups. (D) Box plots of InvSimpson index for each sample under SLE, RA and HC conditions are used to compare the TCR diversity between SLE, RA and HC groups. (E) TRBV gene-usage-frequency stacked histogram showing the distributions of common TRBV gene in the SLE, RA and HC groups, respectively. (F) TRBJ gene-usage-frequency stacked histogram showing the distributions of common TRBJ genes in the SLE, RA and HC groups, respectively.

Figure 3

B-cell receptor repertoire construction and analysis. (A) Box plot showing the top 10 high-frequency B-cell clonotypes for each sample under SLE, RA and HC conditions in this study. (B) Box plots of Shannon–Weiner index for each sample under SLE, RA and HC conditions are used to compare the BCR diversity between SLE, RA and HC groups. (C) Distributions of BCR CDR3 amino acid sequence length in the SLE, RA and HC groups. (D) Box plots of InvSimpson index for each sample under SLE, RA and HC conditions are used to compare the BCR diversity between SLE, RA and HC groups. (E) IGHV gene-usage-frequency stacked histogram showing the distributions of common IGHV genes in the SLE, RA and HC groups, respectively. (F) IGHJ gene-usage-frequency stacked histogram showing the distributions of common IGHJ genes in the SLE, RA and HC groups, respectively.

Figure 4

Single-cell RNA-seq data analysis of peripheral blood mononuclear cells (PBMC) in SLE, RA and HC groups. (A) The results of using canonical correction analysis (CCA) method to remove batch effects after integrating all the SLE, RA and HC scRNA-seq data. (B) TSNE plot representing the nine clusters across 39,446 PBMCs from seven samples (three SLE samples, three HC samples and one RA sample). (C) Representative marker genes define CD4+ T cells, CD8+ T cells, CD14+ monocytes, FCGR3A+ monocytes, natural killer cells, dendritic cells, platelets and B cells, respectively. (D) Bar plots showing the proportion of different immune cell types in each sample.

Figure 5

Differentially expressed gene identification and functional and pathway enrichment analysis of distinct immune cells. (A) Heat map of the top featured marker genes for distinct immune cell types. (FDR < 0.05, logFC > 1) (B) Heat map of the enriched functional and signaling pathways of differentially expressed genes between SLE group and HC group in distinct immune cell types. The heat map cells are colored according to the p-value of the enriched terms, and white cells indicate a lack of enrichment for that term. (C) The Circos plot shows how differentially expressed genes between SLE group and HC group from the given immune cell types overlap. Each arc represents each gene list’s identity. Purple lines link the same gene that is shared by multiple gene lists. Blue lines link the different genes where they fall into the same ontological term (the term needs to be statistically significantly enriched with a size no larger than 100). (D) The Circos plot shows how differentially expressed genes between RA group and HC group from the given immune cell types overlap. (E) Heat map of the enriched functional and signaling pathways of differentially expressed genes between RA group and HC group in distinct immune cell types.

Figure 6

Cell-to-cell communications of MIF family signaling pathways among the distinct immune cells in the PBMC between SLE group and HC group predicted by the CellChat software. (A) Chord plots showing the interactions of ligand/receptor pairs MIF–(CD74 + CXCR4) between SLE group and HC group. (B) Circle plots summarizing the interactions of MIF signaling pathways among individual cell types in both SLE and HC groups. (C) Heat map showing the relative contribution of each cell type based on the computed four-network centrality measures of MIF signaling network between SLE and HC groups. (D) Violin plots showing the expression levels of MIF family cytokines in each immune cell type for both SLE and HC groups.

Figure 7

Cell-to-cell communications of GALECTIN family signaling pathways among the distinct immune cells in the PBMC between RA group and HC group predicted by the CellChat software. (A) Chord plots showing the interactions of ligand/receptor pair LGALS9–CD45 between RA group and HC group. (B) Circle plots summarizing the interactions of GALECTIN signaling pathway among individual cell types in both RA and HC group. (C) Heat map showing the relative contribution of each cell type based on the computed four-network centrality measures of GALECTIN signaling network between RA and HC groups. (D) Violin plots showing the expression levels of GALECTIN family cytokines in each immune cell type for both RA and HC groups.