Identification of Key Genes and Pathways in Myeloma side population cells by Bioinformatics Analysis

Abstract

Background: Multiple myeloma (MM) is the second most common hematological malignancy, which is still incurable and relapses inevitably, highlighting further understanding of the possible mechanisms. Side population (SP) cells are a group of enriched progenitor cells showing stem-like phenotypes with a distinct low-staining pattern with Hoechst 33342. Compared to main population (MP) cells, the underlying molecular characteristics of SP cells remain largely unclear. This bioinformatics analysis aimed to identify key genes and pathways in myeloma SP cells to provide novel biomarkers, predict MM prognosis and advance potential therapeutic targets. Methods: The gene expression profile GSE109651 was obtained from Gene Expression Omnibus database, and then differentially expressed genes (DEGs) with P-value <0.05 and |log2 fold-change (FC)| > 2 were selected by the comparison of myeloma light-chain (LC) restricted SP (LC/SP) cells and MP CD138⁺ cells. Subsequently, gene ontology (GO) and Kyoto encyclopedia of genes and genomes (KEGG) pathway enrichment analysis, protein-protein interaction (PPI) network analysis were performed to identify the functional enrichment analysis of the DEGs and screen hub genes. Cox proportional hazards regression was used to select the potential prognostic DEGs in training dataset (GSE2658). The prognostic value of the potential prognostic genes was evaluated by Kaplan-Meier curve and validated in another external dataset (MMRF-CoMMpass cohort from TCGA). Results: Altogether, 403 up-regulated and 393 down-regulated DEGs were identified. GO analysis showed that the up-regulated DEGs were significantly enriched in innate immune response, inflammatory response, plasma membrane and integral component of membrane, while the down-regulated DEGs were mainly involved in protoporphyrinogen IX and heme biosynthetic process, hemoglobin complex and erythrocyte differentiation. KEGG pathway analysis suggested that the DEGs were significantly enriched in osteoclast differentiation, porphyrin and chlorophyll metabolism and cytokine-cytokine receptor interaction. The top 10 hub genes, identified by the plug-in cytoHubba of the Cytoscape software using maximal clique centrality (MCC) algorithm, were ITGAM, MMP9, ITGB2, FPR2, C3AR1, CXCL1, CYBB, LILRB2, HP and FCER1G. Modules and corresponding GO enrichment analysis indicated that myeloma LC/SP cells were significantly associated with immune system, immune response and cell cycle. The predictive value of the prognostic model including TFF3, EPDR1, MACROD1, ARHGEF12, AMMECR1, NFATC2, HES6, PLEK2 and SNCA was identified, and validated in another external dataset (MMRF-CoMMpass cohort from TCGA). Conclusions: In conclusion, this study provides reliable molecular biomarkers for screening, prognosis, as well as novel therapeutic targets for myeloma LC/SP cells.

Keywords: Bioinformatics analysis; Cancer stem cell; Multiple myeloma; Side population cells; differentially expressed gene.

Figure 1

Volcano plot of DEGs (393 down-regulated genes and 403 up-regulated genes).

Figure 2

Heat map of the top 100 DEGs (52 down-regulated genes and 48 up-regulated genes).

Figure 3

Bubble plots of GO enrichment analysis of DEGs. (A) Bubble plot of GO enrichment analysis of up-regulated DEGs. (B) Bubble plot of GO enrichment analysis of down-regulated DEGs.

Figure 4

KEGG pathway analysis of DEGs.

Figure 5

Distribution of DEGs in myeloma LC/SP cells for the top 5 KEGG enriched pathways.

Figure 6

The top 3 modules with relatively high scores from the protein-protein interaction network. Red: up-regulation; Blue: down-regulation. (A) Module 1 with 53 nodes and 684 edges was significantly enriched in defense response, immune system process and immune response. (B) Module 2 with 21 nodes and 197 edges was significantly enriched in cell cycle phase, cell cycle and M phase of mitotic cell cycle. (C) Module 3 with 13 nodes and 60 edges was significantly enriched in G-protein coupled receptor protein, signaling pathway and chemotaxis.

Figure 7

Kaplan-Meier survival analysis of 9 prognostic genes in MM patients in the training and validation datasets. (P-value < 0.0001 in both GSE2658 and MMRF-COMMPASS). (A) Kaplan-Meier survival analysis of 9 prognostic genes in MM patients in GSE2658. (B) Kaplan-Meier survival analysis of 9 prognostic genes in MM patients in MMRF-COMMPASS.