3.1 Screening of DEGs between the T1DM samples and normal samples
The 216 DEGs between the T1DM samples and normal samples were screened out (Figure 1a), including 92 up-regulated genes (Supplementary Table 1) and 124 down-regulated genes (Supplementary Table 2). The distribution of the DEGs was displayed by the volcano plot (Figure 1b), and the expressions of the DEGs in each sample were shown in Figure 1c.
Fig. 1 The amount of 216 DEGs between the T1DM samples and normal samples. a Histogram; b volcano plot; c heatmap
3.2 The raw identification and functional enrichment analysis of the immune-related DEGs
The 1639 immune-related genes were overlapped with the DEGs using the Venn diagram, getting 50 immune-related DEGs (Figure 2a) that were used for subsequent analysis.
To further investigate the functions of the immune-related DEGs, the GO and KEGG analyses were performed. The GO terms with the top 5 of gene enrichment were protein binding, membrane, extracellular region, extracellular space, and integral component of membrane (Figure 2b and Supplementary Table 3), suggesting that the immune-related DEGs may be involved in the biological activities of the cell membrane. Meanwhile, the analysis of KEGG pathway indicated that these genes were mainly enriched in Type 1 diabetes mellitus pathway (Figure 2c and Supplementary Table 4), further verifying that the immune-related DEGs played a vital role in the occurrence and development of T1DM.
Fig. 2 Raw identification and functional enrichment analysis of the immune-related DEGs. a Venn diagram of overlapping genes; b GO terms; c KEGG pathway
3.3 Identification of the optimal immune-related biomarkers
To identify immune-related hub genes, the LASSO regression analysis for the 50 immune-related DEGs was performed to screen the gene signatures, getting 11 gene signatures (Figure 3a and 3b). Besides, the SVM for the 50 immune-related DEGs was also used to screen gene signatures, getting 6 gene signatures (Figure 3c). Subsequently, the 11 gene signatures identified by LASSO were overlapped with the 6 gene signatures identified by the SVM, and ultimately obtaining 5 hub genes [C-C motif chemokine receptor 3 (CCR3), major histocompatibility complex, class II, DQ alpha 1 (HLA-DQA1), NCR3, toll like receptor 3 (TLR3) and tumor necrosis factor (TNF), Figure 3d], which were considered as the optimal immune-related biomarkers.
Fig. 3 Identification of the optimal immune-related biomarkers. a b LASSO regression analysis; c SVM analysis; d venn diagram of overlapping genes
3.4 Distribution of immune cells between T1DM and normal samples
To further investigate the correlation of the hub genes related to immune with immune cells, the distribution of the 22 immune cells between T1DM and normal samples was analyzed by the CIBERSORT algorithm. Figure 4a showed the percentage of immune cells in each sample, indicating that Monocytes, resting NK cells, resting memory CD4+ T cells, and CD8+ T cells had a larger proportion (Supplementary Table 5). Removing the four types of immune cells that were not in the sample, the proportion of the other immune cells in each sample was displayed by a heatmap (Figure 4b). According to statistical results, the distribution of plasma cells (P = 0.05), resting mast cells (P = 0.013) and neutrophils (P = 0.018) in T1DM samples was significantly increased, while the distribution of resting NK cells (P = 0.009) in T1DM samples was significantly reduced, compared with that in normal samples (Figure 4c). The above results suggested that these differential immune cells might be involved in the immune regulation process of T1DM pathogenesis.
Fig. 4 Distribution of immune cells between T1DM and normal samples. a Percentage of immune cells in each sample; b heatmap; c vioplot
3.5 Correlation analysis of the immune cells and the hub genes
To further analyze the correlation of the hub genes and the immune cells, the Spearman correlation heatmaps were plotted respectively in T1DM and normal samples. As shown in Figure 5a and Figure 5b, the differential immune cells were marked with red. NCR3 was only significantly related to the activated NK cells (P < 0.05) in the normal samples, while was significantly associated with multiple immune cells including M0 macrophages, monocytes, resting NK cells, and resting memory CD4+ T cells (all P < 0.05) in the T1DM samples. Moreover, TNF had no significant correlations to all immune cells (all P > 0.05) in the normal samples, but was significantly related to naive B cell (P < 0.05) and naive CD4+ T cell (P < 0.01) in the T1DM samples. Additionally, there were significant differences between T1DM and normal samples that TLR3, HLA-DQA1 and CCR3 were only related to CD8+ T cells (P < 0.05), resting mast cells (P < 0.05) and Neutrophils (P < 0.05), respectively. These results suggested that NCR3 and TNF might play an important role in the regulation of immune cell-mediated T1DM progression.
Fig. 5 Spearman correlation of the hub genes and the immune cells. a T1DM samples; b normal samples
3.6 Construction of the hub genes and pathways network
To in-depth investigate functions of the hub genes, a network of the 50 immune-related DEGs and GO-BP interaction was constructed using ClueGO Plug-in of Cytoscape software, in which the hub genes were marked as red frame. As shown in Figure 6a, as expected, the pathways that interacted with TNF were the most, mainly cellular response to interferon-gamma-related anti-infection pathways. Next was NCR3, which was mainly related to natural killer cell mediated immunity-related innate immune response pathways. Further, Figure 6b also indicated that proportions of genes in cellular response to interferon-gamma and natural killer cell mediated immunity terms were the largest, 59.52% and 16.67% respectively. To more clearly show the interaction of TNF and NCR3 with pathways, the up-regulated genes and down-regulated genes were separated to construct networks which were shown in Figure 6c-d and Figure 6e-f. Similar results were obtained, TNF participated in antimicrobial humoral response-related anti-infection pathways, NCR3 was still involved in natural killer cell mediated immunity. Besides, the network of the 50 immune-related DEGs and KEGG interaction also indicated similar results (Supplementary Figure 1).
The above results suggested that TNF and NCR3 might regulate T1DM progression respectively by anti-infection pathways and natural killer cell mediated immunity. NK cells recognize and kill virus-infected cells in the absence of antibodies and major histocompatibility complex (MHC), allowing for a much faster immune reaction. The role of NK cells in both the innate and adaptive immune responses is becoming increasingly important in research using NK cell activity as a potential therapy (17).
Fig. 6 The network of hub genes and GO-BP interaction. a All genes; b the pie chart of all genes; c up-regulated genes; d the pie chart of up-regulated genes; e down-regulated genes; f the pie chart of down-regulated genes
3.7 Diagnostic value of TNF and NCR3
To explore the accuracy of the TNF and NCR3 as the diagnostic biomarkers for T1DM, the ROC curves were plotted, respectively. The AUC was 0.763 (TNF, Figure 7a) and 0.918 (NCR3, Figure 7b), suggesting that NCR3 possessed a higher accuracy than TNF in diagnosis of T1DM. And this might be related to NK cell mediated innate immune system which was often ignored in the design of novel immune-based therapies (18).
Fig. 7 The ROC curves. a TNF; b NCR3