Identification of Tumor Suppressive Genes Regulated by miR-31-5p and miR-31-3p in Head and Neck Squamous Cell Carcinoma

We identified the microRNA (miRNA) expression signature of head and neck squamous cell carcinoma (HNSCC) tissues by RNA sequencing, in which 168 miRNAs were significantly upregulated, including both strands of the miR-31 duplex (miR-31-5p and miR-31-3p). The aims of this study were to identify networks of tumor suppressor genes regulated by miR-31-5p and miR-31-3p in HNSCC cells. Our functional assays showed that inhibition of miR-31-5p and miR-31-3p attenuated cancer cell malignant phenotypes (cell proliferation, migration, and invasion), suggesting that they had oncogenic potential in HNSCC cells. Our in silico analysis revealed 146 genes regulated by miR-31 in HNSCC cells. Among these targets, the low expression of seven genes (miR-31-5p targets: CACNB2 and IL34; miR-31-3p targets: CGNL1, CNTN3, GAS7, HOPX, and PBX1) was closely associated with poor prognosis in HNSCC. According to multivariate Cox regression analyses, the expression levels of five of those genes (CACNB2: p = 0.0189; IL34: p = 0.0425; CGNL1: p = 0.0014; CNTN3: p = 0.0304; and GAS7: p = 0.0412) were independent prognostic factors in patients with HNSCC. Our miRNA signature and miRNA-based approach will provide new insights into the molecular pathogenesis of HNSCC.


Introduction
Head and neck squamous cell carcinoma (HNSCC) arises from the oral cavity, larynx, or pharynx and is ranked the sixth most common cancer [1,2]. In 2018, approximately 84,000 cases of HNSCC were newly diagnosed, and more than 43,000 people died of this disease worldwide [3]. Surgery, radiation therapy, and cisplatin-based chemotherapy are the main treatment strategies for head and neck cancers. At the time of the initial diagnosis, most patients have advanced-stage disease and a poor prognosis (5-year survival rate < 60%) due to lymph node metastasis or recurrence [1]. In addition, cancer cells acquire resistance to cisplatin-based treatment, and the prognosis of patients who fail treatment is extremely poor [4]. The therapeutic effects of molecular-targeted drugs and immune checkpoint inhibitors in patients after treatment failure are poorly understood [5,6].
The Human Genome Project revealed that an extremely large number of non-coding RNAs (ncRNAs) are transcribed from the human genome, and these ncRNAs function in both normal and diseased cells [7,8]. Among ncRNAs, microRNAs (miRNAs) are endogenous single-stranded RNA molecules [19][20][21][22] nucleotides long that function as fine tuners of RNA expression in a sequence-dependent manner [9,10]. A unique feature of miRNAs is that a single miRNA negatively regulates a vast number of RNA transcripts (both protein-coding RNAs and ncRNAs) in each cell [10]. Moreover, bioinformatic studies have shown that more than half of protein-coding genes are controlled by miRNAs [11]. Numerous studies have indicated that aberrantly expressed miRNAs disrupt the tightly controlled RNA networks in normal cells, and these events trigger transformation to a disease state [12,13].
Identification of differentially expressed miRNAs in the cancer tissues of interest is the initial step. The latest RNA sequencing technology has successfully resulted in identification of miRNA expression signatures in cancer tissues. Our HNSCC miRNA signature revealed that both strands of the miR-31 duplex (miR-31-5p and miR-31-3p) are upregulated in cancer tissues. Numerous cohort data from The Cancer Genome Atlas (TCGA) confirmed that miR-31-5p and miR-31-3p are upregulated in HNSCC tissues. The aim of this study was to investigate the oncogenic roles of these miRNA strands and to identify their tumor suppressor gene targets in HNSCC cells.
Identification of differentially expressed miRNAs and their regulated molecular networks may be an effective strategy for elucidating the molecular pathogenesis of HNSCC.

Identification of the miRNA Expression Signature of HNSCC by RNA Sequencing
Six cDNA libraries (derived from three HNSCC tissues and three normal oral epithelial tissues) were analyzed by RNA sequencing. After a trimming procedure, 955,347-1,927,436 reads were successfully mapped to the human miRNAs (Table S1). The clinical features of the HNSCC specimens using in this study are summarized in Table S2.
A total of 168 miRNAs were identified as upregulated (log2 fold change > 1.5) in HNSCC tissues ( Figure 1A and Table S3).   were evaluated using TCGA-HNSC data. (C) Kaplan-Meier survival analyses of HNSCC patients using data from TCGA-HNSC. Patients were divided into two groups according to the median miRNA expression level: high and low expression groups. The red and blue lines represent the high and low expression groups, respectively.

Expression Levels and Clinical Significance of miR-31-5p and miR-31-3p in HNSCC
We focused on miRNAs of which both strands (the guide strand and passenger strand) derived from pre-miRNAs were upregulated in this signature. A total of 7 pre-miRNAs (miR-31, miR-223, miR-4655, miR-4781, miR-6753, miR-6830, and miR-6871) were detected in this signature ( Figure 1A and Table S3). From TCGA-HNSC database analysis, it was confirmed that miR-31 is the only miRNA whose expression of both strands were significantly upregulated in HNSCC tissues among 7 pre-miRNAs ( Figure 1B). The expression of neither miRNA was associated with worse overall survival rates in patients with HNSCC according to analysis of TCGA-HNSC data ( Figure 1C).
In this study, we focused on miR-31-5p and miR-31-3p, and continued to investigate the functional aspects of these miRNAs.

Effects of Inhibition of miR-31-5p and miR-31-3p Expression on the Proliferation, Migration, and Invasion of HNSCC Cells
First, we measured the expression levels of miR-31-5p and miR-31-3p in 11 HNSCC cell lines compared with fibroblast cell lines (IMR-90 and MRC-5). Detailed information on the cell lines used is shown in the Table S4. Overexpression of miR-31-5p and miR-31-3p was detected in several HNSCC cell lines, e.g., Ca9-22, HSC-2, HSC-4, and SAS ( Figure S1), relative to fibroblasts. We selected two of these HNSCC cell lines, SAS and HSC-2, to investigate the oncogenic roles of these miRNAs. To suppress the expression of miR-31-5p and miR-31-3p, we used inhibitors (Anti-miR TM miRNA Inhibitor) of these miRNAs. The inhibitors were used at a concentration of 30 nM. To evaluate their effects in functional assays, we confirmed the expression of miR-31-5p and miR-31-3p after transfection of inhibitors into SAS and HSC-2 cells ( Figure S2). Inhibition of miR-31-5p and miR-31-3p attenuated the proliferation (Figure 2A and Figure S3) and markedly decreased the migration and invasion ( Figure 2B,C and Figure S4) of SAS and HSC-2 cells. These data suggest that upregulation of miR-31-5p and miR-31-3p has an oncogenic effect in HNSCC cells.

Direct Regulation of CGNL1 by miR-31-3p in HNSCC Cells
We focused on CGNL1, which has the most significant difference in clinical statistics, from among the five target genes of miR-31-5p and miR-31-3p, and verified the direct regulation of CGNL1 by miR-31-3p. In cells transfected with miR-31-3p, the levels of CGNL1 mRNA and CGNL1 protein were significantly lower than in mock-or miR-controltransfected cells (Figure 8A,B).  [14]. (D) Dual-luciferase reporter assays showed decreased luminescence activity in SAS cells co-transfected with miR-31-3p together with a vector harboring the "wild-type". Normalized data were calculated as Renilla/firefly luciferase activity ratios.
We performed dual-luciferase reporter assays to determine whether CGNL1 was directly regulated by miR-31-3p. We used vectors encoding the partial wild-type sequences of the 3 -UTR of CGNL1 and vector with partially deleted CGNL1 3 -UTR ( Figure 8C). We found that luciferase activity was significantly decreased by cotransfection with miR-31-3p and the vector carrying the wild-type 3 -UTR of CGNL1, whereas transfection with the deletion vector blocked the decrease in luminescence in SAS cells ( Figure 8D). These data showed that miR-31-3p directly bound to CGNL1.

Discussion
The latest RNA-sequencing technologies have enabled identification of genome-wide miRNA expression signatures in human cancers. Our recent studies of miRNA signatures revealed that the passenger strands of some miRNA duplexes, such as miR-99a-3p, miR-145-3p, miR-150-3p and miR-199a-3p, act as tumor suppressors by directly targeting several oncogenes in HNSCC cells [15][16][17][18]. The original theory regarding miRNA biogenesis is that the guide strand of the miRNA duplex is incorporated into the RNA-induced silencing complex and functions as a negative regulator of gene expression, whereas the passenger strand is degraded in the cytoplasm and nonfunctional [10,11]. However, numerous in silico studies (involving over 5200 patients with 14 types of cancers) have shown that both strands (5p and 3p) of some miRNA duplexes (e.g., miR-30a, miR-139, miR-143, and miR-145) function together to regulate pivotal targets and pathways in several types of cancers [19]. Studies on the passenger strands of miRNAs will reveal novel molecular mechanisms of cancer pathogenesis.
In this study, we focused on miR-31-5p and miR-31-3p based on our miRNA signatures. Upregulation of miR-31-5p and miR-31-3p in HNSCC tissues was confirmed by TCGA data analysis. Our functional assays indicated that these miRNAs act as oncogenic miRNAs in HNSCC cells. Previous studies demonstrated that miR-31 has opposing roles (oncogene vs. tumor suppressor) depending on the type of cancer [20].
In esophageal squamous cell carcinoma, miR-31 was upregulated in clinical specimens and acted as an oncogenic miRNA by targeting the tumor suppressor gene LATS2, which is involved in the Hippo pathway [21]. Upregulation of miR-31 was reported in HNSCC tissues, and its expression activated the hypoxia-inducible factor pathway by targeting factor-inhibiting hypoxia-inducible factor [22,23]. Signaling via epidermal growth factor and its receptor is an essential oncogenic pathway in HNSCC and oral squamous cell carcinoma (OSCC), and this signaling pathway enhanced AKT activation and upregulated C/EBPβ expression in OSCC [24]. These events induced upregulation of miR-31 in OSCC cells [24]. Interestingly, a previous study showed that exogenous expression of miR-31 and telomerase reverse transcriptase transformed normal oral keratinocytes into immortalized cells [25]. Those previous and our present results indicate that upregulation of miR-31 downregulates genes/pathways intricately involved in malignant transformation of HNSCC and OSCC.
Our other aim was to clarify the novel molecular pathways regulated by miR-31-5p and miR-31-3p in HNSCC cells. Our in silico analysis revealed that five genes (CACNB2, IL34, CGNL1, CNTN3, and GAS7) are closely associated with HNSCC molecular pathogenesis. Functional analyses of these genes are needed to reveal the molecular mechanisms of HNSCC malignant phenotypes.
Of the five genes, GAS7 was initially cloned from serum-starved mouse NIH3T3 cells, and it consists of a series of different functional domains from the N-to C-termini: Src homology 3 domain, WW domain, and FES-CIP4 homology domain [26]. GAS7 regulates the dynamic activities of the membrane, actin cytoskeleton, and microtubules [26,27]. Downregulation of GAS7 has been reported in several cancer types, and ectopic expression of GAS7 inhibited the migration of lung and breast cancer cells [28]. More recently, it was reported that loss of GAS7 expression accelerated metastasis of neuroblastoma harboring MYCN overexpression or amplification [29]. Previous studies indicated that GAS7 acts as a tumor suppressor in human cancers.
Analysis of TCGA data showed that IL34 is downregulated in HNSCC tissues, and its low expression significantly predicts a poor prognosis in patients with HNSCC. IL34 stimulates the differentiation of monocytes into macrophages via the CSF-1 receptor [30]. IL34 is also a ligand of the macrophage colony stimulating factor receptor [31]. Recent studies showed that IL34 is expressed in various types of cancers and is involved in cancer progression and metastasis [32]. In the future, it is necessary to investigate the functional significance of IL34 in HNSCC.
CGNL1 is a paralogue of cingulin, which is ubiquitously expressed in endothelial cells and localized at tight junctions [33,34]. A previous study showed that cingulin binds to actin filament bundles to bridge tight junctions and actin filaments [35]. CGNL1 is localized on actin filament bundles and has multiple roles depending on its binding partner [35]. Previous reports showed that CGNL1 is an inhibitor of RhoA activity in tight junctions but is also involved in Rac1 activation in Madin-Darby canine kidney epithelial cells [36,37]. CGNL1 likely has various functions by interacting with different types of GEFs and GAPs in each cell. Few detailed functional analyses of CGNL1 have been performed in cancer cells. Expression of CGNL1 was downregulated in HNSCC tissues compared with normal epithelial tissues in numerous TCGA datasets. GEPIA2 database (http://gepia2.cancer-pku.cn/#index, accessed on 20 April 2021) analyses showed that expression of CGNL1 was significantly downregulated in cervical squamous cell carcinoma, esophageal carcinoma, and lung squamous cell carcinoma, suggesting that CGNL1 is downregulated in human squamous cell carcinoma. These findings suggest that CGNL1 plays a tumor suppressor role in HNSCC cells [38]. Sufficient functional analysis of CGNL1 in HNSCC remains unresolved in this study. By clarifying the tumor suppressive function of this gene in the future, a part of the molecular mechanism of HNSCC will be clarified.
We newly created the miRNA expression signature of HNSCC by RNA sequencing. Analysis of the signature revealed that both strands of pre-miR-31 (the guide strand of miR-31-5p and the passenger strand of miR-31-3p) acted as oncogenic miRNAs in HNSCC cells. Our in silico analysis showed that a total of 5 genes (CACNB2, IL34, CGNL1, CNTN3, and GAS7) were independent prognostic factors in patients with HNSCC. Our HNSCC miRNA signature and miRNA-based analyses will provide important insights into the molecular pathogenesis of HNSCC.

Clinical HNSCC and Normal Epithelial Tissue Specimens and HNSCC Cell Lines
Six specimens (three HNSCC tissues and three normal oral epithelial tissues) were analyzed by RNA sequencing to determine the HNSCC miRNA signature. The clinical features of HNSCC patients are summarized in Table S2.
All specimens used were obtained by surgical resection at Chiba University Hospital. All patients provided written informed consent for the use of their specimens. This study was approved by the Bioethics Committee of Chiba University (approval number: 28-65, 10 February 2015).
Two human HNSCC cell lines (HSC-2 and SAS) were obtained from the RIKEN BioResource Center (Tsukuba, Ibaraki, Japan) and used in this study.

Determination of the miRNA Expression Signature in HNSCC by RNA Sequencing
Small RNAs were sequenced to determine the miRNA expression signature of HNSCC. The RNA sequencing procedure was described in our previous studies [39][40][41][42].

Transfection of Mirnas Precursors and Inhibitors into HNSCC Cells
The procedures used for transfection of miRNA precursors and inhibitors into HNSCC cells have been described previously [15][16][17][18]. The reagents used in this study are listed in Table S5.

Functional Assays (Cell Proliferation, Migration, and Invasion) in HNSCC Cells
The procedures used for the functional assays in cancer cells (proliferation, migration, and invasion) have been described in our previous studies [15][16][17][18]. Cells were transfected with 30 nM miRNA inhibitors. Cell proliferation was evaluated by XTT assay. Migration assays were performed using uncoated transwell polycarbonate membrane filters, and invasion assays were conducted using modified Boyden chambers.

Analysis of the Clinical Significance of HNSCC Patients Using TCGA-HNSC Data
The strategy used to identify miRNA target genes is presented in Figure 3. We selected putative target genes with miR-31-5p and miR-31-3p binding sites using TargetScanHuman ver. 7.2 (http://www.targetscan.org/vert_72/; data were downloaded on 10 July 2020). The expression profiles of HNSCC clinical specimens (genes downregulated in HNSCC tissues) were used for screening miRNA target genes [14]. Our expression data were deposited in the GEO database (accession number: GSE172120). Furthermore, we narrowed down the candidate genes by factoring in clinical information from TCGA-HNSC analyses.
Multivariate Cox regression analyses were also performed using TCGA-HNSC clinical data and survival data according to the expression level of each gene from OncoLnc to identify factors associated with HNSCC patient survival [43]. In addition to gene expression, the tumor stage, pathological grade, and age were evaluated as potential independent prognostic factors. The multivariate analyses were performed using JMP Pro 15.0.0 (SAS Institute Inc., Cary, NC, USA).

Western Blotting
Cell lysates were prepared 48 h after transfection with RIPA buffer (Nacalai Tesque, Chukyo-ku, Kyoto, Japan). Then, 20 µg of protein lysates were separated on 4-12% Bis-Tris gel and transferred to nitrocellulose membranes (Invitrogen, Carlsbad, CA, USA) and blocked for 1 h at room temperature with Blocking One (Nacalai Tesque, Inc., Kyoto, Japan). The antibodies used in this study are shown in Table S5.

Plasmid Construction and Dual-Luciferase Reporter Assays
The partial wild-type sequence of the CGNL1 3 -untranslated region (3 -UTR) was inserted between the XhoI-PmeI restriction sites in the 3 -UTR of the hRluc gene in the psiCHECK-2 vector (C8021; Promega, Madison, WI, USA). Alternatively, we used sequences that were missing the miR-31-3p target sites. The synthesized DNA was cloned into the psiCHECK-2 vector. SAS cells were transfected with 50 ng of the vector, 10 nM microRNAs, and 0.5 µL Lipofectamine 2000 in 50 µL Opti-MEM (both from Invitrogen, Carlsbad, CA, USA).

Statistical Analysis
Statistical analyses were performed using GraphPad Prism 7 (GraphPad Software, La Jolla, CA, USA) and JMP Pro 15 (SAS Institute Inc., Cary, NC, USA). Dunnet's test were used for multiple group comparisons. For correlation analyses, Spearman's test was applied. A p value less than 0.05 was considered statistically significant. Bar graphs ( Figure 2, Figure 8A,D, Figure S1 and Figure S2) showed mean value and standard error.

Conclusions
In this study, we focused on miR-31-5p and miR-31-3p based on our miRNA signatures. Our functional assays indicated that these miRNAs play an oncogenic role in HNSCC cells. Using in silico database analysis to identify gene targets regulated by miR-31-5p and miR-31-3p, we rapidly identified candidate tumor suppressor genes in HNSCC. Our HNSCC miRNA signature and miRNA-based analyses will provide important insights into the molecular pathogenesis of HNSCC.