Microarray Expression Profiling of Long Non-Coding RNAs Involved in Nasopharyngeal Carcinoma Metastasis

Increasing evidence has demonstrated a significant role for long non-coding RNAs (lncRNAs) in tumorigenesis. However, their functions in nasopharyngeal carcinoma (NPC) metastasis remain largely unknown. In this study, a model comparing high and low metastatic NPC cell lines (5-8F vs. 6-10B and S18 vs. S26) was constructed to determine the expression profile of lncRNAs using the microarray analysis, and we found 167 lncRNAs and 209 mRNAs were differentially expressed. Bioinformatic analysis indicated that the dysregulated mRNAs participated in important biological regulatory functions in NPC. Validation of 26 significantly dysregulated lncRNAs by qRT-PCR showed the expression patterns of 22 lncRNAs were in accordance with the microarray data. Furthermore, the expression level of ENST00000470135, which was the most upregulated lncRNA in high metastatic cell lines, was significantly higher in NPC cell lines and tissues with lymph node metastasis (LNM) and knocking down ENST00000470135 suppressed the migration, invasion and proliferation of NPC cells in vitro. In conclusion, our study revealed expression patterns of lncRNAs in NPC metastasis. The dysregulated lncRNAs may act as novel biomarkers and therapeutic targets for NPC.


Introduction
Nasopharyngeal carcinoma (NPC) is one of the leading head and neck malignancies and prevalent in Southeast China, Southeast Asia, the Middle East, Northeast Africa and Alaska [1][2][3][4][5]. Although local and regional control of NPC has been improved after the introduction of intensity-modulated radiation therapy and concurrent chemoradiotherapy, the prognosis of patients with distant metastasis still remains very poor [6,7]. Therefore, identification of the molecular mechanisms involved in the metastasis of NPC is urgently required to develop individualized treatment of this disease.
Metastasis is an important characteristic of malignant tumors and the leading cause of cancer-related deaths [8]. Recently, long non-coding RNAs (lncRNAs), which are defined as RNAs longer than 200 nucleotides without protein-coding potential [9], have been found to regulate cell growth, metastasis and apoptosis [10][11][12] in most tumor types via controlling transcriptional, post-transcriptional and/or epigenetic mechanisms [13][14][15]. To date, several studies have shown that a number of lncRNAs are involved in the progression of NPC, including NEAT1, HNF1A-AS, MALAT1, HOTAIR and LINC00312 [16][17][18][19][20]. However, despite this progress, the precise functions and mechanisms of lncRNAs in the metastasis of NPC still remain unclear.
Two sets of cell lines, 5-8F (high tumorigenic and metastatic potential) and 6-10B (low tumorigenic and metastatic potential), which were derived from the NPC cell line SUNE-1 [21], and S18 (high tumorigenic and metastatic potential) and S26 (low tumorigenic and metastatic potential), which were derived from CNE-2 cells [22,23], represent good models for investigating the metastasis of NPC. For the first time, we compared these high metastatic potential and low metastatic potential cell lines (5-8F vs. 6-10B and S18 vs. S26) using microarrays to identify dysregulated lncRNAs that participate in NPC tumorigenesis. Gene Ontology (GO) and pathway analysis was performed to better understand the differentially expressed mRNAs. Subsequently, we validated the results via quantitative reverse transcription polymerase chain reaction (qRT-PCR) and selected the most differentially expressed lncRNA ENST00000470135, which had previously not been reported in order to clarify its biological function. Our research provides a new insight into metastasis of nasopharyngeal carcinoma.

Profiles of Differentially Expressed Long Non-Coding RNAs (lncRNAs) and mRNAs
We conducted lncRNA and mRNA expression profiling of two groups of cell lines (5-8F vs. 6-10B and S18 vs. S26) via microarray analysis. In total, 167 lncRNAs (94 upregulated and 73 downregulated; Table S1) and 209 mRNAs (162 upregulated and 47 downregulated; Table S2) were found to be differentially expressed in high metastatic cells when compared with low metastatic cells (fold change ≥2 in both groups, normalized intensity of at least one cell line of each group ≥5). The lncRNAs are collected from highly respected databases, including Gencode, RefSeq, UCSC Knowngene and five other high quality publications ( Figure 1A). To clarify the expression signatures of the dysregulated lncRNA, we analyzed upregulated or downregulated lncRNAs identified in high metastatic potential cell lines according to their classification and chromosome distribution ( Figure 1B,C).

GO and KEGG Pathway Analysis
To explore the potential function of the differentially expressed mRNAs in high metastatic cells when compared with low metastatic cells, GO analysis was performed to describe biological process (BP), cellular component (CC) and molecular function (MF ; Tables S3-S5). The GO terms were determined by calculating Enrichment Score (p < 0.05) as previously described [24]. In our study, the aberrantly expressed genes were mainly enriched for GO terms related to regulation of cellular component organization, wound healing and cell migration involved in biological process, and cytoplasm, extracellular region and extracellular space linked with cellular component, as well as peptidase inhibitor activity, protein binding and peptidase regulator activity in molecular function. The top ten highest and most significant GO terms are shown in Figure 2A-C.
Pathway analysis based on the Kyoto Encyclopedia of Genes and Genomes (KEGG) database was used, and we identified a total of 26 pathways with significant differences (p < 0.05) in gene expression between the high metastatic potential and low metastatic potential cell lines (Table S6). The pathway terms of top ten highest Enrichment Scores are shown in Figure 2D; a number of these pathways, including the apoptosis pathway and small cell lung cancer pathway, are associated with carcinogenesis.

Validation of Significantly Dysregulated lncRNAs by qRT-PCR
Among the aberrantly expressed lncRNAs, 26 were significantly dysregulated (fold change >5 in both groups; Table 1). In order to verify the microarray data, we selected the 26 most significantly dysregulated lncRNAs (fold change >5 in both groups), which included 15 upregulated lncRNAs and 11 downregulated lncRNAs and then validated their expression level by quantitative RT-PCR (qRT-PCR) in two sets of NPC cells (5-8F vs. 6-10B and S18 vs. S26). The results showed that the expression patterns of 22 lncRNAs were consistent with the microarray data ( Figure 3A,B), which demonstrated the reliability of the microarray data. Among the 22 validated lncRNAs, the most differentially expressed lncRNA was ENST00000470135 (fold change >60 in both groups).

ENST00000470135 Is Upregulated in Nasopharyngeal Carcinoma (NPC) Cells and Tissues with Lymph Node Metastasis
To validate the importance of ENST00000470135 in NPC, we firstly examined the expression levels of ENST00000470135 in the immortalized nasopharyngeal epithelial cell line NP69 and ten NPC cell lines using qRT-PCR. The RNA level of ENST00000470135 was remarkably higher in all of the NPC cell lines tested ( Figure 4A). Moreover, we analyzed the expression of ENST00000470135 in 16 freshly frozen NPC tissues (six without lymph node metastasis (LNM) and 10 with LNM), and found that ENST00000470135 was significantly upregulated in tumors from patients with lymph node metastasis compared to those from patients without lymph node metastasis ( Figure 4B; p = 0.033). These results strongly suggest that ENST00000470135 is upregulated in NPC.

Depletion of ENST00000470135 Has Significant Effect on NPC Cell Migration, Invasion and Proliferation In Vitro
To assess whether aberrant expression of ENST00000470135 affects the motility and invasion ability of NPC cells, 5-8F and HNE-1 cells were transiently transfected with siRNA targeting ENST00000470135 or Ctrl siRNA ( Figure 5A). In the Transwell migration and invasion assays, the migratory and invasive ability of 5-8F and HNE-1 cells transfected with ENST00000470135 siRNA was significantly lower than negative control cells ( Figure 5B,C; * p < 0.05, ** p < 0.001). The results suggest that the knockdown of ENST00000470135 dramatically suppresses the migration and invasion of NPC cells.
Colony formation assay and 3-(4,5)-dimethylthiahiazo (-z-y1)-3,5-di-phenytetrazoliumromide (MTT) assay were performed to further explore whether depletion of ENST00000470135 affects the viability and proliferation of NPC cells. The colony formation rate was significantly lower in 5-8F and HNE-1 cells transiently transfected with ENST00000470135 siRNA than in cells transfected with respective control (Figure 5D; * p < 0.05, ** p < 0.001). Moreover, 5-8F and HNE-1 cells transfected with ENST00000470135 siRNA displayed significant growth inhibition ( Figure 5E; * p < 0.05, ** p < 0.001). These in vitro biological assays indicate that depletion of ENST00000470135 decreases the ability of NPC cells to migrate, invade and proliferate, providing strong evidence that ENST00000470135 may promote carcinogenesis and metastasis in NPC and further confirming the reliability of the microarray data.

Depletion of ENST00000470135 Has Significant Effect on NPC Cell Migration, Invasion and Proliferation In Vitro
To assess whether aberrant expression of ENST00000470135 affects the motility and invasion ability of NPC cells, 5-8F and HNE-1 cells were transiently transfected with siRNA targeting ENST00000470135 or Ctrl siRNA ( Figure 5A). In the Transwell migration and invasion assays, the migratory and invasive ability of 5-8F and HNE-1 cells transfected with ENST00000470135 siRNA was significantly lower than negative control cells ( Figure 5B,C; * p < 0.05, ** p < 0.001). The results suggest that the knockdown of ENST00000470135 dramatically suppresses the migration and invasion of NPC cells.
Colony formation assay and 3-(4,5)-dimethylthiahiazo (-z-y1)-3,5-di-phenytetrazoliumromide (MTT) assay were performed to further explore whether depletion of ENST00000470135 affects the viability and proliferation of NPC cells. The colony formation rate was significantly lower in 5-8F and HNE-1 cells transiently transfected with ENST00000470135 siRNA than in cells transfected with respective control ( Figure. 5D; * p < 0.05, ** p < 0.001). Moreover, 5-8F and HNE-1 cells transfected with ENST00000470135 siRNA displayed significant growth inhibition ( Figure 5E; * p < 0.05, ** p < 0.001). These in vitro biological assays indicate that depletion of ENST00000470135 decreases the ability of NPC cells to migrate, invade and proliferate, providing strong evidence that ENST00000470135 may promote carcinogenesis and metastasis in NPC and further confirming the reliability of the microarray data.

Discussion
The highest incidence of nasopharyngeal carcinoma worldwide is observed in Southern Asia; the age-standardized incidence per 100,000 males varies from 20-50 in southern China to 0.5 in white populations [1,25,26]. The remarkable geographic and racial distribution of NPC indicates that the pathogenesis and development of this cancer may be associated with genetic factors. It is well recognized that distant metastasis is the major pattern of failure in NPC. In this study, we established a model to compare high metastatic potential cell lines with low metastatic potential cell lines (5-8F vs. 6-10B and S18 vs. S26) in order to further explore the mechanisms that regulate metastasis in NPC.
Increasing evidence shows that lncRNAs play important roles in carcinogenesis, tumor progression and metastasis [11,27]. In our study, high-throughput microarray analysis was used to profile aberrantly expressed lncRNAs and mRNAs that might be involved in NPC metastasis and 167 lncRNAs and 209 mRNAs were found to be dysregulated. Recently, several studies confirm that some lncRNAs identified in other cancer types are also involved in NPC. Li et al. [28] showed that H19 (a long non-coding RNA) expression was significantly upregulated in NPC tissues and cell lines. Additionally, they found that H19 promoted invasion of NPC cells via the miR-630/EZH2 pathway. Jin et al. [18] concluded that upregulation of MALAT1 in NPC increased radioresistance by modulating miR-1/slug axis. Bo et al. [29] showed that AFAP1-AS1 promoted NPC cell metastasis via regulation of actin filament integrity. There are also several studies that have clarified the expression patterns of lncRNAs in NPC focusing on paired analyses of different types of tissues [30][31][32]. However, few studies have explored the relationships between dysregulated expression of lncRNAs and metastasis in NPC, especially using models based on cell lines.
GO analysis predicted that the dysregulated mRNAs were linked to biological process, cellular component and molecular function in NPC. The GO terms such as cell migration and cell motility involved in biological process indicated the associated gene product contributes to the metastasis of NPC. The results of pathway analysis showed that aberrantly expressed mRNAs were related to 26 signalling pathways in NPC. The correlation between these pathways and multiple diseases including NPC has been proved by previous studies. Among which, apoptosis [33][34][35], p53 [36,37] and NF-κB [38,39] have been reported to be closely related to NPC pathogenesis.
In validation of the reliability of the microarray data, strong correlations were observed between the expression levels of representative lncRNAs in both the microarray and qRT-PCR analyses. To confirm whether the lncRNAs identified by the microarray contribute to the progression of NPC, we explored the biological function of ENST00000470135, the most aberrantly expressed lncRNA. ENST00000470135 is a non-coding RNA from Ensembl [40] with a length of 673 bp, and it is spliced from RP5-884M6.1. This lncRNA has not been reported in the literature, so its biological function and regulation mechanism in nasopharyngeal carcinoma is still unclear. In this study, we found that ENST00000470135 was expressed at a significantly higher level in the NPC cell lines than the nasopharyngeal epithelial cell line NP69. Additionally, the expression of this lncRNA was upregulated in NPC tissues taken from patients with LNM compared to those from patients without LNM. Furthermore, in vitro functional experiments, including the Transwell, colony formation and MTT assays, indicated that knocking down ENST00000470135 suppresses the migration, invasion and proliferation of NPC cells. Above all, these results suggest that ENST00000470135 plays an important role in NPC metastasis, which is consistent with the findings of the microarray analysis.
Additionally, among the various kinds of biomarkers, circulating lncRNA is ideal due to its accessibility and noninvasiveness. Recent studies have demonstrated that circulating lncRNAs could forecast prognosis and predict therapeutic efficacy in different kinds of cancers [41][42][43]. For example, three circulating lncRNAs (LincRNA-p21, GAS 5, HOTAIR) have been discovered as biomarkers to predict chemoradiotherapy sensitivity in head and neck cancers including nasopharyngeal carcinoma [44], indicating the important value of circulating lncRNAs in clinical outcome predictions.
A total of 16 freshly frozen NPC tissue samples (from 10 patients with LNM and 6 patients without LNM) were collected at Sun Yat-sen University Cancer Center. Written informed consent was obtained from each patient before biopsy. This study was approved by the Institutional Ethical Review Board of our Cancer Center (GZR2016-108, 23 February 2016).

RNA Extraction and Quality Control
Total RNA was extracted from the 16 frozen NPC tissue specimens and 11 cell lines using TRIzol reagent (Invitrogen, Grand Island, NY, USA) according to the manufacturer's instructions. The quality and amount of RNA were assessed using a NanoDrop ND-2000 spectrophotometer (Thermo Scientific, Rockford, IL, USA) and RNA integrity was assessed by standard denaturing agarose gel electrophoresis. Isolated RNAs were stored at −80 • C prior to lncRNAs microarray analysis and quantitative reverse transcription PCR.

Microarray Analysis
Arraystar Human LncRNA Microarray V3.0 (ArrayStar, Rockville, MD, USA) which contains 30,586 lncRNAs and 26,109 coding transcripts was performed to detect the genome-wide profiling compared the high metastatic potential (5-8F and S18) and low metastatic potential cell lines (6-10B and S26). The lncRNAs were collected from the majority of landmark public databases (Gencode, Refseq, ect.) as well as high-quality publications. In order to recognize individual transcripts accurately, we used a splice junction or specific exon probe to represent each transcript. mRNA was purified with the extraction of rRNA according to the manufacturer's instructions (mRNA-ONLYTM eukaryotic mRNA Isolation Kit, Epicentre, Madison, WI, USA). Each of the samples was amplified and transcribed into fluorescent cRNA without 3 bias. After the labeled cRNAs were measured by NanoDrop ND-1000, blocking agent and fragmentation buffer were added. Then, the mixture was heated and the cRNA was diluted using the hybridization buffer. Finally, the labeled cRNAs were assembled to the human lncRNA array version 3.0 (8 × 60 K; Arraystar, Rockville, MD, USA). After washing and fixation, the hybridized arrays were scanned using an Agilent DNA scanner G2505C (Agilent Technologies, Santa Clara, CA, USA). KangChen Bio-tech (Shanghai, China) perfomed the microarray experiments.

Bioinformatics Data Analysis and Data Mining
Acquired array images were analyzed by Agilent Feature Extraction software (v11.0.1.1). GeneSpring GX v12.1 software package (Agilent Technologies, Santa Clara, CA, USA) was used to perform quantile normalization and subsequent data processing. Differentially expressed lncRNAs and mRNAs between the high and low metastatic potential cell lines (5-8F vs. 6-10B and S18 vs. S26) were identified through fold change and normalized intensities filtering (fold change ≥2 in both groups, normalized intensity of at least one cell line of each group ≥5). We have deposited the microarray data in the National Center for Biotechnology Information's Gene Expression Omnibus (GSE89804).
The aberrantly expressed mRNAs in NPC metastasis were selected to perform GO and KEGG pathway anaysis. For GO analysis (Available online: http://geneontology.org/), the corresponding genes were divided into 3 classifications by enrichment analysis, including biological process (BP), cellular component (CC) and molecular function (MF). Using the latest KEGG (Available online: http://www.genome.jp/kegg/), we analyzed the differentially expressed mRNAs and calculated the enrichment of different pathways.

Reverse Transcription and Quantitative Real-Time RT-PCR
To measure the expression of selected lncRNAs, we firstly performed reverse transcription using random primers (Promega, Madison, WI, USA) and Moloney Murine Leukemia Virus (M-MLV) reverse transcriptase (Promega, Madison, WI, USA) and stored the RT products briefly at 4 • C or at −20 • C until use. Then, quantitative RT-PCR was conducted on a CFX96 Touch sequence detection system (Bio-Rad, Hercules, CA, USA) using SYBR Green (Platinum SYBR Green qPCR SuperMix-UDG reagents; Invitrogen, Grand Island, NY, USA). GAPDH was used as the normalization control and the relative expression levels were calculated using the 2 −∆∆Ct method [45]. The real-time RT-PCR primers for the lncRNAs and GAPDH are shown in Table S7.

Transwell Migration and Invasion Assays
Cell migration and invasion ability were assessed using Transwell chambers (8-µm pores; Corning) coated with or without Matrigel (BD Biosciences, San Diego, CA, USA). The plates were placed in the culture incubator for at least 30 min at 37 • C, and then 5-8F or HNE-1 cells (5 × 10 4 or 1 × 10 5 ) were transfected and suspended in 200 µL serum-free medium and placed into the upper chambers. The lower chambers were filled with 500 µL medium supplemented with 10% FBS. After incubation for 12-24 h, the NPC cells that had migrated or invaded through the Matrigel and pores were fixed with 4% paraformaldehyde, stained using 0.5% crystal violet and manually counted using an inverted microscope (100×).

MTT Assay and Colony Formation Assay
For the 3-(4,5)-dimethylthiahiazo(-z-y1)-3,5-di-phenytetrazoliumromide (MTT) assay, 5-8F and HNE-1 cells transfected with siRNAs were seeded at a density of 1000 cells per well in 96-well plates. At 1, 2, 3, 4 and 5 days, an ELX800 spectrophotometric plate reader (Bio-Tek, Winooski, VT, USA) was used to measure the cell viability at 490 nm. For the colony formation assay, the transfected cells were plated at 400 cells per well in 6-well plates and cultured in the cell culture incubator for 7 to 12 days at 37 • C, fixed with 4% paraformaldehyde, stained using 0.5% crystal violet, and then the numbers of colonies were counted.

Statistical Analysis
Significance test was used for GO analysis and Fisher's exact test was used for pathway analysis. For qRT-PCR validation analysis and functional analysis, the Student's t-test was used to determine the significance of the differences between two groups. p-Values < 0.05 were considered significant. Statistical analysis was performed using SPSS software version 16.0 (IBM, Chicago, IL, USA).

Conclusions
In conclusion, we described the expression profile of lncRNAs and mRNAs in NPC metastasis by microarray and found several hundred dysregulated lncRNAs and mRNAs. After qRT-PCR verification, we explored the biological function of the most abnormally expressed lncRNA, ENST00000470135. The results of in vitro functional experiments suggest that the lncRNA promote migration, invasion and proliferation of NPC cells. This investigation to elucidate the roles of lncRNAs in NPC metastasis provides novel insight into the molecular mechanisms of disease progression, which may have significance for NPC therapy. In the near future, we intend to study the value of the dysregulated lncRNAs identified in this study as prognostic biomarkers using a large number of patient samples, and also explore their biological function and regulation mechanism, with the overall objective of providing novel and effective therapeutic strategies for NPC.