This article is mentioned in:

Introduction

Esophageal squamous cell carcinoma (ESCC), which begins in flat cells lining the esophagus, is a common esophageal malignant tumor globally, and particularly in China (1). Alcohol intake and tobacco smoking are major risk factors for ESCC (2,3). It is possible to treat small localized ESCC (<2 cm) with surgery, while larger tumors may be slowed with chemotherapy and radiation therapy; however, as diagnosis is frequently late the outcomes tend to be relatively poor (4,5). Therefore, multiple studies have been performed to study the pathogenesis of ESCC and to identify more effective molecular methods of diagnosis.

Lee et al (6) suggested that hereditary genetic polymorphisms of autophagy related 5 and collagen type IV α3 chain may be prognostic predictors for patients with ESCC. Elevated expression of AXL receptor tyrosine kinase in esophageal tumor tissues has been suggested to be associated with increased risk of death and recurrence of ESCC (7). Previous reports have identified genetic alternations of oncogenes and tumor suppressors in ESCC (8,9), but the genetic and molecular basis for esophageal carcinogenesis remains largely unknown. Altered gene expression, which is the result of gene mutation or dramatic changes in gene regulation, is associated with a variety of diseases. MicroRNAs (miRNAs) act as oncogenes or tumor suppressors, mediating multiple gene expression in human cancers (10,11). Previous studies demonstrated that in ESCC, tumor-suppressive miRNAs of miR-145, miR-133a, miR-133b (12) and miR-143 (13) inhibited tumor cell proliferation and invasion. MiR-155 (14), miR-21 (15) and miR-31 (16) act as oncogenes in ESCC. miRNA expression is correlated with carcinogenesis, and associated with prognosis and therapeutic outcomes. Guo et al (17) suggested that high expression of miR-103/107 was correlated with poor survival of patients with ESCC. miR-142-3p is a potential prognostic biomarker for ESCC (18). However, whether these miRNAs are efficient enough or whether it is possible to use them to evaluate prognosis in patients with ESCC remains unclear. Therefore, further studies are required to improve the survival rates of patients with ESCC.

Using the non-coding RNA (miRNA) profile GSE43732, Chen et al (19) identified four miRNA signatures which had good prognostic value, and investigated the replicability and accuracy of these miRNAs as predictors. The present study used different methods to screen prognosis-associated miRNAs, and the results may provide novel predictable biomarkers and contribute to improved survival in patients with ESCC.

Materials and methods

Microarray data

The non-coding RNA profile GSE43732, which was contributed by Chen et al (19), was downloaded from the Gene Expression Omnibus database (20). Paired frozen tissues (cancer tissue and adjacent normal tissue) taken from 119 patients were assessed by microarray for miRNA expression analysis. This data was obtained based on the platform of GPL16543 Agilent-038166 cbc_human_miR18 (Cancer Institute and Hospital, Chinese Academy of Medical Sciences, Beijing, China).

Data preprocessing

Series Matrix Files containing normalized chip signals were transformed to gene symbols. The aggregate function in R package (bioconductor.org/packages/2.1/bioc/html/impute.html) was used to calculate the mean value when multiple probes were mapped to the same gene. The k-nearest neighbor method (21) from the impute package (22) in R was used with the default k-value of 10. Quantile normalization of miRNA expression profile was performed with the preprocessCore package (23) in R. The distribution of miRNA expression prior to and following normalization was visualized as a box plot.

Screening for differentially expressed miRNAs

Differentially expressed miRNAs were screened using the limma package (24) in R. The difference of mean expression level between two groups of samples was examined using the paired Student's t-test in the limma package. The P-value was adjusted using the Benjamini-Hochberg method (25) to control the false discovery rate. Genes with log2 fold change >2 and adjusted P<0.05 were considered to be differentially expressed miRNAs.

Identification of prognosis-related differentially expressed miRNAs

The high and low critical values of miRNA expression were derived based on the expression values of miRNAs and the survival time of patients using receiver operating characteristic (ROC) curve analysis. The ROC curves were generated and analyzed using Proc (26), KMsurv (27) and survival package (28) in R. Patients with lower and higher values compared with the critical values were divided into two groups. Subsequently, the significant variation of survival time between these two groups was examined using the Kaplan-Meier (KM) method and log rank test, with P<0.05 considered to indicate a statistically significant difference.

Only miRNAs with an area under the curve (AUC) of ≥0.7 were selected. They were subsequently screened again for miRNAs which were highly expressed in tumor tissues and associated with worse prognosis than others. Similarly, those miRNAs expressed at low levels in tumor tissues and associated with a worse prognosis were also identified.

Construction of the microRNA regulatory network for prognostic risk

Databases providing information on miRNA-target genes were used to identify the target genes of differentially expressed miRNAs. Verification databases included miRNecords (29) and MiRWalk (30), and predictive databases included miRanda (31), MirTarget2 (32), PicTar (33), PITA (34) and TargetScan (35). The identified target genes should be predicted by at least three of these databases. Target genes obtained from the two types of database were used to construct the miRNA-regulatory network, which was subsequently visualized using Cytoscape software (36).

Functional enrichment analysis of target genes

Gene Ontology (GO, http://www.geneontology.org/) analysis consisting of molecular function, cell component and biological process (BP), the Kyoto Encyclopedia of Genes and Genomes (KEGG) pathways (Reactome database, http://www.reactome.org/ and KEGG database, http://www.genome.jp/kegg/), and the Disease Ontology (BMC Genomics 10 Suppl 1:S6) enrichment analyses (37) for miRNAs-target genes were performed using the TargetMine (targetmine.nibio.go.jp) (38) online tool. The hypergeometric test within TargetMine for enrichment results was performed and the P-value was adjusted with the Holm-Bonferroni correction. An adjusted P<0.05 was considered to indicate a statistically significant difference.

Results

Data normalization and differentially expressed miRNAs

Following normalization, a total of 726 miRNA expression values were obtained (Fig. 1). A total of 107 differentially expressed miRNAs, including 54 upregulated and 53 downregulated miRNAs, were obtained. The volcano plot of differentially expressed miRNAs is depicted in Fig. 2.

Figure 1. Box plots for each array for (A) adjacent normal tissue and (B) cancer tissue prior to normalization, and (C) adjacent normal tissue and (D) cancer tissue following normalization.

Figure 2. Volcano plot of differentially expressed miRNAs. The horizontal axis represents log2 fold change and the vertical axis represents-log10 (adjusted P-value).

Prognosis related miRNAs

Critical values of the 107 differentially expressed miRNAs were determined. Kaplan-Meier (KM) curve analysis revealed that 44 miRNAs had significant effects on prognosis (P<0.05; data not shown). Among them, 27 miRNAs were upregulated and 17 miRNAs were downregulated. Furthermore, 12 miRNAs, including 9 upregulated and 3 downregulated miRNAs, were screened with AUC ≥0.7. In the end, hsa-miR-143-3p and hsa-miR-145-5p were revealed to be highly expressed in tumor tissues and were associated with prognosis (Fig. 3). Meanwhile, hsa-miR-182-5p and hsa-miR-455-5p, which were also associated with prognosis, were expressed at low levels in tumor tissues (Fig. 3).

Figure 3. ROC curves (left) and Kaplan-Meier curves (right) of four miRs associated with good prognosis. The critical value is on the top left corner of the ROC curves. ROC, receiver operating curve; miR, microRNA.

Regulatory network of prognosis related miRNAs and target genes

A total of 8 target genes of hsa-miR-455-5p were obtained from the verification database, but not the other three miRNAs (hsa-miR-143-3p, hsa-miR-145-5p and miR-182-5p). The target genes were activating transcription factor 6, interleukin 10 (IL10), toll-like receptor 4 (TLR4), activator of heat shock 90 kDa protein ATPase homolog 1, scavenger receptor class B member I, prostaglandin-endoperoxide synthase 2, calreticulin and proopiomelanocortin. From the predicted databases, a total of 245 target genes including ring finger protein (RNF) 44/111, lysine demethylase (KDM) 3B/6A and chromosome 15 open reading frame 26/41, and the miRNA regulatory network is depicted in Fig. 4.

Figure 4. Regulatory network of hsa-miR-455-5p. (A) Regulatory network of eight verified target genes of hsa-miR-455-5p. (B) Regulatory network of 8 verified and 245 predicted genes targeted by hsa-miR-455-5p. hsa-miR-455-5p, human microRNA 455-5p.

Functional enrichment analysis of target genes of miRNAs

No pathway was significantly enriched by the 253 target genes of hsa-miR-455-5p; however, two GO-BP terms were statistically significant (Table I). These were associated with the macromolecule metabolic process. The Disease Ontology enrichment results revealed that 15 target genes of hsa-miR-455-5p were enriched in term of melanoma (Table I).

Table I. GO and DO enrichment analyses of target genes.

Table I.

GO and DO enrichment analyses of target genes.

Category	Term	P-value	Count
GO-BP	Macromolecule metabolic process [GO:0043170]	0.00763	117
GO-BP	Cellular macromolecule metabolic process [GO:0044260]	0.00860	109
DO term	Melanoma [DOID:1909]	0.00215	15

[i] GO, Gene Ontology; DO, Disease Ontology; BP, biological process; count, the number of target genes enriched in this term.

Discussion

The poor prognosis of ESCC frequently results in severe adverse effects to patients. Novel, effective molecular markers for ESCC are, therefore, urgently required. In the present study, several miRNAs associated with prognosis of ESCC were discovered. miRNAs including miR-143-3p and miR-145-5p were upregulated, while miR-182-5p and miR-455-5p were downregulated in ESCC tissues. Previous research has demonstrated the involvement of miR-143-3p (39) and miR-145-5p in ESCC (40,41), which were recorded in the MalaCards database; however, not the downregulated miRNAs miR-182-5p and miR-455-5p. miR-182-5p and miR-455-5p may be novel miRNAs associated with ESCC.

The association between miR-182 and tumors has been widely studied. Segura et al (42) demonstrated that abnormal miR-182 expression represses forkhead box O3 and microphthalmia-associated transcription factors to promote melanoma metastasis, and contribute to renal cell carcinoma proliferation (43). miR-182 may be a prognostic marker for glioma progression and its expression was associated with patient survival (44). In the present study, miR-182-5p was revealed to be under expressed in tumor tissues and was associated with prognosis. Therefore, miR-182-5p may be a novel miRNA associated with the prognosis and survival rate of ESCC.

miR-455-5p was another downregulated miRNA identified in the present study. The expression pattern and its association with prognosis was similar to that of miR-182-5p. Swingler et al (45) hypothesized that miR-455 may regulate apoptosis in certain mouse embryo tissue. Besides, decreased expression of miR-455-5p has been demonstrated to be correlated with vascular invasion and poor overall survival, and to be involved in endometrial serous adenocarcinomas progression (46). Based on these previous findings and the results of the present study, miR-455-5p may be a novel miRNA participating in the development and prognosis of ESCC. In addition, multiple target genes were identified. IL-10 (47) and TLR4 (48), which were overexpressed in the results of the present study, were demonstrated to be involved in ESCC. Among the predicted target genes, no genes recorded in the MalaCards database were associated with ESCC, but the homologous genes including RNF6 and KDM4C were illustrated in this database. Furthermore, the Disease Ontology enrichment results demonstrated that the target genes were tumor-associated. In addition, IL-10 and TLR4, as well as multiple other target genes of miR-455-5p, were revealed to be enriched in macromolecule biosynthesis and metabolism, which is associated with tumor initiation and progression (49,50). Cellular metabolic perturbation is a crucial hallmark of cancer, and is associated with angiogenesis, evasion of apoptosis, metastasis and avoidance of immune detection (51). Macromolecule biosynthesis and metabolism are critical to supplying enough carbohydrates, proteins, lipids and nucleic acids to support the rapid cell division that occurs within a tumor (52). Therefore, miR-455-5p may be another novel miRNA that participates in the progression of ESCC and is associated with its prognosis.

In conclusion, in addition to the two known miRNAs (miR-143-3p and miR-145-5p), two novel tumor suppressor miRNAs, miR-182-5p and miR-455-5p, were demonstrated to be associated with ESCC progression. As multiple target genes of miR-455-5p were identified, it may be involved in complicated regulatory mechanisms and be involved in ESCC. Although there were no identified target genes for miR-182-5p, it may also be important in ESCC and warrants of further research. Whether the two miRNAs are effective targets for ESCC requires confirmation and further study in future experiments.

References