Comparison of the miRNA expression profiles in fresh frozen and formalin-fixed paraffin-embedded tonsillar tumors

MicroRNAs are considered as promising prognostic and diagnostic biomarkers of human cancer since their profiles differ between tumor types. Most of the tumor profiling studies were performed on rarely available fresh frozen (FF) samples. Alternatively, archived formalin-fixed paraffin-embedded (FFPE) tissue samples are also well applicable to larger-scale retrospective miRNA profiling studies. The aim of this study was to perform systematic comparison of the miRNA expression profiles between FF and macrodissected FFPE tonsillar tumors using the TaqMan Low Density Array system, with the data processed by different software programs and two types of normalization methods. We observed a marked correlation between the miRNA expression profiles of paired FF and FFPE samples; however, only 27-38% of the differentially deregulated miRNAs overlapped between the two source systems. The comparison of the results with regard to the distinct modes of data normalization revealed an overlap in 58–67% of differentially expressed miRNAs, with no influence of the choice of software platform. Our study highlights the fact that for an accurate comparison of the miRNA expression profiles from published studies, it is important to use the same type of clinical material and to test and select the best-performing normalization method for data analysis.


Introduction
MicroRNAs (miRNAs) are small non coding RNAs (~21 nucleotides) which play an important role in post-trancriptional regulation of gene expression. Their binding with perfect or imperfect complementarity to the 3´untranslated region of the target mRNA leads either to mRNA cleavage and degradation or to inhibition of mRNA translation [1]. Of interest, Ørom et al. also reported an association of miRNA with the 5´untranslated region which resulted in the activation of gene expression [2,3]. MiRNAs are involved in many biological processes such as PLOS  Hospital and Institute of Hematology and Blood Transfusion. The study set was selected based on the presence of HPV DNA and HPV E6 mRNA we defined in our previous studies [11,15]. The sampling and tissue handling for FF samples has been described before [11,15]. Sections of FF tumor samples were cut on a cryostat, and the number of tumor cells was determined by a pathologist. All FF samples contained more than 50% of tumor cells. FFPE samples were macrodissected as described before [11].

Processing of samples
Total RNA from all FF tissues was isolated by the miRVana kit (Life Technologies, USA) according to the manufacturer's protocol. Total RNA from FFPE samples was isolated simultaneously with DNA from four 10-μm sections enriched for the tumor cells by macrodissection using the Ambion RecoverAll™ Total Nucleic Acid Isolation kit for FFPE according to the manufacturer's protocol (Applied Biosystems, USA). RNA concentration and quality was measured by a Nanodrop™ Spectrophotometer (Thermo Scientific, USA) and Experion chip electrophoresis (Bio-Rad, USA).

TaqMan Low Density Array (TLDA) analysis
The miRNA expression profiling was performed by the TaqMan 1 Array Human MicroRNA A +B Cards Set v3.0 (Life Technologies, USA) containing a total of 384 miRNA probes and controls per card. Overall, we analyzed 15 FF samples and 15 FFPE samples (five HPV-positive tumor samples, five HPV-negative tumor samples, and five non-malignant tissues of each material). For the analysis of FF samples, 1000 ng of total RNA was used for reverse transcription. In the case of FFPE samples, the input was 350 ng due to the smaller concentration of extracted RNA from FFPE samples. Total RNA of each sample was reverse transcribed using the Taq-Man 1 MicroRNA Reverse Transcription Kit and Megaplex™ RT Primers specific for each card (both Life Technologies, USA) according to the manufacturer's protocol. The workflow without preamplification was chosen for both types of material, FF as well as FFPE samples. The TLDA cards were analyzed on the Applied Biosystems 7900HT Real-Time PCR System.

Data analysis
Each set of samples was evaluated separately, and the results were compared. For data processing and evaluation, we used the SDS 2.4 and the ExpressionSuite v1.0.1 software (Life Technologies, USA). Ct values identified by automated thresholding were exported separately for cards A and B, and the RQ (relative quantity) was calculated from the detected Ct values using the 2 -Ct formula. Data were further processed by the GeneSpring GX v13.1 software (Agilent), GenEx v6.1 (MultiD), and MeV v4.9 (TIGR TM4 suite) for comparison. We also applied two types of data normalization-50 th percentile shift and global normalization-to perform within-sample normalization, globally across the data sets. For comparison of two groups in censored analyses, the T-test in combination with Pavlidis Template Matching (PTM) algorithm were used. In the comparison of differential abundance of miRNAs in the groups, the P-value (P<0.05) and fold-change (FC>1.33) were set. Only the miRNAs differentially expressed in at least 3/5 or 6/10 (60%) samples with measured results were considered for further analyses. The correlation was evaluated in the GraphPad InStat 3.0 tool using the Spearman nonparametric correlation. by the individual RT qPCR using the TaqMan 1 MicroRNA Assays (Life Technologies, USA) as described recently [11]. Total RNA from five FFPE samples of non-malignant tonsils were used as a calibrator, and RNU48 served as the endogenous control. The data were analyzed in the GenEx v6.1 (MultiD) software. The 2 -ΔΔCt method was used for calculations of the fold change. The cut-off fold change was set as for the arrays to +/-1.33 (33% fold change). T-test or Mann Whitney nonparametric test was applied depending on the data distribution. All results were statistically significant (P-value 0.05).

Results
The quality and concentration of total RNA was analyzed by a Nanodrop™ Spectrophotometer and Experion chip electrophoresis. The RNA integrity number (RIN) was higher than seven in all FF samples. As documented in Fig 1, total RNA isolated from FFPE samples showed characteristic profiles for degraded RNA fragments (RIN in the range 1.9-2.6), but fragments with a length of around 200 bp were present allowing for a miRNA expression analyses. As illustrated by the principal component analysis (PCA) (Fig 2), the miRNA expression profiles of all samples differ depending on the analyzed type of material (FF vs FFPE) and characteristics of the tissue, tumor vs normal tissue (TT vs NT). The correlation between the miRNA expression profiles of paired FF and FFPE samples was evaluated using the Spearman correlation test. The correlation of the expression level for the well characterized miRNAs located on card A was better (correlation coefficient !0.65 in 7/10 samples (p<0.0001)) in comparison to all analyzed miRNAs (miRNAs on cards A and B). The exact Spearman correlation values of each comparison and the confidence intervals are summarized in Table 1. The correlation of paired samples is also visualized in scatter plots (Fig 3). Our assumption that FF samples with higher percentage of tumor cells would be more likely to give better correlation of the expression miRNA profiles with FFPE macrodissected samples was not proven. The correlation coefficient of the miRNA expression profiles of FF and FFPE samples was not influenced by the number of tumor cells in the particular FF sample (Table 1).
Differentially expressed miRNAs were identified by the P-value (P<0.05) and fold-change >1.33 with the 50 th percentile shift normalization method using the GeneSpring software. In our study, macrodissected FFPE samples revealed less differentially expressed miRNAs (N = 58) than paired fresh frozen tumors (N = 83 miRNAs). When comparing the lists of deregulated miRNAs between FF and FFPE samples, we detected overlaps in 27-38% (Table 2, on the left).
To establish how the method used for normalization influences the final results, we applied, besides the 50 th percentile shift normalization, also global normalization for each group of the analyzed samples. As illustrated in Fig 4 and Table 3, the differentially expressed miRNAs based on each normalization method overlap in 58-67% for FF and/or FFPE samples. Using different software and the identical normalization method, almost 90% of the deregulated miRNAs were concordant between FF and FFPE samples (data not shown). We performed additional analysis and identified 16 miRNAs commonly detected in both type of material and using both type of normalization (Table 4). These miRNAs were screened in database miRSearch V3.0 (Exiqon). Most of them have experimentally proved target genes participating in signaling pathways and cellular processes including not only tumor development. Noteworthy, target genes are transcription factors E2F, MAP kinases, members of RAS protein family, tumor suppressor PTEN, transforming growth factors TGFs, tumor suppressor TP53, proteins participating in cell cycle (cyclins, cyclin-dependent kinases) or regulating cell death.
Finally, from the deregulated miRNAs in FFPE samples, we have selected 11 tumor-specific miRNAs, also with regard to viral or non-viral etiology of the tumor. The expression of these miRNAs was evaluated in a larger set of 64 macrodissected tumor samples by individual Taq-Man assays (Fig 5). MiR-106b# and miR-9 were selected as specific for HPV-positive tonsillar tumors, miR-16, miR-34a, miR-193b, miR-31, miR-221, and miR-21 as specific for HPV-negative tumors, and miR-155, miR-126, and miR-205 as specific for tonsillar tumors of any etiology. The differential expression of six miRNAs (miR-106b#, miR-9, miR-16, miR-34a, miR-155, and miR-126) (P<0.05; FC>1.33) was confirmed in a larger set of 64 tumor samples, the fold change of miR-193b was equal to 1.23, but the trend of expression was also maintained. The trend of expression of four miRNAs (miR-31, miR-221, miR-21, and miR-205) analyzed in the large set of samples was opposite to that derived from the results revealed by arrays.

Discussion
The archives of formalin-fixed paraffin-embedded material in pathology laboratories offer the possibility to analyze the samples retrospectively and provide the option for clinical cancer research to extend the studies, and, moreover, to correlate the results with clinical data obtained during long-term follow-up of patients. The utilization of FFPE material for RNA based studies is still challenging, especially due to RNA fragmentation. However, the small size of miRNAs contributes to their stability during the processing of FFPE samples [9] so that FFPE tissues are applicable for miRNA expression analysis studies as previously shown [16][17][18][19]. Chatterjee et al. have presented a cross comparison of miRNA detection technologies in FFPE samples and revealed that it is possible to obtain high quality sequence reads for miRNA profiling in FFPE samples with an RNA integrity number around 2 [20]. In our study, we performed the miRNA expression profiling by TLDA analysis in a set of macrodissected FFPE samples and compared the results with the miRNA expression profiles of paired fresh frozen samples. Most studies comparing the miRNA expression between FF and FFPE samples were performed using only several selected miRNAs [21][22][23][24][25].
A comparison of the miRNA expression profiles in FF and FFPE samples obtained by the microarray analysis has also been published. Romero-Cordoba et al. have performed a comparative analysis of paired samples of breast cancer using the TaqMan Low Density Array platform, similarly to Hui et al. [12,14]. They revealed a high correlation; however, both studies analyzed small groups of paired samples and did not use dissected samples. Goswami et al. have also observed a high correlation between the samples; however, their workflow differed   in several aspects (assessment of extraction method, platform optimization, etc.) [13]. Additional studies were performed using other microarray platforms [9,26,27] or next generation sequencing platforms [28]. None of this comparative studies has been performed in head and neck tumors. Since our recent study has shown that the tumor tissue homogeneity is important for the robustness of the miRNA profiling studies, we performed macrodissection of FFPE samples to enrich them for tumor cells. To our knowledge, only one study dealt with the issue of FFPE dissection. It tested the possibility of a non-satisfactory correlation of the expression of several miRNAs between FF and dissected FFPE samples due to the enrichment of the tumor cells.   However, the improvement in the correlation when using non-dissected FFPE samples was due to a low number of samples which was not convincing [25]. In our analysis, we compared FF tumor samples with more than 50% of tumor cells and macrodissected FFPE tumor cells.
Although, based on the PCA, the miRNA expression profiles differ between FF and FFPE, we performed the Spearman correlation test which showed a good correlation (Spearman correlation coefficient in the range from 0.65 to 0.84) in 70% of samples. Based on the study of de Biase et al. [25] who have pointed to the fact that dissection can influence the expression results which is necessary to take into account if the miRNA analysis is performed with or without dissection, we evaluated the correlation of the differentially expressed miRNAs in FF and FFPE based on the number of tumor cells in FF samples. However, no relationship between the number of tumor cells and correlation level was observed.
We revealed higher numbers of differentially regulated miRNAs in FF tonsillar tissues compared to the FFPE material which we attributed to the higher amount of RNA used for the miRNA expression analysis in FF samples. Therefore, some less expressed miRNAs could have been missed in FFPE samples. The overlap of differentially expressed miRNAs between FF and FFPE samples was around 30% which should be explained by the fact that unlike the tumor samples, the control non-malignant tissues were not paired between FF and FFPE and could influence the miRNA expression profiles more than we expected. Nevertheless, the analysis of both groups of tested material revealed the deregulation of the known tumor associated miR-NAs such as miR-205, miR-210, or miRNAs from the family mir-8, as has been demonstrated in a number of solid tumors [29][30][31][32].
Further, when comparing the results obtained using two different types of data normalization, we revealed an overlap of between 58-67%. Therefore, the comparison between studies is obviously difficult because of the use of different methods of data analysis. In our analysis, both tissue types and both normalization methods overlap in 16 common miRNAs. All of these miRNAs were shown to be cancer related, nevertheless their role in head and neck cancer has not been elucidated. The role of miR-27a was mentioned by Venkatesh et al. as it regulates the expression of tumor suppressors in oral squamous cell carcinomas [33]. The role of miR-126# was found in endothelial proliferation [34,35]. For miR-140-3p it was shown that it contributes to reduction of proliferation and migration of cancer cells in breast and lung cancer [36,37] and in spinal chordoma, in contradictory, to be a marker of poor prognosis [38]. MiR-141 plays role in AKT signaling and inhibits prometastatic mesenchymal characteristics [39]. MiR-143 is a part of miR cluster 143/145 expressed in many tissues and regarded as tumor suppressing [40]. It´s reduced expression was found also in cervical squamous cell carcinomas [41] and low expression of this miRNA contributes to poor prognosis in oral squamous cell carcinomas [42]. As recently discovered the expression of miR-196b promotes cell migration and invasion in oral cancer [43,44]. Members of mir-200 family, which includes miR-200b, miR-200c and miR-205, are the major regulators of EMT pathway, primarily targeting transcriptional factors ZEB1 and SP1 [45], and play a prognostic role in various malignancies [46]. Deregulated expression of hypoxia-induced miR-210 was found in various tumors and influences cancer cell proliferation, apoptosis, angiogenesis and other cellular processes involved in tumor development [47]. The group of Seki et al. have been interested in the functional significance of mir-29 family, including miR-29a, in cervical cancer as well as in head and neck cancer. They found out that this tumor-suppressive miRNA inhibits cancer cell migration and invasion [48,49]. Although some studies focused on the functional impact of identified miR-NAs in HNC have been published, more extensive research is necessary to more accurately explain the relevance in HNC tumor development.
From the groups of deregulated miRNAs identified in FFPE tumors of either viral or nonviral etiology, we selected 11 miRNAs based on the results and relevance in the literature, the expression of which was evaluated in a larger set of 64 FFPE samples. From the group of miR-NAs specific for HPV-positive tonsillar tumors, we choose upregulated miR-9 and miR-106b# whose deregulation was previously reported in HNC [50][51][52]. Furthermore, miR-9 has been shown to be activated by HPV leading to increased cell motility [53] and to be involved in the pathways regulating metastasis [54]. Upregulated miR-31, miR-221, and miR-21 were selected for confirmation from the group of miRNAs specific for HPV-independent tumors including HNC as playing a role in increasing cell proliferation, invasion, and migration [55][56][57][58], participating in regulation of epithelial-mesenchymal transition (EMT), or having a prognostic impact [59,60]. Further choices were miR-34a, specific for HPV-negative tumors, a tumor suppressor whose downregulation has been shown in a number of tumor types including HNC [61,62], as well as miR-16, the inhibition of which promotes cell proliferation, migration, invasion, and EMT and contributes to tumor progression [63]. The last miRNA selected for confirmation from the group of HPV-negative tumors was miR-193b whose high expression in tissue was identified as an independent prognostic risk factor in patients with ovarian cancer [64]. From the group of tonsillar specific miRNAs regardless of the virus presence, we selected miR-205, miR-155, and miR-126 which participate in oncogenic pathways including cell proliferation, migration, or invasion and serve as prognostic predictors of patients with HNC [31,65,66].
We were able to confirm deregulated expression in the same set of samples as run on array. In this subset, deregulated expression was confirmed for all miRNAs except miR-21 and miR-205. However, their trend of expression was maintained and the fold change was close to the cut-off value. In the larger set of samples, deregulated expression was confirmed in 64% of comparisons. The fold change of miR-193b was equal to 1.23, and the trend of expression of four additional selected miRNAs was opposite to that revealed by arrays. The results suggest possible variability between samples and show that not all deregulated miRNAs detected in a smaller tested cohort might be representative of a broader spectrum of samples.
In conclusion, our study compares the miRNA expression profiles in tonsillar tumors as identified in fresh frozen and formalin-fixed paraffin-embedded samples. Although the correlation between the defined groups was relatively good for all miRNAs, the overlap of the selected differentially expressed miRNAs was rather suboptimal. We concluded that for accurate comparison between studies, the key factors are the use of the same type of clinical material and selection of the normalization method for data analysis. Additionally, we suggest that combining multiple analytical methods rather than a single test alone is advisable for more robust, reliable detection of differential abundance of miRNAs and interpretation of results.
Genetics Services Platform of the International Agency for Research on Cancer for help with additional TLDA card processing.
The paper was supported by grant by the Ministry of Education, Youth and Sports of CR within the National Sustainability Program II (Project BIOCEV-FAR) LQ1604 and by the project "BIOCEV" (CZ.1.05/1.1.00/02.0109).