Identification of key biomarkers associated with cell adhesion in multiple myeloma by integrated bioinformatics analysis

Background Multiple Myeloma (MM) is a hematologic malignant disease whose underlying molecular mechanism has not yet fully understood. Generally, cell adhesion plays an important role in MM progression. In our work, we intended to identify key genes involved in cell adhesion in MM. Methods First, we identified differentially expressed genes (DEGs) from the mRNA expression profiles of GSE6477 dataset using GEO2R with cut-off criterion of p < 0.05 and [logFC] ≥ 1. Then, GO and KEGG analysis were performed to explore the main function of DEGs. Moreover, we screened hub genes from the protein–protein interaction (PPI) network analysis and evaluated their prognostic and diagnostic values by the PrognoScan database and ROC curves. Additionally, a comprehensive analysis including clinical correlation analysis, GSEA and transcription factor (TF) prediction, pan-cancer analysis of candidate genes was performed using both clinical data and mRNA expression data. Results First of all, 1383 DEGs were identified. Functional and pathway enrichment analysis suggested that many DEGs were enriched in cell adhesion. 180 overlapped genes were screened out between the DEGs and genes in GO terms of cell adhesion. Furthermore, 12 genes were identified as hub genes based on a PPI network analysis. ROC curve analysis demonstrated that ITGAM, ITGB2, ITGA5, ITGB5, CDH1, IL4, ITGA9, and LAMB1 were valuable biomarkers for the diagnosis of MM. Further study demonstrated that ITGA9 and LAMB1 revealed prognostic values and clinical correlation in MM patients. GSEA and transcription factor (TF) prediction suggested that MYC may bind to ITGA9 and repress its expression and HIF-1 may bind to LAMB1 to promote its expression in MM. Additionally, pan-cancer analysis showed abnormal expression and clinical outcome associations of LAMB1 and ITGA9 in multiple cancers. Conclusion In conclusion, ITGA9 and LAMB1 were identified as potent biomarkers associated with cell adhesion in MM.

immune cells, etc. [4]. MM cells and the BM microenvironment interacted through a complex network mediated by cytokines and adhesion molecules, including the integrins, the immunoglobulin superfamily, the selectins, the cadherins, and the proteoglycans. The cross-linking activated several known regulatory pathways, which involved in the survival, proliferation, migration, homing, as well as drug resistance of MM cells, playing a crucial role in MM development [5].
Integrins are heterodimeric membrane glycoproteins expressing on the surface of many types of cells, serving as CAMs and the major receptors for the ECM [5,6]. Their expression can vary considerably between normal and cancer tissues. Studies showed that integrin expression levels were correlated with pathological outcomes, including patient survival and metastasis in many cancer types. Furthermore, Integrins regulated a series of cellular functions such as cell proliferation, migration, invasion, and survival which play a crucial role in the cancer progression [7,8]. ITGA9 is one of the least studied integrins among the 24-member integrin family in human. Furthermore, studies have suggested that ITGA9 had the abnormal expression in numerous cancers and has been found to be crucial for a number of biological processes in many types of cancers, such as breast carcinoma, melanoma and lung cancer [9][10][11][12].
Laminins are extracellular heterotrimeric basement membrane glycoproteins [13]. They are composed of three polypeptide chains, named as α, β, and γ [14]. Laminins are involved in diverse physiological and pathological processes, including involvement in basement membrane assembly, neurite outgrowth, and promotion of cell adhesion, migration, protease activity, proliferation and angiogenesis in cancers [15,16]. The link between cancer cells and laminins is vital in tumor invasion and metastasis. Invading tumor cells interact with laminins to acquire more metastatic potential [17]. LAMB1 is expressed in most tissues and is one of the 3 chains constituting laminin 1. Studies have shown that LAMB1 was shown to be a potential biomarker in some cancers [18,19].
However, there is rare research on the role which ITGA9 and LAMB1 play in myeloma development and metastasis. With the development of gene expression profiling array and second-generation high-throughput sequencing technology, people analyze the data on gene expression profiles to screen DEGs to find candidate biomarkers and potential therapeutic targets in MM [20,21]. In this study, we aim to find potential biomarkers in MM. First, we analyzed a microarray profile to obtain DEGs between active MM and normal plasma cells. Then, 12 genes were identified as hub genes based on a PPI network analysis. Finally, through integrative bioinformatic analysis, we identified that ITGA9 and LAMB1 were correlated with cell adhesion in MM.

Microarray data information
Three gene expression profiles (GSE6477, GSE2658 and GSE136324) were from the GEO database (http://www. ncbi.nlm.nih.gov/geo). The GSE6477 dataset, which included samples from 147 patients with different stages of plasma cell neoplasm and 15 NDs, was performed on the Affymetrix Human Genome U133A Array platform. The GSE2658 dataset was composed of 559 samples of MM patients. 426 patients with MM was obtained from the GSE136324 dataset. Both of them were performed on the Affymetrix Human Genome U133 Plus 2.0 Array. All the above-mentioned profiling datasets were accessible online with no conflict.

Identification of DEGs in MM
The online analysis tool GEO2R was used to screen DEGs of the GSE6477 dataset. In this study, we set "p < 0.05 and [logFC] ≥ 1" as the cut-off criterion. A heatmap of the top 100 significantly changed DEGs were drawn using the Multiple Experiment Viewer software.

GO and KEGG enrichment analysis
GO enrichment analysis (http://www.geneo ntolo gy.org/) and KEGG pathway enrichment analysis (http://www. genom e.jp/kegg/pathw ay.html) were used to illuminate the biological functions and pathways associated with DEGs [16], which were both integrated in the DAVID (http://david .abcc.ncifc rf.gov/) program. Enriched GO and KEGG terms with p < 0.05 were considered as statistical significance. The ClueGO and CluePedia plug-in Cytoscape software version 3.7.1 (http://www.cytos cape. org/) were used to analyze the pathways interrelation.

PPI network construction
STRING (http://string-db.org) online database is used to predict the PPI information. The cut-off value for the filtration criteria was set at the confidence score > 0.7. Then, the PPI network was visualized by Cytoscape. The MCODE plugin in Cytoscape was applied to select modules of the PPI network.

ROC analysis
The ROC curves were used to explore the sensitivity and specificity of DEGs for MM diagnosis using GraphPad Prism Software (Version V, La Jolla, CA, USA).

Survival analysis
The PrognoScan database (http://dna00 .bio.kyute ch.ac. jp/Progn oScan /) was used to investigate the prognostic impact of selected genes in different cancers. According to the median expression of a particular gene, the patients were split into high and low expression groups. The OS and DSS of GC patients was evaluated using a KM plot. A Cox p value<0.05 was considered as statistically significant.

Analysis of expression in various cancers
The mRNA expression levels of LAMB1 and ITGA9 in various cancers and their normal tissue counterparts were analyzed using the Oncomine database (https :// www.oncom ine.org/resou rce/login .html). The relationship between two specific gene mRNA expression levels in various cancers were analyzed using the GEPIA (Gene Expression Profiling Interactive Analysis) database (http://gepia .cance r-pku.cn/).
The default weighted enrichment method was applied for enrichment analysis. The random combination was set for 1000 times. The analysis was performed with following settings: FDR<0.25, NOM p-value<0.05 and |NES|>1. The enrichment map plugin in Cytoscape was applied to visualize geneset enrichment results.

Statistical analysis
An independent sample t-test and one-way ANOVA were used to compare the statistical significance between two or more samples, respectively. Two-tailed p-values < 0.05 were considered as statistically significance. The correlations analysis was assessed by the Spearman correlation coefficient and Chi square test. Statistical analysis was performed by the SPSS 22.0 and GraphPad Prism Software.

Identification of DEGs
We acquired mRNA expression profiles of plasma cell samples from 101 patients of new or relapse multiple myeloma and 15 NDs from the GSE6477. Gene expression distribution of profiles was matched via boxplot analysis (Fig. 1a). Statistical analysis was required to validate the comparison between NDs and MM patients. We applied the GEO2R to identify DEGs associated with MM. 1383 DEGs (538 upregulated genes and 845 downregulated genes) were identified in the MM samples compared to NDs' samples. We set the cut-off criterion as p < 0.05 and [logFC] ≥ 1 (Fig. 1b). The top 100 genes were clustered in the heatmap between patients with MM patients and NDs (Fig. 1c).

Functional and pathway enrichment analysis
We performed GO and KEGG enrichment analysis to investigate the functions of DEGs using DAVID. The top GO (Fig. 2a-c, Table 1) and KEGG (Fig. 2d, Table 2) terms for DEGs were shown. For biological process (BP), DEGs were mainly enriched in immune response, inflammatory response, ECM organization, leukocyte migration, cell adhesion (Fig. 2a). DEGs in molecular function (MF) were significantly associated with protein binding, ECM structural constituent, serine type endopeptidase activity, protease binding and receptor activity (Fig. 2b). The cellular components (CC) analysis indicated that proteins encoded by DEGs were mostly located in the extracellular exosome, extracellular space, extracellular region and ECM (Fig. 2c). KEGG enrichment analysis showed that CAMs, proteoglycans in cancer, ECM-receptor interaction and PI3K/Akt signaling pathway were significantly enriched in DEGs. The top KEGG pathways for DEGs were shown (Fig. 2d). Subsequently, we analyzed a pathway interrelation and related genes by examining KEGG enrichment analysis results in ClueGO and CluePedia (Fig. 2e). These results suggested that cell adhesion might play a crucial part in the malignant progression of MM.

PPI network and modular analysis of selected genes
To identify key genes in the cell adhesion processes, 180 overlapped genes in DEGs and genes involved with the cell adhesion were screened, including 123 mRNAs and 57 non-coding RNAs (Fig. 3a). Then, we constructed a PPI network containing 123 protein, with 123 nodes and 429 edges based on the STRING database and Cytoscape (Fig. 3c). 12 hub genes whose degree values ≥ 15 for further analysis were chosen, including ITGB1, FN1, ITGB3, ITGAM, PTPRC, ITGB2, ITGA5, ITGB5, CDH1, IL4, ITGA9 and LAMB1 (Fig. 3b). Additionally, Cytoscape displayed a total of eight modules in the default MCODE settings for modular analysis ( Fig. 3d-k, Table 3).

Diagnostic and prognostic values of hub genes in MM
First of all, we performed the ROC curve analysis among 12 hub genes based on the GSE2658. The results showed that ITGAM, ITGB2, ITGA5, ITGB5, CDH1, IL4, ITGA9, and LAMB1 achieved an AUC value of > 0.7, demonstrating that these eight genes have high sensitivity and specificity for MM, suggesting they can be served as biomarkers for the diagnosis of MM (Fig. 4a-h).
Furthermore, we evaluated the influence of all aforementioned 8 hub genes on clinical prognosis using Pro-gnoScan database based on the GSE2658 through K-M curve and log-rank test to identify whether these genes were concerned with the survival of MM patients. ITGA9 ( Fig. 4i) and LAMB1 (Fig. 4j) were significantly associated with DSS in MM by evaluating the correlation between these gene expressions and survival rates.

The relationship between ITGA9/LAMB1 mRNA expression and clinical characters of patients with MM
First of all, the relationship between ITGA9 or LAMB1 mRNA expression and the clinical parameters was analyzed based on both clinical data and mRNA expression data of GSE136324 (Table 4, n = 426). The results showed that LAMB1 expression was negatively correlated with β 2 -MG (Fig. 5a, p = 0.024, spearman correlation coefficient r = -0.110). Patients with low LAMB1 expression had a tendency of being in high R-ISS stage (Fig. 5c, p = 0.001). No significant relationship was found between LAMB1 mRNA expression and malignant plasma cells (PCs) in BM (Fig. 5b, p = 0.9427) and GEP groups (Fig. 5d, p = 0.1785).

Function and signaling pathways analysis of ITGA9 and LAMB1
To investigate the function of ITGA9 and LAMB1 on MM progression, we performed a comprehensive analysis including GSEA, TF prediction and Spearman correlation analysis. Firstly, an enrichment map was constructed using genesets which related to cell adhesion with p-value < 0.05 (Fig. 6a). Genesets including GO calcium dependent cell cell adhesion via plasma membrane cell adhesion molecules, GO cell adhesion via plasma membrane adhesion molecules, KEGG ECM receptor interaction, NABA basement membranes, NABA core matrisome were enriched in patients with low LAMB1 expression in the GSE2658, suggesting that LAMB1 was involved in the cell adhesion in myeloma cells (Fig. 6b-f ). Furthermore, cancer microenvironment-dn, MANALO hypoxia dn, ELVIDGE hypoxia dn genesets were significantly enriched in patients with low LAMB1 expression in the GSE2658 (Fig. 6g-i). For futher study, we found the positive correlation between LAMB1 and HIF-1 mRNA expression in a varity of cancers using GEPIA database (Fig. 6j, p = 0.00018, spearman correlation coefficient r = 0.038). The result of TF prediction by the online site PROMO (http://algge n.lsi.upc.es/cgi-bin/promo _v3/ promo /promo init.cgi?dirDB =TF_8.3/) strengthened the correlation that HIF-1 may bind to the promoter of LAMB1 (Fig. 6k). Therefore, LAMB1 which may be regulated by HIF-1 played a vital role in myeloma cell adhesion.
ITGA9 also was a cricial cell adhesion molecule since GSEA results showed that genesets including KEGG ECM receptor interaction, GO cell adhesion via plasma membrane adhesion molecules, NABA ECM glycoproteins, GO calcium dependent cell cell adhesion via plasma membrane cell adhesion molecules, GO protein complex involved in cell adhesion, NABA core matrisome, cell-cell adhesion were enriched in patients with low ITGA9 expression in the GSE2658 (Fig. 7b-h). An enrichment map was constructed using genesets which related to cell adhesion with p-value < 0.05 (Fig. 7a). Results demonstrated that SCHLOSSER myc targets repressed by serum and hallmark myc targets v1 genesets were enriched in patients with low ITGA9 expression in the GSE2658 (Fig. 7i-j). There was a negative correlation between the mRNA expression of ITGA9 and MYC in multiple cancers based on the GEPIA database (Fig. 7k, p = 7.7e-19, spearman correlation coefficient r = −0.09). We predicted that MYC may bind to the promoter of ITGA9 and repressed its expression based on TF prediction by PROMO (Fig. 7l). In conclusion, ITGA9 was a vital cell adhesion molecule in myeloma which may be negtive regulated by MYC.

ITGA9 and LAMB1 expression and clinical outcome association in different cancers
We explored the expression and prognosis of ITGA9 and LAMB1 mRNA in different cancers based on the Oncomine database and the PrognoScan database respectively.
LAMB1 is highly expressed in most cancers while low in some cancers including breast cancer, leukemia, ovarian cancer and prostate cancer (Fig. 8a). Furthermore, LAMB1 was significantly associated with OS and DSS respectively in AML, DLCBL, glioma, NSCLC, colorectal cancer (Fig. 8b-g) and colorectal cancer, breast cancer ( Fig. 8h-i).
Compared with ND, the transcript levels of ITGA9 indicated significant low expression in bladder cancer, brain and central nervous system cancer, breast cancer, leukemia, liver cancer, lung cancer and etc., suggesting that the down-regulation of ITGA9 was common in various types of cancer (Fig. 9a). ITGA9 was significantly associated with OS and DSS respectively in AML, B cell lymphoma, lung cancer, colorectal cancer and esophagus cancer (Fig. 9b-f ) and colorectal cancer, breast cancer ( Fig. 9g-h). The details including hazard ratio (HR) with 95% confidence intervals(CI) and p-values were shown in the Table 5.

Discussion
MM is a common hematological malignancy so it is vital to investigate the molecular mechanisms. Microarray has been widely used to analyze the expression changes of genes in MM and predict the potential biomarkers.
Here, we found that LAMB1 was abnormally expressed and was associated with OS and DSS in many cancers, which were consistent with previous results. LAMB1 has a high protein level in high-grade gliomas, suggesting a possible correlation with tumor progression [24]. What's more, LAMB1 was identified to take part in cell attachment and have the capacity to inhibit metastasis. In prostate cancer, LAMB1 was shown to be involved in cell motility and invasion into the surrounding ECM [19]. In our study, low LAMB1 expression was significantly associated with high β2-MG concentration and high R-ISS stage in MM. GSEA results showed that LAMB1 was involoved in cell adhesion and may be induced by hypoxia. TF prediction and the correlation analysis suggested that HIF-1 may bind to LAMB1 promoter to increase its transcription. ITGA9 abnormal expression was found in many cancers and was likely to correlate with higher grade cancers [9]. For example, sequencing analysis also found that ITGA9 was significantly down-regulated in cervical squamous cell carcinoma [22]. Genetic variation and epigenetic modification of ITGA9 are related with the tumorgenicity and progression of colorectal cancer [23]. In our study, we found down-regulated ITGA9 was associated with poor outcome in MM. Myeloma patients with low ITGA9 expression had more tendency of having the higher number of malignant PCs in BM and becoming the higher GEP group. GSEA results showed that ITGA9 was a vital cell adhesion molecule and may be repressed by MYC in myeloma. Based on the TF prediction and the correlation analysis, we predicted MYC may bind to the promoter of ITGA9 for transcriptional repression.

Conclusions
We identified eight hub genes, including ITGAM, ITGB2, ITGA5, ITGB5, CDH1, IL4, ITGA9, and LAMB1 to be potential diagnostic markers in MM. Further study demonstrated ITGA9 and LAMB1 which correlated with clinical characters and prognosis may play important roles in the cell adhesion and can be regulated by different TFs in MM.