Long noncoding RNA expression signatures of bladder cancer revealed by microarray

Dysregulation of long noncoding RNAs (lncRNAs) has been regarded as a primary feature of several human cancers. However, the genome-wide expression and functional significance of lncRNAs in bladder cancer remains unclear. The aim of this study was to identify aberrantly expressed lncRNAs that may play an important role in contributing to bladder cancer pathogenesis. In this study, we described lncRNAs profiles in four pairs of human bladder cancer and matched normal bladder tissues by microarray. We finally determined 3,324 differentially expressed human lncRNAs and 2,120 differentially expressed mRNAs (≥2-fold change). A total of 110 lncRNAs were significantly differentially expressed between the tumor and the control groups (≥8-fold change). Four lncRNAs (TNXA, CTA-134P22.2, CTC-276P9.1 and KRT19P3) were selected for further confirmation of microarray results using quantitative PCR (qPCR), and a strong correlation was identified between the qPCR results and microarray data. We also observed that numerous lncRNA expression levels were significantly correlated with the expression of tens of protein coding genes by construction of the lncRNA-mRNA co-expression network. Kyoto Encyclopedia of Genes and Genomes annotation showed a significant association with p53, bladder cancer, cell cycle and propanoate metabolism pathway gene expression in the bladder cancer group compared with the normal tissue group, indicating that deregulated lncRNAs may act by regulating protein-coding genes in these pathways. We demonstrated the expression profiles of human lncRNAs in bladder cancer by microarray. We identified a collection of aberrantly expressed lncRNAs in bladder cancer compared with matched normal tissue. It is likely that these deregulated lncRNAs play a key or partial role in the development and/or progression of bladder cancer.


Introduction
Bladder cancer is the most common urinary tract malignant tumor and accounts for 5% of all diagnosed cancers (1). Although the majority of cases are not clinically advanced at presentation, the disease eventually recurs or develops metastases. Urothelial carcinoma of the bladder, the most common histopathological type of bladder cancer, has a variety of genetic and phenotypic characteristics. Numerous factors, including chromosomal anomalies, genetic polymorphisms, and genetic and epigenetic alterations, contribute to tumorigenesis and progression of urothelial carcinoma of the bladder (2). Therefore, revealing the molecular mechanism for bladder cancer development is required for developing effective therapy.
With the development of high-throughput DNA sequencing and array-based technologies, various classes of noncoding RNAs (ncRNAs) have recently been shown to function as regulators of protein-coding genes (3). Emerging data strongly suggest that long noncoding RNAs (lncRNAs; length >200 bp) are important in the basal regulation of protein coding genes, at the transcriptional and the posttranscriptional levels (3). Dysregulation of lncRNAs, such as HOTAIR, H19, MALAT-1 and PCA3, has been regarded as a primary feature of several human cancers, including prostate cancer, breast cancer, gastric cancer and hepatocellular carcinoma (4)(5)(6)(7). Certain recent studies have reported that several lncRNAs, including UCA1, MALAT-1 and ncRAN, show marked potential in the field of bladder cancer progression (8)(9)(10). However, the genome-wide expression and functional significance of lncRNAs in bladder cancer remains unclear.
In this study, we present the lncRNA expression profiles in four pairs of bladder cancer samples and matched histologically normal urothelium samples by lncRNA microarray. We observed that a collection of lncRNAs were aberrantly expressed in bladder cancer compared with matched normal tissue, several of which were evaluated by quantitative PCR (qPCR) in a total of 51 pairs of tissues.

Materials and methods
Patient samples. Fifty-six patients with urothelial carcinoma of the bladder who received radical cystectomy at the Fudan University Shanghai Cancer Center (Shanghai, China) were included in the study. Urothelial carcinoma of the bladder was diagnosed histopathologically. Written informed consent was obtained from all patients and the study was approved by the Institutional Review Board of Fudan University Shanghai Cancer Center (Shanghai, China). Of the 56 pairs of samples, five pairs were used in lncRNA microarray analysis [one pair was excluded from analysis according to three dimension principal component analysis (3D-PCA)] and 51 pairs were analyzed by qPCR. The tumor sample and matched normal bladder tissue from each subject were snap-frozen in RNA ladder immediately after resection and stored in the tissue bank.
RNA extraction. If the proportion of cancer cells in a tissue section was >80% then the frozen block was subjected to RNA extraction. Total RNA was extracted from 56 pairs of snap-frozen urothelial carcinoma and matched normal bladder tissues using TRIzol reagent (Invitrogen Life Technologies, Carlsbad, CA, USA) according to the manufacturer's instructions. The RNA integrity was evaluated by a NanoDrop ND-1000 spectrophotometer (NanoDrop products, Wilmington, DE, USA).
Microarray and computational analysis. RNA purified from total RNA following the removal of rRNA was amplified and transcribed into fluorescent cRNA along the entire length of the transcripts without 3' bias utilizing a random priming method and cDNA was labeled and hybridized to the Human LncRNA Array V2.0 (8x60 K, Arraystar Inc., Rockville, MD, USA). A total of 33,045 lncRNAs and 30,215 coding transcripts, which were collected from the most authoritative databases, including RefSeq, UCSC Known Genes, Ensembl and the associated literature, were detected by the microarray.
The microarray work was performed by KangChen Bio-tech (Shanghai, China). In brief, the Arraystar lncRNA array protocol was as follows: i) the RNA sample, kit and reagents were prepared, including TRIzol reagent, Biopulverizer (Biospec, Bartlesville, OK, USA) and Mini-Bead-Beater-16 (Biospec); ii) Total RNA Clean-up and RNA quality control; iii) labeling reaction was prepared; iv) labeled/amplified RNA and labeled cRNA QC were purified; v) hybridization was performed; vi) microarray wash was conducted; vii) scanning was performed; and viii) the data were extracted using Agilent feature extraction software (Agilent Technologies, Santa Clara, CA, USA).
The arrays were scanned by the Agilent Scanner G2505B (Agilent Technologies) and the acquired array images were analyzed by Agilent Feature Extraction software (version 10.7.3.1; Agilent Technologies). Quantile normalization and subsequent data processing were performed using the GeneSpring GX v11.5.1 software package (Agilent Technologies).

Gene ontology (GO) and pathway analysis.
To discover the function and associated pathways of differentially expressed mRNAs, GO and pathway analyses were performed. GO annotations of microarray genes were downloaded from NCBI (http://www.ncbi.nlm.nih.gov/), UniProt (http://www.uniprot.org/) and the Gene Ontology (http://www.geneontology.org/). The elim Fisher algorithm was used to perform a GO enrichment test and GO categories with P<0.05 were reported. Pathway annotations of microarray genes were download from KEGG (http://www.genome.jp/kegg/) and a Fisher exact test was performed in order to locate the significant enrichment pathway. The resulting P values were adjusted using the Benjamini Hochberg false discovery rate (BH FDR) algorithm. Pathway categories with a FDR<0.05 were reported.

Construction of the lncRNA-mRNA co-expression network.
The network construction procedures included the following: i) preprocessing of data: if one coding gene has different transcripts the median value is taken to represent the gene expression values, without special treatment of lncRNA expression values; ii) data were screened and the subset of data were removed according to the lists of the differential expression of lncRNA and mRNA obtained from the GO and pathway analyses; iii) the Pearson correlation coefficient was calculated and the R value was used to calculate the correlation coefficient between lncRNA and coding genes; and iv) Pearson's correlation coefficient was used for screening; RNAs with a Pearson's correlation coefficient of ≥0.99 were considered significant and the lncRNA-mRNA co-expression network was constructed by Cytoscape software (The Cytoscape Consortium, San Diego, CA, USA). qPCR confirmation. To validate the microarray results, qPCR was performed. Total RNA was extracted from frozen tumor specimens using TRIzol reagent and then reverse transcribed      Statistical analysis. Statistical analysis was performed using Student's t-test for the comparison of the two groups in microarray, and analysis of variance for multiple comparisons. P<0.05 was considered to indicate a statistically significant result. The statistical significance of the microarray result was analyzed by fold change and Student's t-test. The FDR was calculated to correct the P-value. The threshold value used to screen differentially expressed mRNAs was a fold change of ≥2.0, and a fold change of ≥8.0 for differentially expressed lncRNAs.

Results
Overview of lncRNA profiles. The expression profiles of lncRNAs in paired samples were shown by calculating the log 2 fold-change tumor/normal (T/N). The agreement was formulated as follows: Fold change (FC) cut-off: 2.0. Positive fold change values indicated upregulation and negative values indicated downregulation. Log fold change means log 2 of absolute fold-change (log 2 FC). The fold change and P-value were calculated from the normalized expression. One pair of samples was excluded from analysis according to 3D-PCA.
A total of 110 lncRNAs were significantly differentially expressed between the tumor and control groups (≥8-fold). Twenty-two lncRNAs were upregulated and 88 lncRNAs were downregulated in the tumor group compared with the controls (Table II). Log 2 FC of upregulated lncRNAs in the tumor group ranged between 3.017291 to 4.581319, and -3.00191 to -6.10723 of downregulated lncRNAs. RP11-58A12.3 (log 2 FC=-6.10723) was the most significantly downregulated lncRNA and RNU12 (log 2 FC=4.581319) was the most significantly upregulated lncRNA. We observed that downregulated lncRNAs were more common than upregulated lncRNAs.
Overview of mRNA profiles. Up to 17,069 coding transcripts were detected in four pairs of samples using 30,215 coding transcripts probes. Among the four pairs of samples, 1,269 mRNAs were upregulated in tumor samples compared with the matched normal tissues, while 851 mRNAs were . qPCR analysis confirmed microarray data. After normalization to β-actin RNA, data were presented as mean ± SD (n=51) and the median value for each lncRNA was used for statistical analysis. The experiment was conducted in triplicate. # P<0.05 versus control. lncRNA, long noncoding RNA. downregulated. Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway analysis showed that the differentially expressed mRNAs may be associated with p53, bladder cancer, cell cycle and propanoate metabolism pathways (Fig. 1). These results suggest that bladder cancer has a variety of genetic and phenotypic characteristics.
Confirmation of microarray results by qPCR. Four lncRNAs were selected for further confirmation of microarray results using qPCR. These lncRNAs were among the most downregulated or upregulated lncRNAs. Data analysis showed that KRT19P3 was upregulated and TNXA, CTA-134P22.2 and CTC-276P9.1 were downregulated in bladder cancer samples compared with matched normal tissues (Fig. 2). These data support a strong correlation between the qPCR result and the microarray data (Fig. 3).

Construction of the lncRNA-mRNA co-expression network.
We constructed a coding-non-coding gene co-expression network based on the correlation analysis between the differentially expressed lncRNA and mRNA. We used lncRNAs and mRNAs with Pearson's correlation coefficients of no less than 0.99 to construct the network. In total, 79 lncRNAs and 103 mRNAs were included in the co-expression network. The co-expression network indicated that one mRNA or lncRNA may correlate with one to tens of lncRNAs. The co-expression network may suggest that the inter-regulation of lncRNAs and mRNAs is involved in bladder cancer.

Discussion
The human transcriptome is more complex than a collection of protein-coding genes and their splice variants (11)(12)(13). With the advent of whole-genome and transcriptome sequencing technologies, it was shown that at least 90% of the genome is actively transcribed (13). Although lncRNA was initially suggested to be transcriptional noise, recent evidence suggests that this part of the genome may play a major biological role in cellular development and human diseases (14,15). The newly discovered lncRNAs demonstrate developmental and tissue-specific expression patterns, and aberrant regulation in a variety of diseases, including cancer (16)(17)(18). However, the function of lncRNAs in tumor pathogenesis and growth is less well characterized. Recent studies have started to reveal the importance of lncRNAs in tumorigenesis in bladder cancer. Yang et al demonstrated that UCA1 regulates cell cycle progression through CREB via PI3K-AKT-dependent signaling pathways and may serve as a new diagnostic and therapeutic target in bladder cancer (8). Ying et al demonstrated that MALAT-1 expression levels were upregulated in bladder cancer that subsequently metastasized, and that increased expression of MALAT-1 activated the Wnt pathway to promote epithelial-mesenchymal transition and human bladder cancer cell metastasis (9). However, the genome-wide expression and functional significance of lncRNAs in bladder cancer remains unclear.
In this study, we evaluated the lncRNA expression profile in the tissue of four bladder cancer patients to reveal the potential role of lncRNAs in the pathogenesis of bladder cancer. Microarray techniques revealed a set of significantly differentially expressed lncRNAs, with 22 upregulated and 88 downregulated lncRNAs in bladder cancer tissue compared with matched normal tissue (fold change ≥8). For further validation of microarray results, qPCR was performed to evaluate the expression patterns of TNXA, CTA-134P22.2, CTC-276P9.1 and KRT19P3 in a total of 51 patients with bladder cancer. The qPCR results matched well with the data from the microarray.
A major function of lncRNAs is to modulate the epigenetic status of proximal or distal protein-coding genes via cis-and trans-acting mechanisms (19,20). We also observed that numerous lncRNA expression levels were significantly correlated with the expression of tens of protein coding genes by construction of the lncRNA-mRNA co-expression network. For example, RUNX1T1 and SLC25A4 were correlated with 21 and 28 mRNAs respectively.
In order to gain insight into the function of targets of lncRNAs, GO term and KEGG pathway annotation were applied to the target gene pool. KEGG annotation showed a significant association with the p53, bladder cancer, cell cycle and propanoate metabolism pathway gene expression in the bladder cancer group compared with the normal tissue group, indicating that deregulated lncRNAs may act by regulating protein-coding genes in these pathways.
We demonstrated for the first time the expression profiles of human lncRNAs in bladder cancer by microarray. We identified a collection of aberrantly expressed lncRNAs in bladder cancer compared to matched normal tissue. It is likely that these deregulated lncRNAs play a key or partial role in the development and/or progression of bladder cancer. Further study is required to determine whether these lncRNAs may serve as new therapeutic targets and diagnostic biomarkers in bladder cancer.