1. Infiltrating immunocytes landscape in periodontitis
To explore the differences of infiltrating immunocytes between healthy and periodontitis tissues, CIBERSORT algorithm was conducted to reveal the distribution of 22 immunocyte fractions and we found no matter in healthy or periodontitis tissue the plasma cell was the most in content (Figure 1A). We compared the fractions of 22 immunocytes between healthy and periodontitis, and 13 immunocytes content were altered in periodontitis tissue including 5 immunocytes fractions increased (plasma cell, neutrophils, T cells CD4 memory activated, dendritic cells activated and NK cells resting) and 8 immunocytes fractions decreased (Tregs, dendritic cells resting, eosinophils, mast cells resting, macrophages M1, B cells memory, T cells follicular helper and T cells CD4 memory resting) in periodontitis tissue (Figure 1B). Among them, Tregs and neutrophils changed most in decreasing and increasing (Figure 1C). Then, to confirm the variation in immunocytes proportions might be an intrinsic feature that may characterize the individual diversities, we performed PCA analysis and found that the proportions of immune cells displayed distinct clustering and individual differences (Figure 1D). Among them, plasma cells, B cells memory and mast cells resting were the most contributed in the healthy and periodontitis tissue variation. Correlation analysis was performed in periodontitis samples and all samples. We found in both periodontitis and all samples that B cells memory and plasma cells had a big and significantly negative correlation, so as to plasma cells and macrophages M1 (Figure 1E). Further, the differences of infiltrating immunocytes in different clinical subgroups were also investigated including periodontitis subtype, gender, and age by taking similar analysis strategies (Figure S1).
To explore and identify the robustness and contributing weight of immunocytes that are related to periodontitis independently, we conducted a series of regression algorithms to adjust the immunocyte interaction effect including univariate logistic regression, LASSO regression, and multivariate logistic regression. 12 periodontitis related immunocytes were first identified by univariate analysis (Figure 2A). Then, LASSO regression was performed to make feature selection and dimension reduction so that the unimportant immunocytes could be excluded. Therefore, we screen out 11 important periodontitis related immunocytes through LASSO regression (Figure 2B-C). Last, multivariate logistic regression was conducted on the 11 immunocytes to get the independently periodontitis related immunocytes. Only the dendritic cells resting, dendritic cells activated and T cells follicular helper were independently related to periodontitis (Figure 2D).
2. Identification of immune-related periodontitis subtypes and their characteristics
To identify immune-related periodontitis subtypes, 1214 immune-related genes were compiled, and 702 of them were significantly related to periodontitis (FDR < 0.05). Using the 702-gene panel, 2 distinct immune-related periodontitis subtypes were identified by performing consensus clustering analysis (Figure 3A, Figure S2A-C). Then, 116 genes were selected from all 702 immune-related genes and the sparse structure of the metagenes appears in the localized patterns of the gene contributions, in which the 116 genes were defined as hub genes for each subtype representing the subtype characteristics (Figure 3B, Figure S2D). There were 45 immune-related genes represent subtype 1 (cluster 1) and 71 genes for subtype 2 (cluster 2). To get a further understanding of the regulation effects among these genes, the PPI networks of the 116 hub genes revealed their relationships and importance in the network., in which CXCR4 of cluster-1 has the most links with other genes and PIK3RI of cluster-2 has the most (Figure 3C).
To understand the characteristics and differences between the two immune-related subtypes, GO enrichment analysis and canonical pathway enrichment analysis was performed. Biological processes of GO were conducted and to simplify the results, we grouped categories according to functional theme, which represents the biological characteristics of the cluster (Figure 4A). The characteristics of cluster-1 mainly concerning FC-receptor and immune-associated receptor-mediated immune regulations such as FC/FC-gamma receptor signaling pathway and cell surface receptors. Another characteristic of cluster-1 is about B cells immune processes regulations such as B cells differentiation and B cell activation. As for cluster-2, the IL-6 cytokine feature is predominant such as regulation of IL-production. Then, canonical pathway analysis revealed the significantly enriched or activated pathways in cluster-1 and cluster-2 (Figure 4B-C). The results were similar as above and we found some famous pathways were shown such as JAK-STAT signaling in cluster-1 and PK3K-AKT in cluster-2.
3. Associations between immune subtypes and immunocytes
To investigated the diversities between the two immune subtypes from cell level, the 22 infiltrating immunocyte fractions were compared and clinical features were also compared. From the results, the two immune subtypes had a large difference in infiltrating immunocytes, in which 13 immunocytes were different (Figure 5A-B). We found cluster-1 had more dendritic cells activated, plasma cells and neutrophils, while B cells memory, T cells CD4 memory resting, T cells follicular helper, NK cells activated, monocytes, macrophages M1, macrophages M2, dendritic cells resting, mast cells activated and mast cells resting were enriched in cluster-2. Besides, we noticed cluster-2 was significantly associated with chronic periodontitis and cluster-1 associated with aggressive periodontitis (Figure 5C), while there is no different distribution between two clusters about gender and age.
4. Differentially expressed genes and gene-sets between healthy and periodontitis
To understand the molecular mechanisms by which genes or gene-sets drive the biological reactions in the periodontitis, DEGs analysis and gene-set enrichment analysis were conducted. Using the criterion of FDR < 0.01 and |logFC| > 1, 171 DEGs were identified and 134 of them were up-regulated in periodontitis while 37 of them were down-regulated. These up-regulated genes may be involved in the progress of periodontitis, and among them, the MZB1, TNFRSF17, and IGLL5 are the top three up-regulated genes. The interaction relationships among these genes were also investigated and presented as a PPI network (Figure S3A). According to outward traversal from a locally dense seed protein and vertex weighting by local neighborhood density, these DEGs were divided into several molecular complex and presented by different shapes. Besides, the top ten hub genes among these genes were also identified with the highest links to nodes, indicating they are in the core position of the network and play an essential role in the regulation of the network.
Not only the individual gene was investigated, but also gene-sets were. We mainly focus on four aspects of gene-sets including immunologic signatures, hallmark gene sets, GO biological processes and canonical pathways. Results showed many enriched gene-sets in periodontitis were associated with immune processes or regulations (Figure S3B-E) such as CD4 T CELL VS B CELL DOWN in immunologic signature, IL6 JAK STAT6 SIGNALING in hallmark, POSITIVE REGULATION OF B CELL DIFFERENTIATION in GO and T HEPLER PATHWAY in the canonical pathway. All these four aspects implied immune-related response is essential in periodontitis.
5 Identification of infiltrating immunocytes related hallmark pathways
To further understand the molecular mechanisms by which hallmark pathways are involved in immunocyte content, we examined the correlation between the fractions of 22 immunocytes and the activity of 50 hallmark pathways. In total, 463 statistically significant correlations between 22 immunocytes and 50 hallmark pathways were observed, in which the T cells CD4 naive and B cells memory had the most negative correlated pathways as 24 while Plasma cells had the most positive correlated pathways of 21 (Figure 6A). The links with high correlation coefficients were presented (Figure 6B) such as dendritic cells resting is negatively correlated with HALLMARK COAGULATION with -0.56 correlation and plasma cells is positively correlated with HALLMARK UNFOLDED PROTEIN RESPONSE with 0.57 correlation.
6. Infiltrating immunocytes related gene module identified by WGCNA
To understand the molecular mechanisms of alteration in infiltrating immunocytes in detail, the genes or gene modules regulating immunocytes or affected by immunocytes need to be uncovered. Therefore, WGCNA, an advanced bioinformatic algorithm, was conducted to find clusters (modules) of highly correlated genes related to external sample characteristics (immunocytes fractions). The analysis of network topology for thresholding powers was calculated from 1 to 20 and identified 13 was used as network topology for thresholding power considering its relatively mean connectivity power and balanced scale independence (Figure S4A-B). Clustering of module eigengenes was performed in the top 25% variance 4329 genes and 27 modules were identified (Figure 7A-B). Then, the immunocyte fractions and clinical features related gene modules were identified presenting in a heatmap (Figure 7C). For example, plasma cells, we can find the turquoise module was the most positive associated gene module to plasma cells fractions from the heat map and scatter plot (Figure 8A), suggesting the genes in turquoise module is involved in regulation for plasma cell population or affected by increased plasma cell fractions. The genes in the turquoise module with high gene significance and module membership were presented as PPI network so that we could locate the central genes according to links degree (Figure S4C), and we noticed BTK is the hub gene with the most links in the network. GO analysis was also conducted to assess the biological function of the turquoise module and it involved in B cell regulation (Figure 8B). Furthermore, plasma cell-related genes were also investigated (Figure 8C-E, Figure S4D).
7. An easy-to-use and user-friendly web-tool for periodontitis related immunity analysis
To help periodontitis researchers used the results described and facilitate broad access to results in this work to any aspect of interest about periodontitis, we established an easy-to-use and user-friendly interactive web-tool for researchers working on periodontitis (http://118.24.100.193:3838/tool-PIA/). We named this platform as Periodontitis Immune Atlas (PIA). PIA provides 6 modules for users: “Genes”, “Pathways”, “WGCNA”, “Correlations”, “Immunocytes” and “Functions”, and all the results in these six modules can be downloaded for usage. The “Genes” module (Figure 9A) provides the DEGs between healthy and periodontitis in several forms. Users can browse the DEG list in a table and sort them by the way you want. Besides a gene list of immune-related genes with their category is also presented. The specific gene of interest can be visualized by volcano-plot or box-plot. Two genes’ correlation analysis are also provided. The “Pathway” module (Figure 9B) allows users to inquire about the activity score of specific pathways and presented them by box-plot. Three commonly used pathway sets are involved including hallmarks, KEGG and GO-BP. The correlation between the two pathways can also be investigated. The “WGCNA” module (Figure 9C) provides analysis for phenotype and gene module. Users can find gene modules related to interesting phenotypes and the gene modules function can also be investigated. In the “Correlation” module (Figure 9D), users can fully explore the relationships among genes, pathways, and immunocytes. “Immunocytes” module (Figure 9E) allows users to visualize the status of immunocytes in different clinical status and the correlation between two immunocytes is also included. The “Functions” module (Figure 9F) provides GSEA analysis between different clinical status. Besides, single-gene GSEA can be used as a function prediction for a specific gene in periodontitis. All figures generated in this platform can be downloaded for further analysis and the query results are demonstrated in a comfortable form. The functions provided in the web-tool, which will be continuously updated, may serve as a guide for researchers interested in periodontitis related immunity