Identification of a novel prognosis-associated ceRNA network in lung adenocarcinoma via bioinformatics analysis

Abstract

Background: Lung adenocarcinoma (LUAD) is the most common subtype of nonsmall-cell lung cancer (NSCLC) and has a high incidence rate and mortality. The survival of LUAD patients has increased with the development of targeted therapeutics, but the prognosis of these patients is still poor. Long noncoding RNAs (lncRNAs) play an important role in the occurrence and development of LUAD. The purpose of this study was to identify novel abnormally regulated lncRNA-microRNA (miRNA)-messenger RNA (mRNA) competing endogenous RNA (ceRNA) networks that may suggest new therapeutic targets for LUAD or relate to LUAD prognosis.

Methods: We used the SBC human ceRNA array V1.0 to screen for differentially expressed (DE) lncRNAs and mRNAs in four paired LUAD samples. Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway analyses were performed to annotate the DE lncRNAs and mRNAs. R bioinformatics packages, The Cancer Genome Atlas (TCGA) LUAD database, and Kaplan-Meier (KM) survival analysis tools were used to validate the microarray data and construct the lncRNA-miRNA-mRNA ceRNA regulatory network. Then, quantitative real-time PCR (qRT-PCR) was used to validate the DE lncRNAs in 7 LUAD cell lines.

Results: A total of 2819 DE lncRNAs and 2396 DE mRNAs (P < 0.05 and fold change ≥ 2 or ≤ 0.5) were identified in four paired LUAD tissue samples. In total, 255 of the DE lncRNAs were also identified in TCGA. The GO and KEGG analysis results suggested that the DE genes were most enriched in angiogenesis and cell proliferation, and were closely related to human cancers. Moreover, the differential expression of ENST00000609697, ENST00000602992, and NR_024321 was consistent with the microarray data, as determined by qRT-PCR validation in 7 LUAD cell lines; however, only ENST00000609697 was associated with the overall survival of LUAD patients (log-rank P = 0.029). Finally, through analysis of ENST00000609697 target genes, we identified the ENST00000609697-hsa-miR-6791-5p-RASL12 ceRNA network, which may play a tumor-suppressive role in LUAD.

Conclusion: ENST00000609697 was abnormally expressed in LUAD. Furthermore, downregulation of ENST00000609697 and its target gene RASL12 was associated with poor prognosis in LUAD. The ENST00000609697-hsa-miR-6791-5p-RASL12 axis may play a tumor-suppressive role. These results suggest new potential prognostic and therapeutic biomarkers for LUAD.

Keywords: Biomarker; CeRNA; LncRNAs; Lung adenocarcinoma; Prognosis.

Fig. 1

The flowchart of this study

Fig. 2

Screening of DE lncRNAs and mRNAs between LUAD tissues and normal tissues. A Heatmap of DE lncRNAs; B Heatmap of DE mRNAs; the x axis shows samples and the y axis shows DE genes. Red and blue represent upregulated and downregulated genes, respectively. C Scatter and volcano plots showing the expression profiles of DE lncRNAs based on the lncRNAs expression values detected by microarray. D Scatter and volcano plots showing expression profiles of DE mRNAs based on the expression values of mRNAs detected by microarray; |FC|≥ 2.0, P < 0.05. E Heatmap of DE lncRNAs in 512 LUAD tissue samples and 59 adjacent tissue samples. Green and red represent downregulated and upregulated lncRNAs, respectively. F Volcano plot of the 1271 upregulated lncRNAs and 645 downregulated lncRNAs that were screened (|log₂FC|> 1, P < 0.05). G Venn diagram of DE lncRNAs identified by the TCGA and microarray analyses

Fig. 3

GO and KEGG enrichment analysis of DE lncRNAs and mRNAs. Barplot of the top 30 enriched GO classification terms for DE lncRNAs (A) and DE mRNAs (B). Bubble plot of the top 30 GO level 2 terms enriched by the DE lncRNAs (C) and DE mRNAs (D). Bubble plot of the top 30 KEGG pathways enriched by the DE lncRNAs (E) and DE mRNAs (F). G, H Barplot of the top 30 KEGG classifications for the DE lncRNAs (G) and DE mRNAs (H). GeneRatio ≥ 2, P < 0.05

Fig. 4

GO and KEGG enrichment analysis of the intersection between DE lncRNA target genes and DE mRNAs. A Venn diagram of DE mRNAs and DE lncRNA target genes. B Cell component (CC). C Biological process (BP). D Molecular function (MF). E KEGG pathway. The y‐axis shows significantly enriched categories for the targets, and the x axis shows the enrichment scores of these terms. The bar plot height indicates the number of genes in the functional area

Fig. 5

Validation of candidate lncRNAs. A The expression of the downregulation candidate lncRNAs in 7 LUAD cell lines and BEAS-2B cells was determined by qRT-PCR. B The expression of upregulation candidate lncRNAs in 7 LUAD cell lines and BEAS-2B cells was determined by qRT-PCR. The data are presented as the mean ± standard error of three independent experiments. *P < 0.05; **P < 0.01; ***P < 0.001. C The relative expression of the candidate DE lncRNA ENST00000609697 in the TCGA dataset. D Kaplan–Meier (KM) survival analysis based on the candidate DE lncRNA ENST00000609697. x axis: overall survival (years); y axis: survival rate. Green and red represent the low and high lncRNA expression groups, respectively

Fig. 6

Target gene functional annotation of the ENST00000609697 ceRNA regulatory network. A CeRNA network of ENST00000609697. B Bar plot of the top 30 GO classification terms enriched by the target genes, P < 0.05. C BP cnetplot of the target genes, P < 0.001. D BP emapplot of the target genes, P < 0.001. E Relative expression of RASL12 in the TCGA dataset, P < 0.0001. F Survival analysis based on RASL12 in the TCGA dataset, P = 0.034