Transcriptomic signatures of classical monocytes reveal pro-inflammatory modules and heterogeneity in polyarticular juvenile idiopathic arthritis

Abstract

Introduction: Polyarticular juvenile idiopathic arthritis (pJIA) is a childhood-onset autoimmune disease. Immune cells contribute to persistent inflammation observed in pJIA. Despite the crucial role of monocytes in arthritis, the precise involvement of classical monocytes in the pathogenesis of pJIA remains uncertain. Here, we aimed to uncover the transcriptomic patterns of classical monocytes in pJIA, focusing on their involvement in disease mechanism and heterogeneity.

Methods: A total of 17 healthy subjects and 18 premenopausal women with pJIA according to ILAR criteria were included. Classical monocytes were isolated, and RNA sequencing was performed. Differential expression analysis was used to compare pJIA patients and healthy control group. Differentially expressed genes (DEGs) were identified, and gene set enrichment analysis (GSEA) was performed. Using unsupervised learning approach, patients were clustered in two groups based on their similarities at transcriptomic level. Subsequently, these clusters underwent a comparative analysis to reveal differences at the transcriptomic level.

Results: We identified 440 DEGs in pJIA patients of which 360 were upregulated and 80 downregulated. GSEA highlighted TNF-α and IFN-γ response. Importantly, this analysis not only detected genes targeted by pJIA therapy but also identified new modulators of immuno-inflammation. PLAUR, IL1B, IL6, CDKN1A, PIM1, and ICAM1 were pointed as drivers of chronic hyperinflammation. Unsupervised learning approach revealed two clusters within pJIA, each exhibiting varying inflammation levels.

Conclusion: These findings indicate the pivotal role of immuno-inflammation driven by classical monocytes in pJIA and reveals the existence of two subclusters within pJIA, regardless the positivity of rheumatoid factor and anti-CCP, paving the way to precision medicine.

Keywords: autoimmunity; classical monocytes; inflammation; polyarticular juvenile idiopathic arthritis; transcriptomic.

Figure 1

Differential expression analysis of classical monocytes’ genes in pJIA reveals similarities with RA. (A) Volcano plot showing DEGs identified from the comparison of pJIA vs. control group. A total of 440 DEGs were highlighted, of which 360 were upregulated (red dots) and 80 downregulated (blue dots). DEGs were identified based on a fold change of 2 and a Benjamini–Hochberg false discovery rate (FDR) using a cutoff set at <0.05. (B) Pearson’s correlation between the fold changes of pJIA and RA, both compared with control group, revealed high positive correlation. Genes involved in inflammatory mechanisms are highlighted in the plot. (C) Heatmap showing the expression profiles of the top 50 DEGs and unsupervised clustering of sample. This reveals a clear stratification of pJIA and controls in different clusters and a heterogeneity within the pJIA group. Hierarchical clustering of samples was performed based on the Euclidean distance calculated from the normalized and scaled expression. (D) Forest plot showing genes previously associated with pJIA in Open Target database. Log2 of fold change and standard errors are plotted. DEGs, differentially expressed genes; pJIA, polyarticular juvenile idiopathic arthritis; padj, adjusted p-value; RA, rheumatoid arthritis.

Figure 2

Functional gene set enrichment analysis and pJIA-associated gene modules. (A) Dot plot of functional analysis performed using the hallmark gene sets (GSEA) shows the enrichment of pathways associated with the activation of immuno-inflammation. Dot size represents the normalized GSEA enrichment score (NES), and the X-axis indicates the FDR. Pathways with FDR of <0.1 was considered significant. (B) Correlation analysis between clinical and laboratorial characteristics with normalized gene expression. Only genes associated with TNF-α and IFN-γ were included. Red color indicates positive correlation and blue, negative correlation. *p-value. (C) Protein–protein interactions network reconstructed with upregulated genes. Four highly connected modules are represented. Colors indicate each module (cluster), and the dot size indicates the page rank score of genes belonging to these modules. (D) Heatmap showing the expression profiles of genes belonging to modules identified in panel (C). Unsupervised clustering of sample shows a stratification of pJIA and controls in three different clusters. Hierarchical clustering of samples was performed based on the Euclidean distance calculated from the normalized and scaled expression. Numbers (1, 2, and 3) on top of the dendrogram indicate the clusters. GSEA, gene set enrichment analysis; FDR, false discovery rate; NES, normalized enrichment score.

Figure 3

Unsupervised learning and clustering unveil the heterogeneity in pJIA. (A) Optimal number of clusters indicated by silhouette analysis. Based on average silhouette score, the vertical dashed line shows two clusters as the optimal cluster number for pJIA accurate stratification. (B) Forest plot showing the differential pattern of cytokines and their receptors in the identified pJIA clusters. Log2 fold change of cluster 2 vs. cluster 1 and standard errors are plotted. (C) Boxplot showing varying activation levels of inflammatory pathways between control group (green), cluster 1 (blue), and cluster 2 (red). Points indicate the single sample GSEA score of each sample. p-values from Kruskal–Wallis test are shown. (D) Boxplot showing the comparison of clinical characteristics between cluster 1 and cluster 2. Mann–Whitney test was performed. (E) Heatmap revealing the pattern of expression of informative features selected with recursive feature elimination approach using SVM with a linear kernel. A total of 33 genes were identified, of which four were upregulated and 29 downregulated in cluster 2. GSEA, gene set enrichment analysis; SVM, support vector machine.