Differentially Expressed miRNAs in Tumor, Adjacent, and Normal Tissues of Lung Adenocarcinoma

Lung cancer is the leading cause of cancer deaths. Non-small-cell lung cancer (NSCLC) is the major type of lung cancer. The aim of this study was to characterize the expression profiles of miRNAs in adenocarcinoma (AC), one major subtype of NSCLC. In this study, the miRNAs were detected in normal, adjacent, and tumor tissues by next-generation sequencing. Then the expression levels of differential miRNAs were quantified by quantitative reverse transcription-polymerase chain reaction (qRT-PCR). In the results, 259, 401, and 389 miRNAs were detected in tumor, adjacent, and normal tissues of pooled AC samples, respectively. In addition, for the first time we have found that miR-21-5p and miR-196a-5p were gradually upregulated from normal to adjacent to tumor tissues; miR-218-5p was gradually downregulated with 2-fold or greater change in AC tissues. These 3 miRNAs were validated by qRT-PCR. Lastly, we predicted target genes of these 3 miRNAs and enriched the potential functions and regulatory pathways. The aberrant miR-21-5p, miR-196a-5p, and miR-218-5p may become biomarkers for diagnosis and prognosis of lung adenocarcinoma. This research may be useful for lung adenocarcinoma diagnosis and the study of pathology in lung cancer.


Introduction
Lung cancer is the leading cause of cancer-related deaths in the world [1]. Non-small-cell lung cancer (NSCLC) is the main type of lung cancer, which accounts for 80-85%. NSCLC mainly consists of squamous cell carcinoma (SCC) and adenocarcinoma (AC) [2]. Patients with NSCLC were frequently diagnosed at advanced stages, resulting in an overall 5-year survival rate of about 14% [3]. In order to improve the outcome of these patients, the early diagnosis of NSCLC becomes necessary. In recent years, studies showed that molecular biomarkers were discovered in various cancers and could be applied in the diagnosis and prognosis of NSCLC [4][5][6].
MicroRNAs (miRNAs) are endogenous noncoding small RNAs of about 18-25 nucleotides [7]. Mature miRNAs are highly conserved RNA molecules that can regulate the expression of genes by hybridizing to complementary sequences in the 3 -untranslated (UTR) region of target mRNAs [8]. Studies have shown that miRNAs participate in a lot of biological processes, including cell proliferation, cell differentiation, and cell apoptosis [9]. In addition, a large number of studies have shown that miRNAs may play roles in the tumorigenesis and development of lung cancer as tumor oncogenes or suppressor genes [10][11][12][13][14]. However, few researchers studied the expression profiles of miRNAs in paired tumor, adjacent, and normal tissues of lung cancer.
Three principal methods can be adopted to detect the expression levels of miRNAs: quantitative reverse transcription-polymerase chain reaction (qRT-PCR), microarray, and next-generation sequencing (NGS). Quantitative RT-PCR and microarray can only detect relative expression levels of known miRNAs, and the number of analyzed miRNAs is limited. Compared to qRT-PCR and microarray, NGS has an advantage in global miRNA analysis and it provides the possibility of detecting novel miRNAs. What is more, the analysis of miRNAs becomes faster and cheaper based on NGS [15].
In this study, the miRNAs were detected in tumor, adjacent, and normal tissues from AC pooled samples by NGS. The differential miR-21-5p, miR-196a-5p, and miR-218-5p were validated by qRT-PCR. The purpose of this study was to systematically characterize the expression of miRNAs in AC tumor, adjacent, and normal tissues. miR-21-5p, miR-196a-5p, and miR-218-5p can be biomarkers for diagnosis and prognosis of AC. The gradually changed miRNAs may be useful for AC diagnosis and the study of pathology in lung cancer.

Material and Method
2.1. Sample Preparation and RNA Isolation. We collected 18 surgically resected samples including paired tumor, adjacent, and normal lung tissues from 18 AC patients. Tissue samples were all obtained from Jiangsu Province Hospital with informed consent. Tissues were stored at −80 ∘ C and then ground for subsequent RNA isolation.
Total RNA was isolated from 50 mg tissue of each sample by mirVana6 miRNA Isolation Kit (Ambion, USA). The quantity and quality of obtained RNA were measured by Qubit5 2.0 Fluorometer (Life Technologies, USA). The total RNA isolated from 9 AC tumor tissue samples was pooled with equimolar amounts, respectively, and the same operation was applied to 9 AC adjacent and 9 normal tissues. The miRNAs of the 3 pooled samples were analyzed by nextgeneration sequencing. Another 9 paired AC samples were used for following qRT-PCR validation.

Small RNA Library Preparation for Ion Torrent Sequenc-
ing. The small RNA libraries for Ion Torrent sequencing were prepared according to the manufacturer's protocols of Ion Total RNA-Seq Kit v2 (Life Technologies, USA) and Ion Xpress6 RNA-Seq Barcode 1-16 Kit (Life Technologies, USA). The libraries were sequenced by Ion Torrent Proton6 sequencing platform after quality control at Nanjing Percare Biotechnology Co., Ltd.

Bioinformation Analysis.
The raw expression values (read counts) of Ion Torrent were filtered by quality control. The effective reads were analyzed through four steps. Firstly, the basic results were obtained, including the length distribution, mapping to reference genome, mapping to rfam, and the repeat distribution. Secondly, miRNAs (including novel miRNAs) were predicted. Thirdly, miRNA target genes were predicted and annotated. Lastly, other ncRNAs were annotated, including tRNA, rRNA, snoRNA, snRNA, and piRNA. Since low copy number was less reliable, the miRNAs with more than 10 copies were reserved for subsequent analysis. In order to compare miRNA expressions across databases, the total copy number of each sample was normalized to 1,000,000 [16,17]. According to some reports on multiple isomiRs from a given miRNA locus [18,19], isomiRs were also comprehensively surveyed. Sequences that matched the pre-miRNAs in the mature miRNA region ±4 nt (no more than 1 mismatch) were defined as isomiRs.

Quantitative Validation of miRNAs by qRT-PCR.
We performed quantitative RT-PCR (qRT-PCR) to validate the expression levels of miRNAs analyzed by NGS. The reverse transcription reaction was carried out with AMV reverse transcriptase (Takara, Japan) at 42 ∘ C for 60 min and 75 ∘ C for 15 min in a total volume of 20 L. The universal RT primer was 5 -GTCGTATCCAGTGCAGGGTCCGAGGTATTC-GCACTGGATACGACNNNNNN-3 . The quantitative PCR amplification was performed with SYBR Premix (Takara, Japan). The forward primers were miRNAs-specific, and the universal reverse primer was GTGCAGGGT CCGAGGT.U6 snRNA was used as internal reference. The reaction was incubated in ABI 7500 PCR system at 95 ∘ C for 5 min, followed by 40 cycles of 95 ∘ C for 15 s and 55 ∘ C for 45 s in a total volume of 10 L. The cycle threshold (Ct) values were calculated by automatic Ct setting of 7500 system SDS software v1.4 (Applied Biosystems, USA).

Prediction and Functional Analysis of Target mRNA.
We used the starbase v2.0 to predict the target gene of miRNAs verified by qRT-PCR [20,21]. The target genes were obtained by the intersection of three prediction software programs, which were TargetScan, PicTar, and miRanda based on the starbase. Then we used the Gene Ontology (GO) database for function enrichment [22,23]. In addition, DNA intelligent analysis (DIANA) software was used for pathway identification based on Kyoto Encyclopedia of Genes and Genomes (KEGG) database [24].
2.6. Statistical Analysis. We adopted the log 2 method to transform the normalized reads, and the expression patterns of the sequencing data were visualized by hierarchical clustering. The value of differential expression of miRNAs was calculated based on Poisson's distribution [25] and the threshold of value was determined by false discovery rate (FDR) [26]. The average Ct for each triplicate from qRT-PCR was used for analysis. Fold change in miRNA expression was calculated by ΔΔCt, normalized with ΔCt = AvgCt miRNA − AvgCt U6 [27]. A comparison between two groups was performed by using the -test. values of <0.05 were considered to indicate statistically significant differences.

Results of miRNA Sequencing.
Based on Ion Torrent sequencing platform, the total RNAs of the 3 pooled samples from 9 paired AC tissues were characterized. In each pooled sample, more than 50% reads were annotated to miRNAs and ncRNAs database (data not shown). In this study, we detected 259, 401, and 389 known miRNAs in AC tumor, adjacent, and normal tissues, with 65, 95, and 34 novel mature miRNAs detected in the 3 pooled samples, respectively (Table 1). In addition, the sequences mapped with reference genome, on Figure 1: Predicted stem-loop structures of precursors. These 3 precursors had the ability to be processed to form the 3 novel miRNAs in Table 2. The red curve indicated the location that could be processed to miRNAs. (a), (b), and (c) represent tumor, adjacent, and normal control tissues, respectively. average 89.72% reads, were mapped to the reference human genome (data not shown).
To confirm the existence of miRNA among our sequencing miRNAs, we performed a prediction of the miRNA precursor hairpin structure. The novel miRNAs which had the highest expression levels in each pooled sample were shown in Table 2. The secondary structure of their precursor ( Figure 1) had the characteristic of stem-loop hairpin structure, and mature miRNA was far from loops and bulges.
The biological feature of the reads annotated to miRNAs was analyzed afterwards. The length distributions of all miR-NAs in the 3 pooled samples focused on 21, 22, and 23 nt. As expected, reads with a length of 22 nt were the most abundant (data not shown). The peak distributions were consistent with the biological feature of mature miRNAs, indicating that mature miRNAs were well enriched during the preparation of miRNA sequencing library. Simultaneously, all the isomiRs were obtained in each pooled sample. Half of the most dominant isomiR sequences had lengths that were different from those of canonical miRNA sequences. Generally, the most dominant isomiR sequence may be longer or shorter than registered miRNA sequence through altering 5 or 3 ends, especially for the 3 ends (data not shown). Figure 2 illustrated all the isoform profiles of the differential miR-21 which were validated in subsequent experiment (see Sections 3.3 and 3.4).

The Differentially Expressed miRNAs in Tumor, Adjacent,
and Normal Tissues. The differential miRNAs were chosen after three procedures. Firstly, we removed the unreliable miRNAs which had low copy number (≤10 reads) in paired samples at the same time. Then the value and FDR value of miRNA reads were calculated, and the miRNAs with FDR < 0.01 were chosen. At last we obtained the applicable miRNAs with 2-fold or more change. So 86, 50, and 76 miRNAs were differentially expressed with 2-fold or greater change in tumor versus normal, adjacent versus normal, tumor versus adjacent tissues, respectively. As shown in Figure 3, with stricter cutoff criteria (fold change > 3), 10 (6 upregulated, 4 downregulated), 12 (8 upregulated, 4 downregulated), and 6 (1 upregulated, 5 downregulated) miRNAs were differentially expressed in tumor versus normal, adjacent versus normal, tumor versus adjacent tissues.   More 5p isomiRs and 3 isomiRs were found in miR-21 compared to 3p isomiRs and 5 isomiRs, respectively. The -axis meant the value of counts as log 10 .

The Gradually Changed miRNAs in Tumor, Adjacent, and
Normal Tissues. In this study, we attempted to find out the miRNAs which gradually changed from normal to adjacent to tumor tissues with 2-fold or greater change and FDR < 0.01 (the normal as control). As shown in Figure 4(a), there were 61 differential miRNAs in tumor but not in adjacent tissues compared to normal control, 25 miRNAs were expressed differentially in adjacent but not in tumor tissues, 51 miRNAs were differentially expressed in both tumor and adjacent tissues but not in control. In addition, 25 miRNAs were differentially expressed in tumor and adjacent tissues at the same time compared to normal control (Figure 4(b)). Among them, miR-21-5p and miR-196a-5p were gradually upregulated from normal to adjacent to tumor tissues and miR-218-5p was gradually downregulated from normal to adjacent to tumor tissues significantly (Figure 4(c)).

qRT-PCR Validation of Differential Expression of Gradually Changed miRNAs.
In order to validate the differential expression of these 3 gradually changed miRNAs from NGS results, we performed qRT-PCR to quantify the expression levels of miR-21-5p, miR-196a-5p, and miR-218-5p in 18 paired AC samples (including 9 paired AC samples used in miRNA sequencing). The results suggested that miR-21-5p, miR-196a-5p, and miR-218-5p had consistent gradual changes in normal, adjacent, and tumor samples with sequencing results ( < 0.001, Figure 5), and the magnitude of changes differed between qPCR and NGS. The expression levels of these 3 miRNAs quantified by qPCR were shown in Table 3.

Discussion
Non-small-cell lung cancer (NSCLC) accounts for 80-85% of all kinds of lung cancer, the leading cause of cancer deaths worldwide. Adenocarcinoma is one of the major subtypes of NSCLC. miRNAs have been shown to play roles in the tumorigenesis and the development of lung cancer. However, few studies characterize the profiles of miRNAs in paired tumor, adjacent, and normal tissue samples of adenocarcinoma at the same time.
In this study, we found 25 differential miRNAs coexpressed in tumor, adjacent, and normal tissue samples from AC patients based on next-generation sequencing. Among them, 3 miRNAs displayed gradual changes from normal to adjacent to tumor tissues. Then qRT-PCR was performed to confirm that miR-21-5p and miR-196a-5p were gradually upregulated and miR-218-5p was gradually downregulated from normal to adjacent to tumor tissues significantly. The predicted target genes of these 3 miRNAs may regulate in several tumor-related KEGG pathways, including NSCLC, colorectal cancer, pancreatic cancer, and chronic myeloid leukemia. Therefore, pathways in cancer ranked first in all of the enriched KEGG pathways. In particular, miR-21-5p might   Figure 5: The qRT-PCR results of 3 gradually changed miRNAs. miR-21-5p and miR-196a-5p were gradually upregulated from normal to adjacent to tumor tissues, and miR-218-5p was gradually downregulated from normal to adjacent to tumor tissues significantly ( < 0.01, ΔCt = AvgCt miRNA − AvgCt U6 ).
adenocarcinoma patients at the same time. The results showed and validated that miR-21-5p and miR-196a-5p were gradually upregulated, and miR-218-5p was gradually downregulated from normal to adjacent to tumor tissues. These indicated that miR-21-5p and miR-196a-5p might serve as oncogene, and miR-218-5p might act as tumor suppressor gene in lung cancer. However, we need to increase the depth of sequencing to validate the results and gain more data of miRNAs in the future study. miR-21-5p was comprehensively studied in lung cancer by many researchers. Seike et al. found that aberrantly increased expression of miR-21 was enhanced by activated EGFR signaling pathway and played a significant role in lung carcinogenesis in never-smokers [28]. Markou et al. revealed that mature miRNA-21 was significantly upregulated in NSCLC patients and was an independent negative prognostic factor for overall survival in NSCLC [29]. Xu et al. showed that miR-21 regulated cellular proliferation, invasion, migration, and apoptosis by targeting PTEN, RECK, and Bcl-2 in lung squamous carcinoma [30]. However, miR-218-5p and miR-196a-5p were little studied in lung cancer. miR-218 was deleted and downregulated in lung squamous cell carcinoma and may be a strong candidate tumor suppressing miRNA potentially involved in lung cancer [31]. Song et al. found that miR-218 suppressed the growth of lung carcinoma by reducing MEF2D expression [32]. Liu et al. showed that miR-196a promoted NSCLC cell proliferation and invasion through targeting HOXA5 [33]. Although the differential expression levels of the 3 miRNAs were approximate to our results in the previous studies, we validated their gradual changes in paired tumor, adjacent, and normal tissues of lung adenocarcinoma for the first time. The gradual changes may verify the function of miR-21-5p, miR-196a-5p, and miR-218-5p as oncogene or suppressor gene more accurately. As little existing research has studied the profiles and functions of miR-196a-5p and miR-218-5p in lung cancer, we will further validate their regulatory roles in lung cancer in future research.
In conclusion, we characterized the differential expression of miRNAs in paired tumor, adjacent, and normal tissue samples of lung adenocarcinoma based on NGS, and we found and validated the gradual changes of miR-21-5p, miR-196a-5p, and miR-218-5p from normal to adjacent to tumor tissues. These 3 miRNAs may regulate as oncogene or suppressor gene in tumorigenesis and progression of lung adenocarcinoma. Our research provided several miRNA biomarkers for lung adenocarcinoma and may be useful for lung adenocarcinoma diagnosis and the study of pathology in lung cancer.  Figure 7: The heatmap of significant KEGG pathways that 3 miRNAs may be involved in. The column showed each pathway and candidate target genes of miRNAs. The color showed red to yellow as log ( value) from low to high.