Microarray data re-annotation reveals specific lncRNAs and their potential functions in non-small cell lung cancer subtypes.

Non‑small‑cell lung cancer (NSCLC) is a leading cause of cancer mortality worldwide. The most common subtypes of NSCLC are adenocarcinoma (AC) and squamous cell carcinoma (SCC). However, the pathophysiological mechanisms contributing to AC and SCC are still largely unknown, especially the roles of long non‑coding RNAs (lncRNAs). The present study identified differentially expressed lncRNAs between lung AC and SCC by re‑annotation of NSCLC microarray data analysis profiling. The potential functions of lncRNAs were predicted by using coding‑non‑coding gene co‑expressing network. Reverse transcription-quantitative polymerase chain reaction (RT‑qPCR) was used to investigate lncRNA expression levels in AC cell lines (A549 and L78), SCC cell lines (H226 and H520) and normal cells (NL‑20). Western blotting analysis was used to investigate the protein expression levels in these cell lines. A total of 65 lncRNAs were differentially expressed between AC and SCC including 28 lncRNAs that were downregulated in SCC subtypes compared with those in AC ones, and 37 upregulated lncRNAs in SCC subtypes compared with AC subtypes. Three lncRNAs, sex determining region Y‑box 2 overlapping transcript (SOX2‑OT), NCBP2 antisense RNA 2 (NCBP2‑AS2) and ubiquitin like with PHD and ring finger domains 1 (UHRF1), were predicted to be associated with lung cancer; RT‑qPCR confirmed that SOX2‑OT and NCBP2‑AS2 were associated with lung cancer. Finally, western blot assays demonstrated that there was no difference in β‑catenin and glycogen synthase kinase 3β (GSK‑3β) expression in cancer cells compared with NL‑20, but increased phosphorylated (p‑)β‑catenin and p‑GSK‑3β was detected in lung cancer cell lines compared with NL‑20, particularly in A549 cells. Although these results require further experimental verification, the analysis of lncRNA signatures between AC and SCC has provided insights into the regulatory mechanism of NSCLC development.

Introduction

Non-small-cell lung cancer (NSCLC) is the leading cause of cancer death worldwide (1). The most common subtypes of NSCLC are adenocarcinoma (AC) and squamous cell carcinoma (SCC) (2). Lung AC is the most common histological type of NSCLC, resulting in >500,000 deaths globally every year (3). In addition, SCC comprises 25–30% of pulmonary cancers (4). Despite advances in surgery, molecular subtyping and targeted therapy, prognosis of AC and NSCLC remains poor (5). However, the pathophysiological mechanisms contributing to AC and SCC remain largely unknown. Thus, there is an urgent need for the identification of novel agents that will provide clinicians with useful assistance in improving the survival rate of patients with AC and SCC.

Long non-coding RNAs (lncRNAs) are a large and important class of heterogeneous ncRNAs with a length >200 nucleotides (6,7). LncRNAs are emerging contributors to the cancer paradigm, demonstrating potential roles in both oncogenic and tumor suppressive pathways (8). The human genome is replete with lncRNA, many of which are transcribed and likely to have a functional role. A number of dysregulated lncRNAs have been detected in multiple human cancers (9–11), and their expression is associated with cancer metastasis and prognosis. For example, the lncRNA Hox transcript antisense intergenic RNA (HOTAIR) is upregulated in primary breast tumors and metastases, and its overexpression is associated with enhanced breast cancer metastasis (12). Notably, lncRNA dysregulation was proposed as a hallmark feature in cancer (10). Recently, several lncRNA signatures were developed as novel predictors of survival in patients with cancer (13,14), displaying a similar prognostic power to that of protein-coding RNA and microRNAs, thus providing a new molecular option for cancer diagnosis and prognosis. However, the survival rate of lncRNA signatures for lung cancer patients has not yet been investigated and the lncRNA-based pathophysiological mechanism contributing to AC and SCC have not been well characterized.

To investigate the underlying role of lncRNA in lung AC and SCC, the present study re-annotated a NSCLC microarray data analysis profiling to identify differentially expressed lncRNAs, then predicted the potential functions of these lncRNAs. Finally, some lncRNAs were validated using reverse transcription-quantitative polymerase chain reaction (RT-qPCR). A total of 65 lncRNAs were differentially expressed between AC and SCC. Of these lncRNAs, 3 lncRNAs were potentially related to cancer-associated pathways. Two of these lncRNAs were verified using RT-qPCR. Although these results require further experimental verification, the analysis of lncRNA signatures between AC and SCC provided interesting insights into the regulated mechanism of development of NSCLC.

Materials and methods

Cell culture and agents

The original microarray data were obtained from tumor tissues of 40 AC and 18 SCC samples in NSCLC (https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE10245) (15). NSCLC cell lines A549, L78, H226, H520 and normal lung cell line NL-20 were obtained from Shanghai Institutes for Biological Sciences, Chinese Academy of Sciences (Shanghai, China). The cells were routinely cultured in RPMI-1640 (Invitrogen; Thermo Fisher Scientific, Inc., Waltham, MA, USA, Carlsbad USA) supplemented with 10% heat-inactivated fetal bovine serum (FBS; Invitrogen; Thermo Fisher Scientific, Inc.), 100 U/ml penicillin and 100 µg/ml streptomycin in a humidified cell incubator with an atmosphere of 5% CO2 at 37°C. Cells growing at an exponential rate were used for the experiments. Rabbit monoclonal glycogen synthase kinase-3β (GSK-3β; ab32391), phosphorylated (p-)GSK-3β (ab131097), β-catenin (ab16051), p-β-catenin (ab73153) and β-actin (ab8227), were obtained from Abcam (Cambridge, UK). Infrared (IR) fluorescent dye-labeled secondary antibody (A32731) was obtained from Thermo Fisher Scientific, Inc. Protease inhibitor mixture was obtained from Roche Diagnostics (Basel, Switzerland).

Microarray and Functional Reannotation

The microarray data set GSE10245 was downloaded from the Gene Expression Omnibus (GEO) database (http://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE10245). The dataset was composed of the histological subtype AC samples (n=40) and SCC samples (n=18) using Human Genome U133 (HG-U133) Plus 2.0 array (Affymetrix, Inc., Santa Clara, CA, USA).

Functional Reannotation of long non-coding RNAs

The probes on the HG-U133 Plus 2.0 array were reannotated for human lncRNA using non-coding RNA function annotation server (ncFANs) as described in its website (http://www.ebiomed.org/ncFANs/) (16). Differentially expressed lncRNAs were selected using the fold change method: lncRNAs with a fold change value >2 or <0.5 were considered as differentially expressed lncRNAs. In the present study, AC was considered as ‘control’.

Obtaining co-expressing genes and functional enrichment for differentially expressed lncRNA

To explore potential functions of differentially expressed lncRNAs, co-expressing genes for each lncRNA were identified by calculating Pearson's correlation coefficient (PCC) between each differentially expressed lncRNA and all genes across all 18 samples. The genes with a strict cut-off of PCC>0.75 or <-0.75) were identified as co-expressing genes. Pathway enrichment was then implemented to identify the affected pathways of lncRNA co-expressing genes by using the Database for Annotation, Visualization and Integrated Discovery (DAVID) version 6.7 (17,18).

RT-qPCR

Total RNA was extracted using TRIzol (Invitrogen; Thermo Fisher Scientific, Inc.) according to the manufacturer's protocol and diluted to 200 ng/ml. RT-qPCR was then performed using One Step SYBR® Prime Script™ RT-PCR Kit II (Takara Biotechnology Co., Ltd., Dalian, China) according to the manufacturer's protocol. To be brief, two stages were included: Stage I for reverse transcription, 42°C 5 min, 95°C 10 sec; stage II for qPCR reaction, denaturation: 95°C 5 sec, annealing and extension: 60°C 34 sec, repeated 35 cycles. All data were collected from three independent experiments. GAPDH was used as an internal control and relative gene expression was calculated using the 2−ΔΔCq method (19). Primers used in RT-qPCR were as follows: Sex determining region Y-box 2 overlapping transcript (SOX2-OT), forward 5′-CGAAATGGATTCACGGTGCC-3′ and reverse 5′-TGCCAGATCAGGGTGTTGTC-3′; NCBP2 antisense RNA 2 (NCBP2-AS2), forward 5′-TTTAATCCTGGGCTGTGCGG-3′ and reverse 5′-CACAAGACTTCAACGGAGCG-3′; and GAP DH, forward 5′-GGACCAATACGACCAAATCCG-3′ and reverse 5′-AGCCACATCGCTCAGACAC-3′.

Protein extraction

Whole cell proteins were extracted using RIPA buffer (P0013C; Beyotime Institute of Biotechnology, Shanghai, China) according to the manufacturer's protocol. In brief, cells were washed with PBS and scraped with RIPA buffer. Cells were collected in microfuge tubes and agitated for 30 min at 4°C. Then, centrifugation was performed at 16,000 × g for 20 min at 4°C. The supernatants were collected in a fresh tube for western blot analysis. Protein concentration was measured using Bradford assay (P0006; Beyotime Institute of Biotechnology).

Western blot analysis

A total of 50 µg proteins from the cells were separated by 10–12% SDS-PAGE and transferred to nitrocellulose membrane (Pall Life Sciences, Port Washington, NY, USA). Membranes were blocked with 5% milk at room temperature for 30 min, then incubated overnight at 4°C with β-catenin (1:1,000 dilution), p-β-catenin (1:1,000 dilution), p-GSK-3β (1:500 dilution), GSK-3β (1:500 dilution) or β-actin (1:1,000 dilution) antibodies. Infrared (IR) fluorescent dye-labeled secondary antibody (1:10,000) was incubated with the membrane at room temperature for 1 h. Western blot bands were visualized using an IR Imaging System (LI-COR Biosciences, Lincoln, NE, USA). Band densities for each group were quantified using Odyssey software version 3.0 (LI-COR Biosciences) and normalized to β-actin expression.

Statistical analysis

Data were obtained from 3–6 independent experiments and are presented as the mean ± standard deviation. Data were evaluated by analysis of variance followed by Student-Newman-Keuls post hoc tests using the SPSS software v.19.0 (IBM Corp., Armonk, NY, USA) and illustrated using GraphPad Prism version 5.0 (GraphPad Software, Inc., La Jolla, CA, USA). P<0.05 was considered to indicate a statistically significant difference.

Results and Discussion

Differentially expressed lncRNAs between AC and SCC subtypes

Following functional annotation using ncFANs and calculating fold change value, mRNA and lncRNA profiles were obtained. As a result of the bidirectional hierarchical cluster of lncRNA profiles, it was observed that the samples were clustered into classes approximately corresponding to AC and SCC, suggesting that the re-annotated profile could be used for subsequent analysis (Fig. 1).

Figure 1.

Bidirectional hierarchical cluster of long non-coding RNA profiles.

Following analysis, a total of 65 lncRNAs were considered differentially expressed (Table I). Of these 65 lncRNA, 28 were downregulated in the SCC subtypes compared with AC subtypes. By contrast, 37 lncRNAs were upregulated in SCC compared with AC. The most upregulated lncRNA was SOX2-OT (ENSG00000242808; Table I). SOX2-OT is located on human chromosome 3q26.33, which is frequently amplified in lung SCCs (20). Recent studies have demonstrated a positive role for it in transcriptional regulation of SOX2 (21,22). Similarly to SOX2, SOX2OT is highly expressed in embryonic stem cells and downregulated upon induction of differentiation (23). SOX2OT is dynamically regulated during the embryogenesis of vertebrates, and delimited to the brain in adult mice and human (23). Shahryari et al (23) investigated the role of SOX2OT in lung cancer and found that the expression level of SOX2-OT is significantly higher in SCCs than in ACs of the lung, which is in agreement with a previous study by this group (20). Another study indicated that SOX2-OT is a novel prognostic factor (24). In brief, SOX2-OT is important in regulating lung cancer cell proliferation and may be a potential novel prognostic indicator (20,23,25).

Table I.

Differentially expressed lncRNAs.

lncRNA	log (fold change)
ENSG00000253898	−4.152017003
ENSG00000235123	−3.540259841
ENSG00000262117	−2.746729646
ENSG00000245750	−2.61388764
ENSG00000272620	−2.564460843
ENSG00000225383	−2.358113042
ENSG00000260943	−2.270111418
ENSG00000264868	−2.21214258
ENSG00000266968	−2.202088651
ENSG00000204832	−1.931326052
ENSG00000256969	−1.767283812
ENSG00000259974	−1.678282345
ENSG00000235584	−1.633014477
ENSG00000260997	−1.563708541
ENSG00000240476	−1.555238143
ENSG00000232018	−1.515156717
ENSG00000223414	−1.454123693
ENSG00000233834	−1.366598222
ENSG00000237499	−1.278009264
ENSG00000250328	−1.221423553
ENSG00000251629	−1.205091444
ENSG00000230487	−1.177625469
ENSG00000247134	−1.155801947
ENSG00000254429	−1.074388165
ENSG00000241295	−1.053131154
ENSG00000272732	−1.042988456
ENSG00000253563	−1.015484492
ENSG00000214049	−1.000246779
ENSG00000266401	1.003846483
ENSG00000246859	1.005233819
ENSG00000231971	1.021838364
ENSG00000215808	1.034891987
ENSG00000034063	1.042316963
ENSG00000267034	1.093965748
ENSG00000245526	1.106449907
ENSG00000225484	1.124075352
ENSG00000233864	1.136818147
ENSG00000255135	1.169779461
ENSG00000270170	1.195208999
ENSG00000261668	1.208732331
ENSG00000271856	1.217547617
ENSG00000226363	1.245074712
ENSG00000273190	1.26525531
ENSG00000234715	1.278314773
ENSG00000227640	1.304396623
ENSG00000228630	1.32037165
ENSG00000261409	1.373814297
ENSG00000254560	1.389245877
ENSG00000110347	1.461512858
ENSG00000261521	1.470414154
ENSG00000269994	1.499538599
ENSG00000230937	1.584953003
ENSG00000225826	1.638461788
ENSG00000213963	1.851214762
ENSG00000249082	1.942238906
ENSG00000224259	1.961743949
ENSG00000233901	1.962040139
ENSG00000241732	1.968491957
ENSG00000104725	2.438796861
ENSG00000261116	2.754710633
ENSG00000245648	2.758743567
ENSG00000260581	2.93322784
ENSG00000265992	3.470736871
ENSG00000242512	3.647947373
ENSG00000242808	5.237044135

lncRNA, long non-coding RNA.

The most downregulated lncRNA was LINC01419 (ENSG00000253898; Table I). Although, to the best of our knowledge, there have been no reports about this lncRNA in lung cancer, Zhang et al (26) validated that LINC01419 was significantly overexpressed in hepatitis B- and C-related HCC when compared with matched non-tumor liver tissues. Furthermore, functional predictions suggested that LINC01419 regulates cell cycle genes (26).

The method of ncFANs has been applied to different cancer types and other diseases. For example, in the studies of Gao et al (27,28), using the functional annotation algorithm from ncFANs, some differentially expressed lncRNAs were identified in tongue SCC tissues (e.g. SPRR2D-1) and in nasopharyngeal carcinoma (e.g. C22orf32-1). These lncRNAs were both different from the results of the present study, indicating that this is not a common phenomenon or specific to a certain cancer type.

LncRNA-gene co-expression network and functional enrichment

Following calculation of the PCC between differentially expressed lncRNAs and all genes, significantly co-expressed genes were identified for these lncRNAs. A lncRNA-gene co-expression network was then constructed, in which nodes were lncRNAs and genes, and were connected if they are co-expressed (Fig. 2A); in this network, there were 31 lncRNA nodes. The potential biological functions of these differentially expressed lncRNAs were then explored using DAVID v.6.7 for corresponding co-expressed genes, which established that there was enrichment in some cancer-related biological processes, including the cell cycle, p53 signaling pathway and Wnt signaling pathway (Figs. 2A and 3; Table II). To better illustrate the regulated roles of lncRNAs, the co-expressed genes of lncRNA-UHRF1 were annotated into pathways of the cell cycle (Fig. 2B); potential regulation of some key proteins in this pathway was identified, including minichromosome maintenance complex component 2 (MCM2), cyclin A2 (CCNA2, also known as CycA), cyclin B1 (CCNB1, also known as CycB), checkpoint kinase 1 (CHEK1, also known as Chk1,2) and checkpoint kinase 2 (CHEK2, also known as Chk2). Zhang et al (29) investigated whether lovastatin inhibits proliferation due to MCM2 in NSCLCs. Using RT-qPCR and western blot analysis, they demonstrated that lovastatin treatment markedly and consistently inhibits the expression of MCM2. MCM2 knockdown inhibited expression of retinoblastoma protein, cyclin D1 and cyclin dependent kinase 4, but increased expression of cyclin dependent kinase inhibitor 1A (p21) and tumor protein p53, suggesting that siMCM2 triggered cell cycle arrest (29).

Figure 2.

(A) The co-expressing network of differentially expressed lncRNAs and mRNAs. (B) Cell cycle pathway annotated with the co-expressing genes of lncRNA-UHRF1. (C) The co-expressing genes of three lncRNAs. lncRNA, long non-coding RNA.

Figure 3.

The enriched pathways of differentially expressed lncRNAs. The dashed line indicates the cutoff value.

Table II.

Enriched pathways of corresponding co-expressed genes of differentially expressed long non-coding RNAs.

difLnc	Pathway-ID	Pathway name	P-value	Gene
ENSG00000034063	Path:04110	Cell cycle	2.82×10−14	1111;4173;4175;4176;5347; 699;701;890;891;9133;9232;9700
ENSG00000034063	Path:03030	DNA replication	3.09×10−9	10535;2237;4173;4175;4176;5985
ENSG00000034063	Path:04914	Progesterone-mediated oocyte maturation	1.60×10−5	5347;699;890;891;9133
ENSG00000034063	Path:04114	Oocyte meiosis	4.91×10−5	5347;6790;699;9232;9700
ENSG00000034063	Path:04115	p53 signaling pathway	9.80×10−5	1111;51512;891;9133
ENSG00000242808	Path:05217	Basal cell carcinoma	2.95×10−3	7482;8324
ENSG00000242808	Path:04916	Melanogenesis	9.45×10−3	7482;8324
ENSG00000270170	Path:04310	Wnt signaling pathway	2.44×10−3	1857;59343;79718
ENSG00000270170	Path:04144	Endocytosis	5.69×10−3	23527;5337;7037

difLnc, differentially expressed long non-coding RNAs.

CCNA2, a secreted matrix-associated molecule, is involved in multiple cellular processes (30). Previous studies have indicated that expression of CCNA2 is inversely correlated with the aggressiveness of NSCLC (31). Jim et al (32) used 3 NSCLC cell line systems to demonstrate that long-term treatment of cells with recombinant CCNA2 leads to a permanent cell cycle arrest in G1 phase, in addition to apoptosis assays to demonstrate that cells remain viable. Furthermore, CCNA2-treated NSCLCs acquire phenotypic characteristics of senescent cells, including an enlarged and flattened cell shape and expression of senescence-associated β-galactosidase (32).

ATR serine/threonine kinase, a DNA damage response kinase, and its effector kinase, CHEK1, are required for cancer cells to survive oncogene-induced replication stress (33). High expression of CHEK1 in lung tumors is associated with poor overall survival (34); CHEK1 has previously been demonstrated to be a direct target of microRNA-195, which decreases CHEK1 expression in lung cancer cells.

Notably, the present study demonstrated that 3 differentially expressed lncRNAs, ENSG00000241732, ENSG00000260581 and ENSG00000230937, co-regulated some mRNAs: The co-expressed genes of these lncRNAs are presented in solid boxes in Fig. 2C. Of these genes, late cornified envelope 3D (LCE3D) has previously been identified as significantly upregulated in head and neck cancer cells compared with normal mucosa (35). Wright et al (36) investigated the possible role of interferon regulatory factor 6 (IRF6) in asbestos lung tumorigenicity, its potential utility as a marker of asbestos related lung cancer for purposes of causal attribution and its potential as a treatment target for lung cancers arising in asbestos exposed persons. Overexpression of small proline rich protein 1B (SPRR1B) in bronchial epithelial cells is a marker for early metaplastic changes induced by various toxicants/carcinogens (37). Numerous studies have established the role of SPRR1B during squamous differentiation of skin and respiratory epithelial cells (38). Woenckhaus et al (39) investigated the gene expression in surgically resected and microdissected samples of NSCLCs, matched normal bronchial epithelium, and peripheral lung tissue from both smokers and non-smokers using an Affymetrix U133A array: They demonstrated that SPRR1B is differentially expressed and is related to cell differentiation (39).

Expression of lncRNA1 and LncRNA2 in AC and SCC cells

Finally, 2 lncRNAs were selected for further investigation: SOX2-OT and NCBP2-AS2. These lncRNAs were selected since SOX2-OT was the most upregulated lncRNA and NCBP2-AS2 potentially regulates the cancer-associated Wnt signaling pathway (Table II) (40). To explore SOX2-OT and NCBP2-AS2 expression in NSCLC cells, RT-qPCR was used to determine SOX2-OT levels in normal lung cell line NL-20 and cancer cell lines A549, L78, H226 and H520. As demonstrated in Fig. 4A, SOX2-OT was significantly upregulated in AC cells (A549 and L78) compared with the SCC cells (H226 and H520). NCBP2-AS2 levels were also significantly upregulated in AC cells compared with the SCC cells (P<0.05; Fig. 4B). In addition, levels of SOX2-OT and NCBP2-AS2 were both significantly lower in NL-20 cell lines compared with AC and SCC cells (Fig. 4A and B).

Figure 4.

Expression of long non-coding RNAs (A) SOX2-OT and (B) NCBP2-AS2 in A549, L78, H226 and H520 cells, respectively. Data are expressed as the mean ± standard error of the mean (n=3). *P<0.05 vs. normal group. SOX2-OT, sex determining region Y-box 2 overlapping transcript; NCBP2-AS2, NCBP2 antisense RNA 2.

Expression levels of Wnt signaling pathway proteins following treatment with lncRNA1 and lncRNA2 in AC and SCC cells

Deregulation of the Wnt signaling pathway is important in the initiation and development of lung cancers. The present study examined whether the Wnt signaling pathway was involved in different lung cancer cells lines. Western blot assays demonstrated that the expression of β-catenin and GSK-3β was not altered in cancer cell lines compared with NL-20, but increased p-β-catenin and p-GSK-3β was detected in lung cancer cell lines compared with NL-20, particularly A549 cells (Fig. 5). These results demonstrated the potentially different involvement of the Wnt signaling pathway in the two types of NSCLC cancer (AC and SCC).

Figure 5.

Activation of Wnt signaling pathway proteins in normal and lung cancer cells following transfection with long non-coding RNAs. (A) Phosphorylation of β-catenin following transfection with SOX2-OT. (B) Phosphorylation of GSK-3β following transfection with NCBP2-AS2. Data are expressed as the mean ± standard error of the mean (n=3). *P<0.05 vs. normal group. p-, phosphorylated; SOX2-OT, sex determining region Y-box 2 overlapping transcript; GSK-3β, glycogen synthase kinase 3β; NCBP2-AS2, NCBP2 antisense RNA 2.

In conclusion, in the present study, differentially expressed lncRNAs between lung AC and SCC were identified by re-annotation of NSCLC microarray data analysis profiling. A total of 65 lncRNAs were differentially expressed between AC and SCC. The potential functions of these lncRNAs were then predicted. Finally, 2 of these lncRNAs were investigated using RT-qPCR. Although these results require further experimental verification, the analysis of lncRNA signatures between AC and SCC still provided insights into the regulated mechanism of development of NSCLC.