Introduction
Each year ~990,000 people are diagnosed with gastric cancer (GC) and ~738,000 die from this disease worldwide, making GC the fourth most common cancer and the second most common cause of cancer death (1). The highest incidence rates are observed in South America, East Europe, and East Asia (particularly in Korea, Mongolia, Japan and China) (2). Detection of GC in the early stage is important to reduce mortality. Many Japanese series have consistently reported early diagnosis and treatment of GC, with a 5-year survival rate >90%. The key factor for preventing GC is early diagnosis, as confirmed by the good results obtained in surgical series in which early cases are frequent (3). Recent developments in the field of gastrointestinal (GI) endoscopy have been remarkable (4). However, endoscopic biopsy and pathological morphological observation cannot find all precancerous lesions and early GC. To improve this situation, there is an urgent need for novel and reliable biomarkers for early diagnosis of GC.
In recent years, several studies have confirmed that long non-coding RNAs (lncRNAs) have emerged as essential regulators in almost all aspects of biology (5). They are a group of RNA transcripts that are >200 nucleotides in length and lack significant open reading frames (6). Accumulating evidence suggests that lncRNAs play an important role in tumorigenesis (7). lncRNAs may function as oncogenes or tumor suppressors by altering the chromatin structure or by regulating the transcription of protein-coding genes (8). lncRNA is the functional end-product, and the level of lncRNA expression correlates directly with the level of the active molecule. Thus, the use of lncRNAs in diagnostics has intrinsic advantages over the use of protein-coding RNAs. In addition, lncRNAs show greater tissue specificity compared to protein-coding mRNAs and miRNAs, making them attractive in the search for novel diagnostic and prognostic cancer biomarkers (9). It is known that gastric adenocarcinoma accounts for 95% of gastric malignancies. However, there is little information to clarify the relationship between lncRNAs, mRNAs and validation in advanced and early stage GC (10). Therefore, further study of endoscopic biopsy of GC tissues may help to establish inner relationships between lncRNA, mRNA and GC.
In the present study, differential expression profiles of lncRNAs and mRNAs were detected from advanced GC tissues and adjacent non-tumor tissues by microarray analysis. We investigated the targeted and regulated genes of lncRNAs through construction of a co-expression network and Gene Ontology (GO) analysis. Then, we used real-time PCR validation of microarray differential expression lncRNAs and mRNAs in tumor tissues and adjacent non-tumor tissues from 30 newly diagnosed advanced GC patients and 20 newly diagnosed early stage GC patients.
Materials and methods
Specimen collection
The procedures used in the present study were approved by the Wuwei Tumor Hospital of Gansu (Wuwei, China) and conformed to the Helsinki Declaration, and current legislation. A total of 10 human advanced GC specimens and their paired adjacent non-cancerous tissue specimens were obtained for lncRNA and mRNA microarray analysis from the Wuwei Tumor Hospital of Gansu (Table I). In addition, 30 advanced GC patients and 20 newly diagnosed early stage GC patients were recruited from the Wuwei Tumor Hospital of Gansu, aged 45–70 years, for quantitative RT-PCR analysis between 2014–2015. All of these patients were assigned a diagnosis of GC based on histopathology and clinical history. Clinical information that was recorded for each specimen included age, tumor grade, cancer stage, tissue dimensions and date of resection. The pathologist assessed the tumor content by microscopic examination of cases where the percentage of tumor tissue was estimated to be ≥80%. None of the patients had received preoperative chemoradiation. Adjacent non-cancerous tissues were located ≥5 cm from the tumor edge. Tissue samples were immersed in RNAlater (Ambion, Inc., Austin, TX, USA) and stored at −80°C until use.
Isolation of RNA
Total RNA was isolated from 10 advanced gastric adenocarcinoma specimens and their paired adjacent non-cancerous tissues using TRIzol reagent (Invitrogen, Carlsbad, CA, USA). The quality of RNA was measured by Agilent Bioanalyzer, which produces an RNA integrity number (RIN) between 1 and 10, with 10 being the highest quality samples, showing the least degradation. The RIN estimates sample integrity using gel electrophoresis and analysis of the ratios of 28S to 18S ribosomal bands.
lncRNA and mRNA microarray analysis
Three pairs of RNA samples were obtained from 10 patients by pooling RNA from two, four and four patients within one group as 'one sample'. Thus, we generated three pairs of pooled RNA samples from GC specimens and their paired adjacent non-cancerous tissues from 10 patients. These three pairs of RNA samples were subjected to microarray analysis. RiboArray Custom Array 1*90K combined microarray was used, which could detect 32,987 lncRNAs from authoritative data sources, including NCBI RefSeq, H-invDB, UCSC, LncRNAdb, John Rinn and GENCODE LncRNA. In addition, the microarray also included Human Gene Expression, which could detect 65,535 mRNAs. Signals were normalized using the median center tool for genes. ANOVA was used to compare the differentially expressed lncRNAs and mRNAs.
Differential expression of mRNA function
Differentially expressed mRNAs were entered into the Database for Annotation, Visualization and Integrated Discovery (DAVID), which utilized GO to identify the molecular function represented in the gene profile and the Kyoto Encyclopedia of Genes and Genomes (KEGG) to analyze the potential functions of these genes in the pathways (11).
Bioinformatics analysis of co-expression of lncRNAs and mRNAs
According to the relationship among miRNA, lncRNA and mRNA, the post-transcriptional regulation of mRNA transcripts bound by single-stranded miRNAs is basically established. In the present study we built lncRNA-mRNA co-expression network which is based on the theory that lncRNA can regulate miRNA abundance by sequestering and binding them, acting as so-called miRNA sponges. Firstly, the differentially expressed lncRNAs and mRNAs were selected from GC specimens and their paired adjacent non-cancerous tissues. Standard selection criteria to identify differentially expressed lncRNAs and mRNAs were established at P<0.05 and fold change >2. Then, lncRNA and mRNA co-expression network was constructed based on the correlation between the differentially expressed lncRNAs and mRNAs. The correlation of expression profiles between biological replicates and group conditions was demonstrated by unsupervised hierarchical cluster analysis (Figs. 1 and 2).
Quantitative real-time PCR validation
To confirm the microarray results, we performed qRT-PCR analysis on larger samples, including 30 newly diagnosed advanced GC and 20 early stage GC patients. Early stage GC patients were all diagnosed as high level intraepithelial neoplasia by pathological morphological observation. GC specimens and their paired adjacent non-cancerous tissues were collected and total RNA was isolated, then reverse-transcribed using AMV reverse transcriptase (Promega, Madison, WI, USA). Real-time PCR was done with GoTaq® qPCR Master Mix (Promega) on StepOnePlus™ PCR System (Applied Biosystems, Waltham, MA, USA). The relative fold change was calculated using the 2−ΔΔCt method (12), where ΔCt = (CtRNAs − CtGAPDH) and ΔΔCt = ΔCttumor tissues − ΔCtadjacent non-tumor tissues. GAPDH was chosen as the endogenous standard. The threshold cycle indicated the fractional cycle number at which the amount of amplified target reached a fixed threshold and the Ct value was negatively correlated with copy numbers. In addition, conditional logistic regression analysis was used to evaluate the association between differentially expressed RNAs and the risk of GC (13).
Statistical analysis
All the results were expressed as mean ± SD. Statistical analysis was carried out with the Student's t-test for comparison of two groups in microarray analysis and ANOVA for multiple comparisons. In both cases, differences with P<0.05 were considered to indicate a statistically significant result. The statistical significance of microarray analysis results was analyzed by fold change and Student's t-test. False discovery rate was calculated to correct the P-value. We used fold change to screen differentially expressed LncRNAs and mRNAs, the threshold values were (fold change ≥2.0 or ≤0.50) and (P<0.05).
Results
Quality of microarray sample RNAs
Microarray sample RNAs were assessed by electrophoresis on a denaturing agarose gel. Total RNA run on a denaturing gel has sharp 28S and 18S rRNA bands, and the 28S rRNA band should be approximately twice as intense as the 18S rRNA band. This intensity ratio was observed for the RNA in the present study, indicating that the RNA was intact. The concentration of RNAs (OD260), protein contamination of RNAs (ratio OD260/OD280), and organic compound contamination of RNAs (ratio OD260/OD230) were measured with the NanoDrop ND-1000. All samples had OD260/OD280 ratios of total RNA higher than 1.8, indicating adequate RNA concentration.
Overview of lncRNA and mRNA expression profiles
The lncRNA expression levels of GC tissues and their paired adjacent non-cancerous tissues were compared. We found that 427 upregulated and 619 downregulated lncRNAs were significantly differentially expressed (fold change ≥2.0). Cluster analysis was used to arrange specimens into groups according to their expression levels (Fig. 2). The top 40 differentially expressed lncRNAs are listed (Fig. 3). RP11-706O15 (fold change, 68.16) was the most upregulated lncRNA, and MID1 (fold change, 140.85) was the most downregulated. There were more downregulated than upregulated lncRNAs according to the microarray data. In addition, in the GC tissues and their paired adjacent non-cancerous samples, mRNA expression profiling data showed that 647 upregulated and 2,349 downregulated mRNAs were significantly differentially expressed (fold change ≥2.0). Hierarchical cluster analysis was used to arrange specimens into groups according to their expression levels (Fig. 2). The top 40 differentially expressed mRNAs are listed (Fig. 4). ZNF737 (fold change, 81.73) was the most upregulated mRNA, and MGME1 (fold change, 227.27) was the most downregulated. We found that there were more downregulated than upregulated mRNAs in this microarray analysis.
Function analysis of differentially expressed genes
Differential expression of mRNAs in GC tissues and adjacent non-cancerous tissues was analyzed by related (GO) analysis performed by the DAVID functional annotation chart. We analyzed the enrichment of these differentially expressed mRNAs. Enrichment provides a measure of the significance of the function, and as the enrichment increases, the corresponding function is more specific, which helps us to identify GO with a more definitive functional description (14). The results showed that the highest enriched GOs targeted by upregulated transcripts were DNA-dependent transcription (GO:0006351), proteolysis (GO:0006508) and mitotic cell cycle (GO:0000278). The highest enriched GOs targeted by downregulated transcripts were DNA-dependent transcription (GO:0006351), transmembrane transport (GO:0055085) and regulation of DNA-dependent transcription (GO:0006355) (Fig. 5).
Pathway analysis indicated that 30 pathways corresponded to upregulated transcripts and that the most enriched network was small cell lung cancer, composed of 86 targeted genes. Moreover, pathway analysis also showed that 37 pathways corresponded to downregulated transcripts and that the most enriched network was metabolic pathways, composed of 1,193 targeted genes. Among these pathways, the gene category 'mTOR signaling pathway', is involved in the carcinogenesis of gastric cardia adenocarcinoma and gastric adenocarcinoma (15). The gene category 'gastric acid secretion' has been investigated as a cause of adenocarcinoma at the gastroesophageal junction and the distal esophagus (16). The gene category 'pathways in cancer' is involved in the development of gastric cancer (17) (Fig. 6).
Construction of the lncRNA-mRNA co-expression network
We used a coding and non-coding gene co-expression (CNC) network to evaluate the interactions among lncRNAs and identify the core regulatory lncRNAs in the network. The principle of co-expression network studies is that if the expression pattern of lncRNAs and mRNAs shows an identical or opposite curve to the other, there is likely to be some interaction between the two genes. By bioinformatics analysis, the co-expression network in GC and their paired adjacent non-tumor tissues was constructed, and by comparing the differences between them, potential core genes could be identified (18). Based on our previous results, we built CNC networks among the differentially expressed lncRNAs and mRNAs in GC specimens and their paired adjacent non-cancerous tissue specimens. The results showed network connections between 184 lncRNAs and 164 mRNAs. We combined the co-expression network, mRNAs function and microarray profile differential expression to select the key lncRNAs. We found that 14 lncRNAs and 21 mRNAs were strongly associated (Table II and Fig. 7). Most of the genes mentioned here are reported to be linked to GC.
Validation of key lncRNAs and mRNAs by quantitative RT-PCR in advanced GC
To confirm the reliability and validity of the microarray data, we selected 14 differentially expressed key lncRNAs (RP5-919F19, RP11-54O7, RP11-20I23, CTD-2541M15, AC010731, UCA1, AP000459, LOC101928316, RP11-167N4, LINC01071, RP11-111K18, TTC28-AS1, MTOR-AS1 and MEG3) and randomly selected four differentially expressed mRNAs (COL5A1, BCL2L1, COL11A2 and PLCB2), and analyzed their actual expression levels in 30 advanced GC samples and paired adjacent non-tumor tissue samples with qRT-PCR. The relative expression levels of lncRNAs and mRNAs were given as ratios of RNA to GAPDH transcript levels in the same RNA sample. We applied the paired t-test to evaluate the difference in gene expression between the tumor tissues and their adjacent non-tumor tissues. Results showed that, RP5-919F19, CTD-2541M15 and UCA1 expression was significantly higher in carcinoma than in adjacent non-tumor tissues (P<0.05) (Table III; Figs. 8 and 9). AP000459, LOC101928316, RP11-167N4 and LINC01071 expression was significantly lower in carcinoma than in adjacent non-tumor tissues (P<0.05) (Figs. 8 and 9). Expression levels of these 14 lncRNAs and mRNAs were approximately consistent with the microarray results (Table III; Figs. 10 and 11).
Validation of candidate lncRNAs in early stage GC
Then, we chose significant differential expression candidate lncRNAs (RP5-919F19, CTD-2541M15, UCA1, AP000459, LOC101928316, RP11-167N4, LINC01071 and MEG3) to further validate 20 early stage GC patient samples by qRT-PCR. Results showed that, the expression levels of CTD-2541M15 and UCA1 were significantly higher in 20 early stage GC patient tumor tissues than in their adjacent non-tumor tissues (P<0.050). AP000459, LINC01071 and MEG3 expression was significantly lower in carcinoma than in adjacent non-tumor tissues P<0.05 (Table IV). In addition, we also analyzed the change regularity of the significant differential lncRNAs both on advanced GC and early stage GC. Based on our results, compared with adjacent non-tumor tissues, the expression of CTD-2541M15 and UCA1 showed a significant gradual upward trend from early stage GC to advanced GC, AP000459, LINC01071 and MEG3 showed a significant gradual downward trend from early stage GC to advanced GC (Fig. 12).
Association between verified lncRNAs and the GC
Conditional logistic regression analysis was used to evaluate the association between differentially expressed lncRNAs and the risk of GC. As shown in Table V, significantly increased risk for GC was associated with increased expression of CTD-2541M15 and UCA1 (OR=2.860 and 30193, respectively) and reduced expression of MEG3 (OR=0.132). This suggested that CTD-2541M15 and UCA1 may function as oncogenes, while MEG3 might act as tumor suppressors.
Discussion
Despite advances in therapeutic options, early diagnosis and treatment remains the most effective way to reduce GC mortality. Even if gastrectomy is performed well, most patients with GC experience recurrence or metastasis within 1–2 years of surgery (19). The aims of cancer studies are to characterize systematically the cellular and molecular mechanisms involved in disease progression, and consequently, to identify potential biomarkers for predicting high-risk populations (20). Therefore, finding novel molecular biomarkers of malignancy has always been important for cancer prevention and treatment. Like proteins, mRNAs and miRNAs, lncRNAs show potential as novel biomarkers and therapeutic targets in different cancer types. These RNAs are >200 nucleotides and do not encode any protein. Some of them are strongly correlated with poor prognosis, suggesting a potential role in cancer progression (21). lncRNAs in GC have predominantly been reported from western countries and Japan, with few studies from China (22). Although the molecular mechanism of GC has been extensively investigated, the exact pathogenesis of this disease is still unclear (23,24). Molecular pathology of GC may vary among populations, which is likely due to differential exposure to disease risk factors including customs and habits, Helicobacter pylori variants and medical conditions. In the present study, our aim was to establish a comprehensive lncRNA and mRNA expression profile for GC in the Chinese population, a known high-risk population in Wuwei, and to investigate the mechanisms underlying GC carcinogenesis.
The present study revealed differential lncRNA and mRNA expression profiles in GC tissues and adjacent non-tumor tissues samples using microarray analysis. Here, we generated three pairs of pooled RNA samples from 10 advanced GC patients to microarray analysis. In this way, we can obtain more information of microarray detection and avoid separate specimen unique genes. In addition, using tissue samples to microarray analysis, we could not guarantee that each piece of specimen contains the same cell types and quantity, especially in tumor tissues and their adjacent non-tumor tissues. Thus, pooled RNA samples will improve the proportion of the same cell type tissues in the total RNA and make the microarray results more meaningful to population based study. In the present study, microarray results showed that lncRNA and mRNA expression levels in GC samples differed from those in adjacent non-tumor tissues. The microarray expression profile showed significantly differential expression (fold change ≥2.0) of 427 upregulated and 619 downregulated lncRNAs in 10 GC samples. In addition, there was significantly differential expression (fold change ≥2.0) of 647 upregulated and 2,349 downregulated mRNAs. GO and pathway analysis, and a lncRNA and mRNA co-expression network were used to study the biological function and potential mechanism of these genes in GC. We selected key lncRNAs and mRNAs for validation by quantitative real-time PCR in 30 pairs advanced GC and 20 early stage GC tissues and their paired adjacent non-tumor tissues. Finally, we used logistic regression analysis to find the association between candidate RNAs and the GC.
Several association studies have identified a large number of lncRNAs that are differentially expressed in disease states, including oncogenesis (25). An increasing number of studies have documented a biological link between aberrant expression of lncRNAs and GC (26). Dysregulation of lncRNAs, such as CCAT1, CHET1, H19, HOTAIR, MEG3, PVTI and SPRY4-IT1, is regarded as an important feature of GC (27–31). However, the expression and functional significance of lncRNAs in GC tumorigenesis have not been completely characterized. We also used different databases, including LncRNAdb, Gene Expression Atlas, Gene Expression Omnibus (GEO), and Array Express Databases, to analyze the expression level of the most recently cited lncRNAs in studies on GC (32,33). We found that differentially expressed lncRNAs and mRNAs, such as UCA1, MEG3, HOTAIR, LINC01071, COL5A1, BCL2L1, COL11A2 and PLCB2, were the same as in other studies of GC (34). In the present study, the GO results showed that the functions of differentially expressed mRNAs in GC and their paired adjacent non-tumor tissues were DNA-dependent transcription, proteolysis, mitotic cell cycle, transmembrane transport and regulation of transcription. Abnormal expression of these regulatory genes is closely related to tumor development (35,36). Furthermore, KEGG pathway analysis revealed 67 pathways that may play key roles in the different core epigenetic mechanisms of GC, including the mTOR signaling pathway, gastric acid secretion, pathways in cancer, apoptosis and transcriptional regulation in cancer. Li et al (17) and Fan et al (37) have characterized differentially expressed genes that are involved in pathways associated with GC. They found that BCL2A1, LAMA4, ICM1 and LAMC2 participated in cancer-related signaling pathways such as mTOR and gastric acid secretion.
To gain insight into the function of lncRNAs, we used a lncRNA-mRNA co-expression network, which may be used for predicting target genes of lncRNAs. The principle of co-expression networks is that, by comparing the similarity or difference in the expression pattern of lncRNAs and mRNAs, any correlations between them can be determined (38,39). With this method, we compiled a list of lncRNAs and mRNAs that contained some key transcripts closely related to the pathogenesis of GC. We combined the results of microarray analysis, gene function and the co-expression network. The differential expression of mRNAs that participated in cancer-related pathways and related lncRNAs were RP5-919F19, RP11-54O7, RP11-20I23, CTD-2541M15, AC010731, UCA1, AP000459, LOC101928316, RP11-167N4, LINC01071, RP11-111K18, TTC28-AS1, MTOR-AS1 and MEG33. The LncRNAdb and LncRNA Diseases database showed that UCA1 is associated with urothelial cancer (40) and MEG3 with GC (41). The function of other lncRNAs was not reported.
Then, we analyzed the expression of 14 lncRNAs and four mRNAs in 30 pairs of tumor tissues and adjacent non-tumor tissues by using qRT-PCR. After qRT-PCR verification, we demonstrated for the first time that RP5-919F19, CTD-2541M15, AC010731 and UCA1 were upregulated in tumor tissues compared with adjacent non-tumor tissues, whereas AP000459, LOC101928316, RP11-167N4 and LINC01071 were downregulated in tumor tissues. Moreover, expression of seven lncRNAs was significantly correlated with GC-related mRNAs, suggesting these lncRNAs may play important roles in GC by regulating their target genes. Although the other indicators of the lncRNA and mRNA expression results did not differ significantly, the increases and decreases were approximately consistent with the microarray results. In the next step, we chose significant differential expression candidate lncRNAs to further validate 20 early stage GC patients samples. Results suggest that CTD-2541M15 and UCA1 was significantly higher expressed, AP000459, LINC01071 and MEG3 expression was significantly lower in the 20 early stage GC patient tumor tissues than in their adjacent non-tumor tissues. In addition, these lncRNAs expression levels appear with gradual upward or downward trends, respectively. Other indicators did not differ significantly in qRT-PCR, the reason may be due to the tumor tissues in this study, which were analysed by matching with adjacent non-tumor tissues instead of normal tissues, some early lesions might have been confused with normal tissues, which led to some differences in lncRNA and mRNA expression levels between the tumor and non-tumor tissues. Furthermore, we found that expression of CTD-2541M15, UCA1, AP000459, LINC01071 and MEG3 were sensitive between tumor and adjacent non-tumor tissues. Zheng et al (42) also found that the levels of UCA1 in gastric juice from GC patients were significantly higher than in normal subjects. In addition, Sun et al (30) found that MEG3 levels were markedly decreased in GC tissues compared with adjacent normal tissues. MEG3 expression level was significantly correlated with TNM stage, depth of invasion and tumor size. However, functions of CTD-2541M15, AP000459 and LINC01071 have not been reported. The conditional logistic regression analysis results showed that CTD-2541M15 and UCA1 may function as carcinogenic genes, while MEG3 might act as tumor suppressor, but the mechanisms of gene regulation are unclear. In the next step we will launch cell validation and gene intervention experiments to further study related functions.
In conclusion, the present study provided preliminary data that have helped to understand the potential mechanisms of GC, via differentially expressed lncRNAs and mRNAs. There was significant differential expression of lncRNAs between the GC and paired non-cancerous tissues, suggesting that lncRNAs may play important roles in the tumorigenesis of GC. Further validating these significant differential expression candidate RNAs we found that CTD-2541M15, UCA1, AP000459, LINC01071 and MEG3 were closely related to the progression and development of GC. Our results suggest that these differentially expressed lncRNAs may be potential biomarkers for early diagnosis of GC, which is an interesting direction for further research.