Three-microRNA signature identified by bioinformatics analysis predicts prognosis of gastric cancer patients

Abstract

AIM

To identify multiple microRNAs (miRNAs) for predicting the prognosis of gastric cancer (GC) patients by bioinformatics analysis.

METHODS

The original microarray dataset GSE93415, which included 20 GC and 20 tumor adjacent normal gastric mucosal tissues, was downloaded from the Gene Expression Omnibus database and used for screening differentially expressed miRNAs (DEMs). The cut-off criteria were P < 0.05 and fold change > 2.0. In addition, we acquired the miRNA expression profiles and clinical information of 361 GC patients from The Cancer Genome Atlas database to assess the prognostic role of the DEMs. The target genes of miRNAs were predicted using TargetScan, miRDB, miRWalk, and DIANA, and then the common target genes were selected for functional enrichment analysis.

RESULTS

A total of 110 DEMs including 19 up-regulated and 91 down-regulated miRNAs were identified between 20 pairs of GC and tumor adjacent normal tissues, and the Kaplan-Meier survival analysis found that a three-miRNA signature (miR-145-3p, miR-125b-5p, and miR-99a-5p) had an obvious correlation with the survival of GC patients. Furthermore, univariate and multivariate Cox regression analyses indicated that the three-miRNA signature could be a significant prognostic marker in GC patients. The common target genes of the three miRNAs are added up to 108 and used for Gene Functional Enrichment analysis. Biological Process and Molecular Function analyses showed that the target genes are involved in cell recognition, gene silencing and nucleic acid binding, transcription factor activity, and transmembrane receptor activity. Cellular Component analysis revealed that the genes are portion of nucleus, chromatin silencing complex, and TORC1/2 complex. Biological Pathway analysis indicated that the genes participate in several cancer-related pathways, such as the focal adhesion, PI3K, and mTOR signaling pathways.

CONCLUSION

This study justified that a three-miRNA signature could play a role in predicting the survival of GC patients.

Key Words: Gene functional enrichment; Prognosis; Bioinformatic analysis; Differentially expressed miRNAs; Gastric cancer

Core tip: We identified 110 differentially expressed miRNAs through mining the datasets of Gene Expression Omnibus database and acquired the miRNA expression profiles and clinical information of 361 gastric cancer (GC) patients from The Cancer Genome Atlas database. Multiple miRNAs together acting as biomarkers may have a stronger reliability in survival prediction. Our study found that a novel three-miRNA signature could be used for predicting the prognosis of GC patients.

INTRODUCTION

Gastric cancer (GC) is the fourth most common cancer in incidence and the second in mortality among all cancers worldwide[1]. In 2008, a total of 989600 individuals were newly diagnosed with GC and 738,000 deaths occurred, therefore, this disease is a serious public health issue worldwide[2]. Research studies that explore the cellular and molecular mechanisms of GC development and the validation of novel biomarkers are urgently needed to achieve early diagnosis and treatment.

MicroRNAs (miRNAs), which are endogenous small noncoding RNAs (20-22 nt), have been identified as the key regulators of genes at the post-transcriptional level[3]. Increasing studies have found that miRNAs are associated with the development and progression of GC, and can act as important biomarkers in diagnosis[4,5], therapy[6], and prognosis[7,8]. Thus, the identification of differentially expressed miRNAs (DEMs) may contribute to the early diagnosis and the prediction of survival prognosis in GC.

Several studies have found that a number of miRNAs are differentially expressed in GC and are associated with survival prognosis. However, these studies lack a large sample size or an appropriate proportion of samples. A reliable survival prediction requires large-scale samples that include detailed clinical characteristics. The Gene Expression Omnibus (GEO) database is a public functional genomics data repository that includes array- and sequence-based data and allows users to query and download experiments or curated gene expression profiles[9]. The Cancer Genome Atlas (TCGA, http://cancergenome.nih.gov) project is one of the most useful cancer genomics programs and has generated, analyzed, and made available genomic sequence, expression, methylation, and copy number variation data on over 11000 individuals who represent over 30 different types of cancer[10]. In the present study, we identified DEMs between GC and adjacent normal tissues by analyzing the miRNA data of GSE93415 from GEO. In addition, the associations between DEMs and survival prognosis were analyzed using the expression profiles and clinical features downloaded from TCGA.

MATERIALS AND METHODS

Microarray data processing and DEMs identification

The microarray data of GSE93415 were downloaded from the GEO database (https://http://www.ncbi.nlm.nih.gov/geo/) and the miRNA expression data were processed with the limma package in R. Statistically significant DEMs between GC and adjacent normal samples were identified with the cut-off criterion P < 0.05 and fold change > 2.0.

Association analysis between DEMs and GC patients’ survival

TCGA (https://cancergenome.nih.gov/) stomach adenocarcinoma and adjacent normal tissue miRNA sequencing data and clinical information were downloaded for analysis. The inclusion criteria included: (1) samples with completed data for analysis; (2) patients had not received preoperative chemoradiation; and (3) overall survival time less than 80 mo. Consequently, 361 GC samples were included in the present study. The Kaplan-Meier method and log-rank test were conducted to test the prognostic value of DEMs. When P < 0.05, miRNAs were considered significantly associated with the prognosis of patients. Then, we ranked prognosis-related miRNAs according to the median expression level. Subsequently, we scored each GC patient in accordance with a high or low level of expression, and a risk grade was defined by the total scores. Finally, GC patients were sorted into high and low risk groups by the risk-score rank. The prognosis-related miRNA signature was used to analyze overall survival between high and low risk group patients using a Kaplan-Meier curve.

Target genes prediction of prognostic DEMs

We used four online tools to predict the potential target genes of the prognostic related DEMs, including TargetScan (http://www.targetscan.org/vert_71/), miRDB (http://www.mirdb.org/), miRWalk (http://zmf.umm.uni-heidelberg.de/apps/zmf/mirwalk2/index.html), and DIANA (http://www.microrna.gr/microT-CDS). In order to obtain the more reliable target genes, the Venn plot was performed to acquire the consensus genes of the four online tools.

Function analysis of target genes

FunRich [Functional Enrichment analysis tool (http://www.funrich.org/)] is a stand-alone software used for functional enrichment and interaction network analysis of genes and proteins[11]. Enrichment analysis was conducted on the consensus genes using the FunRich tool in the following categories: Biological Process, Cellular Component, Molecular Function, and Biological Pathways. P < 0.05 was considered statistically significant.

Statistical analysis

The data of miRNA expression in GC and adjacent normal samples were performed by unpaired t-test. The association between DEMs expression and clinical characteristics was analyzed by the chi-square and t-tests. Kaplan-Meier survival analysis and the univariate/multivariate Cox regression analysis were used to assess the expression levels of DEMs and prognostic features. All the statistical analyses were performed with IBM SPSS version 19.0 and P < 0.05 was considered statistically significant.

RESULTS

Identification of DEMs in GC

The microarray data of GSE93415, including 20 pairs of GC and adjacent normal tissue samples, were obtained from the NCBI-GEO database. After applying cut-off criteria of P < 0.05 and fold change > 2.0, a total of 110 DEMs were identified between GC and adjacent normal tissues (Table 1). The results of 19 downregulated miRNAs and 91 upregulated miRNAs are displayed in the volcano plot (Figure 1). A heat map of hierarchic cluster analysis showed that DEMs could be discriminated between GC and normal tissues (Figure 2).

Table 1 The differentially expressed miRNAs identified between gastric cancer and adjacent normal tissues.

Upregulated DEMs1	P value	Downregulated DEMs	P value
hsa-miR-199a-3p/hsa-miR-199b-3p	7.10E-05	hsa-miR-652-5p	0.000135
hsa-miR-125b-5p	2.99E-05	hsa-miR-1269b	1.13E-09
hsa-miR-199a-5p	7.65E-07	hsa-miR-665	2.86E-06
hsa-miR-223-3p	7.49E-06	hsa-miR-375	0.003304
hsa-miR-196a-5p	3.22E-07	hsa-miR-4501	1.64E-06
hsa-miR-27a-3p	0.000112	hsa-miR-4279	1.94E-10
hsa-miR-23b-3p	4.71E-05	hsa-miR-943	2.46E-05
hsa-miR-21-5p	1.36E-05	hsa-miR-148a-3p	1.53E-05
hsa-miR-100-5p	0.000197	hsa-miR-1275	0.000349
hsa-miR-20a-5p	0.000102	hsa-miR-4290	8.91E-11
hsa-miR-23a-3p	5.15E-10	hsa-miR-4268	9.16E-11
hsa-miR-1	0.002694	hsa-miR-891a	7.29E-10
hsa-miR-214-3p	1.23E-08	hsa-miR-4795-3p	3.07E-05
hsa-miR-10a-5p	2.33E-05	hsa-miR-1298	3.00E-06
hsa-miR-135b-5p	1.04E-05	hsa-miR-660-3p	1.75E-08
hsa-miR-99a-5p	0.001109	hsa-miR-4661-5p	9.22E-08
hsa-miR-20b-5p	0.000365	hsa-miR-4539	2.37E-09
hsa-miR-199b-5p	0.000291	hsa-let-7d-3p	1.10E-08
hsa-miR-10b-5p	1.83E-05	hsa-miR-4636	1.94E-06
hsa-miR-27b-3p	1.95E-05
hsa-miR-126-3p	0.004079
hsa-miR-130a-3p	0.000442
hsa-miR-142-3p	0.002661
hsa-miR-4291	1.12E-09
hsa-miR-24-3p	6.14E-06
hsa-let-7a-5p	0.0059
hsa-miR-145-5p	0.00066
hsa-miR-17-5p	3.84E-05
hsa-miR-143-5p	6.07E-05
hsa-let-7f-5p	0.001364
hsa-miR-4328	0.001281
hsa-miR-4324	3.91E-05
hsa-miR-145-3p	0.000389
hsa-miR-143-3p	0.006402
hsa-miR-95	0.000106

¹Upregulated DEMs are listed according to the rank of fold changes. DEMs: Differentially expressed miRNAs.

Figure 1 The volcano plot of the differentially expressed miRNAs. A total of 110 DEMs were identified between 20 pairs of GC patients and adjacent normal tissues (cut-off criteria are P < 0.05 and fold change > 2.0). The green and red spots represent downregulated and upregulated miRNAs, respectively. DEMs: Differentially expressed miRNAs; GC: Gastric cancer.

Figure 2 Hierarchical clustering analysis of the differentially expressed miRNAs between 20 pairs of gastric cancer and adjacent normal sample. Each row represents the expression of an miRNA and each column represents a sample: orange color for gastric cancer, blue for normal.

Identification of DEMs related with overall survival in GC

To identify the DEMs which could be used to predict the overall survival of GC patients, we collected 361 samples from TCGA to assess the relationship between DEMs and the overall survival of GC patients. The patients’ clinical characteristics including age at diagnosis, gender, race, TNM stage, and histologic grade are shown in Table 2. By using a log-rank test and Kaplan-Meier curve, we found that three DEMs (miR-145-3p, miR-125b-5p, and miR-99a-5p) were negatively associated with overall survival (Figure 3). The association analysis between the three DEMs and clinical characteristics indicated that miR-145-3p, miR-125b-5p, and miR-99a-5p were all significantly associated with histologic grade (P < 0.05). The detailed results are shown in Table 3.

Table 2 Clinical features of gastric cancer patients.

Variables	Case, n (%)
Age at diagnosis (yr)
< 60	113 (31.3)
≥ 60	248 (68.7)
Gender
Male	241 (66.8)
Female	120 (33.2)
T stage
T1 + T2	88 (24.4)
T3 + T4	273 (75.6)
Histologic grade
G1 + G2	134 (37.1)
G3 + G4	227 (62.9)
Race
White	234 (64.8)
Asian	83 (23.0)
Black or African American	12 (3.3)
NA	32 (8.9)
Pathologic stage
I	44 (12.2)
II	117 (32.4)
III	162 (44.9)
IV	30 (8.3)
NA	8 (2.2)
Node status
N0	111 (30.7)
N1-3	249 (69.0)
NA	1 (0.3)
Metastasis
M0	325 (90.0)
M1	22 (6.1)
Mx	14 (3.9)

NA: Not available.

Table 3 Association between the three differentially expressed miRNAs and clinical features.

Variables	miR-145-3p expression		P value	miR-125b-5p expression		P value	miR-99a-5p expression		P value
	Low	High		Low	High		Low	High
Age at diagnosis (yr)
< 60	51	62	0.225	42	71	0.001a	49	64	0.096
≥ 60	129	119		138	110		131	117
T stage
T1 + T2	40	48	0.342	51	37	0.081	46	42	0.603
T3 + T4	140	133		129	144		134	139
N stage
N0	58	53	0.567	57	54	0.731	57	54	0.731
N1-3	119	124		120	123		120	123
M stage
M0	163	162	0.670	165	160	0.191	164	161	0.386
M1	10	12		8	14		9	13
Histologic grade
G1 + G2	77	57	0.026a	87	47	< 0.001a	80	54	0.004a
G3 + G4	103	124		93	134		100	127
Pathologic stage
I + II	81	80	0.876	86	75	0.221	80	81	0.954
III + IV	95	97		90	102		96	96

^aP < 0.05, statistically significant.

Figure 3 Three differentially expressed miRNAs are negatively associated with overall survival in gastric cancer. A: Patients with a high level of miR-145-3p (n = 180) had a poorer prognosis than those with a low level of miR-145-3p (n = 181) (P = 0.0049). B: Patients with a high level of miR-125b-5p (n = 180) had a poorer prognosis than those with a low level of miR-125b-5p (n = 181) (P = 0.0063). C: Patients with a high level of miR-99a-5p (n = 180) had a poorer prognosis than those with a low level of miR-99a-5p (n = 181) (P = 0.047).

Prognostic role of a three-DEM signature in GC patients

We ranked the three DEMs by the median of expression and then scored each GC patient in accordance with high or low-level expression. A risk grade was defined by the total scores. As a result, all the 361 GC patients were sorted into a high or low risk group. A survival analysis with the Kaplan-Meier method and log-rank test was conducted. The results indicated that the overall survival between the high risk and low risk groups was significantly different (P = 0.045). Interestingly, compared to patients in the high risk group, the low risk patients tended to have a better prognosis (Figure 4). Furthermore, we performed univariate and multivariate Cox regression analyses to verify the prognostic role of the three-DEM signature according to clinical features. The univariate analysis showed that pathologic stage (HR = 1.825, P < 0.001), T stage (HR = 1.864, P = 0.006), N stage (HR = 2.005, P = 0.001), and the three-DEM signature (HR = 1.422, P = 0.039) were significantly associated with the prognostic outcome of GC patients. The multivariate analysis revealed that T stage (HR = 1.623, P = 0.044) and the three-DEM signature (HR = 1.451, P = 0.032) were all independent factors in predicting the prognosis of GC patients (Table 4).

Table 4 Univariate and multivariate Cox regression analyses of the association between the three differentially expressed miRNAs and clinical features.

Variables	Univariate analysis		Multivariate analysis
	HR (95%CI)	P value	HR (95%CI)	P value
Age at diagnosis (≥ 60 vs < 60)	1.373 (0.948-1.988)	0.094	1.642 (1.122-2.401)	0.011a
Pathologic stage (III + IV vs I + II)	1.825 (1.332-2.499)	< 0.001a	1.252 (0.812-1.929)	0.309
T stage (T3 + T4 vs T1 + T2)	1.864 (1.197-2.902)	0.006a	1.623 (1.012-2.603)	0.044a
N stage (N1-2 vs N0)	2.005 (1.328-3.026)	0.001a	1.602 (0.935-2.744)	0.086
M stage (M1 vs M0)	1.368 (0.964-1.941)	0.080	1.313 (0.919-1.875)	0.134
Three-DEM signature (high vs low risk)	1.442 (1.018-1.988)	0.039a	1.451 (1.033-2.040)	0.032a

^aP < 0.05, statistically significant. DEMs: Differentially expressed miRNAs.

Figure 4 The Kaplan-Meier curve for the three-miRNA signature in gastric cancer. The three differentially expressed miRNAs were ranked by the median of expression and then scored for each gastric cancer patient in accordance with high or low-level expression. The low risk group (n = 175) and high risk group (n = 186) were defined by the total scores. Compared to the low risk group, patients in the high risk group had a poorer prognosis (P = 0.045).

Target prediction of three DEMs and gene function analysis

The online target prediction tools TargetScan, miRDB, miRWalk, and DIANA were used to predict the targets genes of miR-145-3p, miR-125b-5p, and miR-99a-5p. We then obtained the consensus genes of each DEM from the four online predictions (Figure 5). As a result, we identified a 108 consensus target genes. Furthermore, we conducted gene enrichment analysis to identify the biological function of common target genes (Figure 6). The Biological Process analysis indicated that the genes were mostly enriched in cell recognition, regulation of nucleic acid, and gene silencing. Cellular Component analysis indicated that genes were enriched in the nucleus, RNA-induced silencing complex, chromatin silencing complex, and phosphoinositide 3-kinase complex. Molecular Function analysis showed that genes were enriched in transmembrane receptor, transcription regulator, and transcription factor activity. Biological Pathways were mainly enriched in the VEGFR, PI3K/Akt, and mTOR signaling pathways.

Figure 5 The consensus target genes of each differentially expressed miRNA. The target genes of the three miRNAs were predicted with four online tools (TargetScan, miRDB, miRWalk, and DIANA).

Figure 6 Genetic functional enrichment analysis. A total of 108 consensus genes were used for functional enrichment analysis with the tool of FunRich. A: Biological process analysis showed that the genes are involved in cell recognition and gene silencing; B: Cellular component analysis indicated that the genes are portion of nucleus, chromatin silencing complex, and TORC1/2 complex; C: Molecular function results showed that the genes are involved in nucleic acid binding, transcription factor activity, and transmembrane receptor activity; D: Biological pathway analysis indicated that the genes participate in focal adhesion, PI3K, and mTOR cancer-related signaling pathways. TORC1/2: Target of rapamycin 1/2; PI3K: Phosphatidylinositide 3-kinase; mTOR: Mammalian target of rapamycin.

DISCUSSION

Due to the reduction in chronic Helicobacter pylori infection and improvement of sanitation, the incidence and mortality rates of GC have declined in recent years[12]. However, there are still almost 460,000 new GC cases and 350000 GC deaths each year in China[13]. The prognosis of GC patients is poor and the five-year survival rate is 5%-20% despite advances in GC therapy[14]. Thus, to improve the clinical treatment and management of GC patients, it is urgent to identify reliable prognostic biomarkers. In this study, we identified a total of 110 DEMs by analyzing the GSE93415 data and discovered that three miRNAs (miR-145-3p, miR-125b-5p, and miR-99a-5p) were negatively associated with overall survival. Additionally, we constructed a three-miRNA signature to predict the prognosis of GC patients.

For decades, a large number of studies have reported that miRNAs can play oncogenic or tumor-suppressing roles in regulating cell biological behavior of cells[15-18]. At present, several miRNAs are known to be useful in the early diagnosis of cancers, including miR-21[19], miR-486[20], miR-24[21], and miR-125a-5p[22]. In addition, miR-191[23], miR-1908[24], miR-200c[25], and miR-217[26] were found to be potential prognostic indictors in cancer. However, these studies only used a single indictor or a limited number of patients for survival analysis. In this study, we identified DEMs by analyzing the array data from the GEO database and found that three highly expressed miRNAs (miR-145-3p, miR-125b-5p, and miR-99a-5p) may be potential prognostic indictors in GC. Kaplan-Meier and Log-rank test survival analysis indicated that the three-miRNA signature can be used to predict the prognosis of GC patients.

We then searched present publications online to compare and test our findings. Chang et al[27] showed that miR-125b-5p was overexpressed in GC patients and promoted invasion and metastasis of GC by targeting STARD13 and NEU1. It was also indicated that miR-125b-5p could be a potential biomarker for predicting prognoses and clinical outcomes in patients with HER2-positive GC that receive trastuzumab treatment[28]. Wu et al[29] found that miR-125b-5p promotes cell migration and invasion by targeting PPP1CA-Rb signal pathways and acts as an independent prognostic factor in GC. Furthermore, Zhang et al[30] demonstrated that miR-99a-5p might function as a novel molecule to regulate cisplatin resistance by directly targeting the calpain small subunit 1 (CAPNS1)-associated pathway in GC. Interestingly, there are no studies describing relations between miR-145-3p and GC. However, in non-small cell lung cancer, miR-145-3p was found to inhibit cancer cell migration and invasion by targeting PDK1 via the mTOR signaling pathway[31]. Moreover, miR-145-3p was also identified to be down-regulated in metastatic castration-resistant prostate cancer and target four molecules which can significantly predict survival in prostate cancer[32]. These results may suggest that miR-145-3p has a complicated effect in different cancers and, in future research, we will investigate the role of miR-145-3p in GC.

Dysregulated genes may participate in tumorigenesis and progression by aberrant signaling pathways. In this study, we predicted the target genes of the three miRNAs and performed gene functional enrichment analysis. The results showed that these target genes were associated with the process of gene silencing and cell recognition, as well as focal adhesion, EGFR, PI3K/Akt, and mTOR signaling pathways. Xu et al[33] suggested that the EGFR-Akt signaling pathway regulates drug resistance in GC patients. The PI3K/Akt pathway was demonstrated to be associated with poor prognosis, tumor progression, and resistance to systematic therapy in many cancers including GC[34,35]. In addition, the PI3K/Akt/mTOR pathway is a key signaling pathway that is reported to be involved in GC[36]. Thus, the results of our functional enrichment analysis are in accordance with present studies.

Above all, we identified a three-miRNA signature for predicting the prognosis of patients with GC and analyzed potential signaling pathways in the development and progression of GC. However, to determine the genesis and development mechanism of GC, more large-scale and systematic investigations are required.

ARTICLE HIGHLIGHTS

Research background

Increasing studies have reported that microRNAs (miRNAs) play an important role in the development and progression of cancers, including gastric cancer (GC). Furthermore, miRNAs can also act as accurate biomarkers in diagnosis and prognosis prediction. In this study, we found that a three-miRNA signature could be used for predicting the prognosis of GC patients and multiple miRNAs together acting as biomarkers may have a stronger reliability in survival prediction.

Research motivation

The worldwide incidence and mortality rates of GC are fairly high. Most of GC patients have been in the advanced stage when diagnosed and endure a poor prognosis. Identifying accurate biomarkers in predicting prognosis of patients is an urgent issue to be solved, so that patients could have an individualized treatment and an improvement in prognosis.

Research objectives

We aimed to identify multiple miRNAs for predicting the prognosis of patients with gastric cancer. In the present study, we found that a three-miRNA (miR-145-3p, miR-125b-5p, and miR-99a-5p) signature could be used for predicting the prognosis of patients with gastric cancer. This objective could be applied to clinical practice and have a guidance role in improving the prognosis of patients with gastric cancer.

Research methods

We obtained the differentially expressed miRNAs by analyzing a microarray dataset from the Gene Expression Omnibus database with the limma package in R. The Kaplan-Meier method and log-rank test were used for describing the survival curve. The target genes of the three miRNAs (miR-145-3p, miR-125b-5p, and miR-99a-5p) were predicted with the online tools of TargetScan, miRDB, miRWalk, and DIANA. Venn plot was performed to obtain the common target genes from these four online tools. Enrichment analysis was conducted on the consensus genes using the FunRich tool.

Research results

In the present study, we found that a three-miRNA (miR-145-3p, miR-125b-5p, and miR-99a-5p) signature could be used for predicting the prognosis of patients with gastric cancer. Multiple miRNAs together acting as prognosis-related biomarkers may have a stronger reliability and this finding could be useful in clinical treatment according to gastric cancer patients with different prognoses. However, the role of miR-145-3p in the tumorigenesis and progression of gastric cancer remains unclear. Thus, this problem remains to be solved in our future study.

Research conclusions

Our study identified that the three miRNAs (miR-145-3p, miR-125b-5p, and miR-99a-5p) were up-regulated in gastric cancer patients by analyzing a microarray dataset. Besides, the novel three-miRNA signature could be used for predicting the prognosis of patients with gastric cancer. Multiple miRNAs together acting as prognosis-related biomarkers may have a stronger reliability, so that our finding could be useful in clinical treatment according to gastric cancer patients with different prognoses.

Research perspectives

This study provides us with a new insight that multiple miRNAs can be together used for predicting the prognosis of patients with gastric cancer. In order to further confirm the prognostic value of the three-miRNA signature, our future research may focus on exploring the relation between miR-145-3p and gastric cancer.