Mir-139-5p inhibits glioma cell proliferation and progression by targeting GABRA1

Glioma is an extremely aggressive malignant neoplasm of the central nervous system. MicroRNA (miRNA) are known to bind to specific target mRNA to regulate post-transcriptional gene expression and are, therefore, currently regarded as promising biomarkers for glioma diagnosis and prognosis. The aim of the present study was to examine the pathogenesis and potential molecular markers of glioma by comparing the differential expression of miRNA and mRNA between glioma tissue and peritumor brain tissue. We explored the impact of screened core miRNA and mRNA on cell proliferation, invasion, and migration of glioma. An miRNA expression profile dataset (GSE90603) and a transcriptome profile dataset (GSE90598) were downloaded from combined miRNA-mRNA microarray chips in the Gene Expression Omnibus (GEO) database. Overall, 59 differentially expressed miRNAs (DEMs) and 419 differentially expressed genes (DEGs) were identified using the R limma software package. FunRich software was used to predict DEM target genes and miRNA-gene pairs, and Perl software was used to find overlapping genes between DEGs and DEM target genes. There were 129 overlapping genes regulated by nine miRNAs between target genes of the DEMs and DEGs. The Chinese Glioma Genome Atlas(CGGA) was analyzed in order to identify miRNAs with diagnostic and prognostic significance. MiR-139-5p, miR-137, and miR-338-3p were validated to be significantly linked to prognosis in glioma patients. Finally, we validated that miR-139-5p affected glioma malignant biological behavior via targeting gamma-aminobutyric acid A receptor alpha 1(GABRA1) through rescue experiments. Low miR-139-5p expression was correlated with survival probability and World Health Organization (WHO) grade. MiR-139-5p overexpression inhibited cell proliferation, migration, and invasion of glioma in vitro. GABRA1 was identified as a functional downstream target of miR-139-5p. Decreased GABRA1 expression was related to similar biological roles as miR-139-5p overexpression while upregulation of GABRA1 effectively reversed the inhibition effects of miR-139-5p. These results demonstrate a novel axis for miR-139-5p/GABRA1 in glioma progression and provide potential prognostic predictors and therapeutic target for glioma patients.


Introduction
Glioma is the most prevalent primary malignant neoplasm in the central nervous system, with high recurrence and mortality rates. According to the World Health Organization (WHO) grading system, glioma has been classified into four grades [1]. The higher the grade of the glioma, the worse the prognosis of the patient. Glioblastoma multiforme(GBM), represents the most aggressive glioma subtype, with a five-year survival rate of less than 3% and average survival of less than 12 months [2]. Treatment involves surgical resection followed by adjuvant chemotherapy or concurrent chemoradiotherapy [3]. However, the 3-year survival rate remains at a low level, especially for GBM patients, with the 3-year survival rate only slightly increasing from 2.00%-5.00% to 7.31% [4]. Therefore, researchers have turned their attention to study the potential molecular mechanisms of GBM, and to find novel prognostic biomarkers for the early diagnosis and monitoring of tumorigenesis and for evaluating prognosis.
MicroRNAs (miRNA) are commonly defined as noncoding RNAs consisting of 18-25 nucleotides that regulate post-transcriptional gene expression via specific binding to target mRNA. MiRNA are known to be involved in the regulation of tumor cellular physiological processes, and the aberrant expression and crossregulation of miRNA and target mRNA can be applied as a potential diagnostic and prognostic biomarker in all grades of glioma. Researchers have begun identifying miRNAs that may act as glioma oncogenes or suppressor genes, which could then be applied as diagnostic and predictive biomarkers, and even act as molecular therapeutic targets against glioma [5,6].
To identify the novel diagnostic and predictive value of miRNAs as biomarkers, bioinformatics tools were used to integrate analysis of differential expression of miRNAs and mRNAs from paired 16 fresh-frozen GBM samples and seven peritumoral tissues. Via integrating the miRNA and mRNA expression in the GBM pairs, we identified a number of predictive miRNAs that can be used as biomarkers for diagnosis and prognosis of GBM, and which could then be used as potential clinical therapeutic targets (flow chart see Fig. 1). Specifically, we identified that miR-139-5p and GABRA1 served as inverse agents in the regulation of malignant phenotypes of glioma cells, we doubted whether the function of miR-139-5p on glioma cells was mediated through its inhibitory impression on GABRA1 expression. We, performed rescue experiments to dispose it. In this study, we firstly provide a series of potential targets for future investigation into the molecular mechanisms involved in glioma and offer insight into miRNAs as biomarkers for disease progression. Furthermore, we explore the possible role of miR-139-5p and GABRA1 in the development of glioma and the underlying molecular regulation mechanism. Our findings may provide new insights for the diagnosis and treatment of glioma.

Microarray data
GSE90603 and GSE90598 were downloaded from the GEO database. A total of 23 fresh-frozen samples were analyzed, including 16 individual samples of GBM and seven individual samples of healthy brain tissues. The expression of miRNA (GSE90603) was detected using the platform GPL21582 (Affymetrix Multispecies miRNA-4 Array). The mRNA dataset(GSE90598) was analyzed by using the platform GPL17692 (Affymetrix Human Gene 2.1 ST Array).

Identification of DEMs and DEGs
The DEMs and DEGs between GBM samples and healthy brain tissues was processed using R software. The Affy R package was used to process the downloaded raw data and the limma R package was then used to identify DEMs and DEGs. The adjusted p values were used to decrease the false positive rate using Benjamini and Hochberg false discovery rate method by default. |log2FC|≥ 2 and adjust p < 0.05 were identified as the DEMs or DEGs between GBM sample and healthy brain tissues.

Gene, function and pathway enrichment analysis
FunRich is a stand-alone software tool used mainly for functional enrichment and interaction network analysis of genes and proteins (http:// www. funri ch. org) [7]. FunRich was used to perform DEM downstream transcription factor analysis and DEM target genes prediction. In addition, GO function and KEGG pathway were analyzed using the FunRich software. The core miRNAs and their target genes regulatory network was constructed and visualized by Cytoscape 3.6.1 software.

Survival analysis and glioma grading analysis
The Chinese Glioma Genome Altas (CGGA) database is a user-friendly web application for data storage and analysis to explore brain tumors datasets over 2,000 samples from Chinese cohorts. This database includes miRNA microarray (198 samples) and matched clinical data and could used for assessing the effect of certain miRNA on glioma patients survival. In addition, the miRNA expression level between grade II and grade IV glioma was analyzed via CGGA database. The median expression values of every DEMs in glioma samples were calculated, based on which they were divided into high (above median)-and low-expression(below median) groups. Integrated with the prognostic outcomes, including the overall survival, obtained from CGGA database, the Kaplan-Meier survival curves of these two groups were plotted and the log-rank tests were used to assess the relationship of gene expression to survival. A p-value < 0.05 was set as the threshold for statistical significance.

Quantitative Real time-PCR (qRT-PCR)
TRIzol Reagent (Thermo, USA) was applied to extract the total RNA from samples according to the manufacture's recommendation. The reverse transcription of 200 ng total RNA was utilizing SuperRT RT reagent Kit (CWBio, China) and quantitative real-time PCR with SYBR PCR Master Mix(CWBio, China). The PCR was set at the initial denaturation of 10 min at 95℃, following with 15 s at 95℃, and 30 s at 60 °C in a total of 40 cycles. All experiments were carried out in triplicate. The miR-139-5p expression was normalized to U6 while the GABRA1 was normalized to actin. The relative expression ratios of genes were calculated by the 2-ΔΔCT method. The primers involved in this assay were shown as follow: hsa-miR-139-5p: TCT ACA GTG CAC GTG TCT CC; U6, F: CTC GCT TCG GCA GCACA, R:AAC GCT TCA CGA ATT TGC GT; actin, F:ACC CTG AAG TAC CCC ATC GAG, R:AGC ACA GCC TGG ATA GCA AC; GABRA1, F:ATG ATG GAG CTC GAG GCA AA, R:AGC TCT GAA TTG TGC TGG GT.

Western blot assay
Total protein in cells was collected, and the concentration of protein wasquantified with Pierce BCA Protein Assay Kit. After segregated by 10% sodium dodecyl sulphate-polyacrylamide gel electrophoresis (SDS-PAGE) gel, extracted proteins were transferred to polyvinylidene difluoride (PVDF) membranes and blocked in Tris-buffered saline-Tween with 5% skim milk at room temperature for 60 min. Western blot analysis was performed according to standard procedures. Primary antibodies were GABRA1 (1:1000, proteintech, USA), β-actin (1:5000, proteintech, USA), and corresponding secondary antibodies were anti-mouse and anti-rabbit(proteintech, USA). X-ray film was pressed in dark. Finally, the proteins were detected by enhanced chemiluminescence (ECL), followed by expose about 60 s for scanning and measuring.

EdU assay
The transfected glioma cells were seeded onto the sixwell tissue-culturing plates (1 × 10 5 cells per well). Following the indicated treatments an EdU Apollo-567 assay kit (RiboBio, China) was utilized to test cell proliferation, with nuclear EdU and DAPI staining visualized under a fluorescent microscope.

Cell apoptosis assay by flow cytometry
Apoptotic cells were tested by apoptosis detection kit (KeyGen BioTECH, China) according to the manufacturer's instructions. Cells under log phase were obtained. 500 μl binding buffer was used to resuspend cells after cells were collected. Cells were added with isometric 5 μl Annexin V-APC and 5 μl propidium iodide and incubated for 10 min at room temperature in the dark. Cells were tested using the flow cytometry.

Wound-healing assay
At 6 h post-transfection, the cells were digested, centrifuged and re-suspended in FBS-free culture. The concentration of cells was adjusted to 5 × 10 5 cells/ml. Then cell layers were scratched with a 100 μl sterile pipette tip. After removing cell culture medium and suspension cells and cell debris, each well was added with serumfree medium and stored in incubator for 24 and 48 h. Cell migration area was then viewed and photographed after incubation for 24 h. Meanwhile, scratch test was performed to evaluate the difference in cell healing ability according to the migration area.

Transwell assay
The Matrigel (BD, USA) melted at 4℃ overnight. 100 μl of diluted matrigel was then added in the chamber. Afterwards, 200 μl of serum-free medium was added to the upper chamber; meanwhile, 500 μl of 10% FBS DMEM was added to the lower chambers, and 2 × 10 5 collected cells in total were planted in the upper ones and cultivated in the incubator for another 48 h. Subsequently, the invading chamber was taken out, and cells on the polycarbonate membrane were fixed with 4% paraformaldehydel, followed by staining with 0.1% crystal violet. Three random fields were selected, and invaded cells were counted under a microscope. The experiments were carried out in triplicate.

Statistical analysis
All data are presented as mean ± SD of three independent experiments. Comparisons between the quantitative data were made using Student's t test, with p < 0.05 considered statistically significant. Survival rate was delineated using Kaplan-Meier method. Statistical and graphical analyses were performed with the use of GraphPad Prism 5 software.

Identification of differentially expressed miRNAs (DEMs) and differentially expressed genes (DEGs) in GBM
Using the GSE90603 dataset, 59 DEMs were identified in the GBM samples compared with the peritumoral tissues, of which 37 were upregulated and 22 were downregulated (Additional file 1). The 20 DEMs with the lowest p-value are presented in Table 1. In total, 419 DEGs were obtained in the GSE90598 dataset, and among these, 77 were upregulated and 342 were downregulated (Additional file 2). The top 20 DEGs with the lowest p-value are presented in Table 2. Heatmaps of DEMs and DEGs are shown in Fig. 2.

Transcription factor (TF) and gene ontology (GO) analyses of DEMs
TF and GO analyses of DEMs were performed by Fun-Rich software. The top ten most commonly expressed TFs associated with DEMS were RORA, RREB1, ZFP161, FOXA1, MEF2A, NKX6-1, POU2F1, SP4, SP1, and EGR1 (Fig. 3a). GO enrichment analysis of the top ten DEMs was also performed using FunRich software. The result revealed that the majority of DEMs in the biological process (BP) category were enriched in "signal transduction", "cell communication", "nucleoside, nucleotide and nucleic acid metabolism", and "transport" (Fig. 3b and Additional file 3). The GO cellular component (CC) enrichment analysis revealed that the majority of DEMs in this category were enriched in "nucleus", "cytoplasm", "lysosome", and "golgi apparatus" (Fig. 3c and Additional file 4). In the molecular function (MF) category, the majority of DEMs were associated with "transcription factor activity", "transcription regulator activity", "protein serine/threonine kinase activity", "ubiquitin-specific protease activity", "receptor signaling, and complex scaffold activity" (Fig. 3d and Additional file 5).

DEMs target gene prediction
The potential target genes (PTGs) of each DEM were obtained by using FunRich software. All 59 DEMs were selected for PTG prediction and 34 DEMs were identified. The FunRich tool generated 5147 PTGs for 34 DEMs. We then took the overlapping gene intersection of PTGs and DEGs, and analyzed the correlation of downregulated miRNAs with upregulated genes, upregulated miRNAs, and downregulated genes. There were nine core DEMs and 129 overlapping genes with a targeted relationship (Table 3).

MiRNA-gene regulatory network construction
According to the information obtained from the nine core DEMs and 129 overlapping target genes, core DEMs and target gene regulatory networks were constructed and visualized using Cytoscape software (Fig. 4). Overall, the data indicate that the nine core DEMs and 129 overlapping target genes may play important roles in the diagnosis and prognosis of glioma.

GO and Kyoto Encyclopedia of Genes and Genomes(KEGG) pathway enrichment analysis of overlapping target genes
To explore roles that overlapping target genes paly in biological functions, GO and KEGG pathway analysis were performed. The overlapping genes were significantly enriched in BP, including "synaptic vesicle cycle", "regulation of membrane potential", "glutamate receptor signaling pathway", "modulation of chemical synaptic transmission" and so on ( Fig. 5a and Additional file 6). For MF, The overlapping target genes were mainly related to "ion gate channel activity", "gated channel activity", "ion channel activity", "substrate-specific channel activity" and so on ( Fig. 5a and Additional file 6). In addition, CC analysis indicated that the overlapping target genes were involved in "presynapse", "glutamatergic synapse", "synaptic membrane", "presynaptic membrane" and so on (Fig. 5a Additional file 6). KEGG analysis revealed that the overlapping target genes were mostly enriched in "glutamatergic synapse", "dopaminergic synapse", "retrograde endocannabinoid signaling", "nicotine addiction" and "synaptic vesicle cycle" (Fig. 5b and Additional file 7).
To reveal the correlation of core miRNA expression profiles with glioma grade, we compared core miR-NAs distribution data of grade II-IV cases obtained in the CGGA database. The expression of miR-137 and    miR-338-3p was positively correlated with glioma grade (Fig. 6d, e), whereas the expression of miR-139-5p was negatively correlated with glioma grade (Fig. 6f ).

GABRA1 as the functional target of miR-139-5p
Combining the literature search and the bioinformatics analyses, miR-139-5p and its potential target, GABRA1, were selected for further analysis. We next aimed to ascertain the underlying molecular mechanisms by which miR-139-5p exerted its tumor suppressing roles in glioma. Hence, an miR-139-5p mimic and inhibitor were respectively transfected into the U251 and U87 cell lines. The efficacy of transfection was determined by RT-PCR, with results indicating that the relative expression of miR-139-5p was significantly increased in U251 and U87 cells after transfection of the miR-139-5p mimic (Fig. 7a).
Conversely, miR-139-5p expression was dramatically downregulated after transfection with the miR-139-5p inhibitor (Fig. 7a). Both the previous bioinformatics analyses and online miRNA database analysis (Target Scan) indicated that GABRA1 could be a target of miR-139-5p, according to the putative target sequence in the GABRA1 3'Untranslated Regions (UTR) (Fig. 7b). RT-PCR and western blot analyses were carried out to further identify whether GABRA1 expression was indeed regulated by miR-139-5p. The results indicated that the mRNA and protein levels of GABRA1 were significantly increased after transfection with the miR-139-5p inhibitor, while the levels were reduced after transfection with the miR-139-5p mimics, in U251 and U87 cells (Fig. 8a-e). In addition, the results of RT-PCR assays indicated that the overexpression of GABRA1 had no effect on the expression of miR-139-5p (Fig. 8f ). These results suggested that miR-139-5p could directly target GABRA1 and negatively regulate its expression.

MiR-139-5p regulates glioma cell migration and invasion through targeting of GABRA1.
We next performed transwell assays and wound healing assays to evaluate the effects of miR-139-3p on the migratory and invasive capabilities of U251 and U87 cells. Both transwell and wound healing assays indicated that miR-139-5p significantly repressed U251 and U87 cell migration and invasion (Figs. 9 and 10). Lossof-function experiments suggested that the miR-139-5p inhibitor elevated U251 and U87 cell migration and invasion (Figs. 9 and 10).
Next, we explored whether GABRA1 mediated the effects of miR-139-5p on the cell migration and invasion of U251 and U87 cells. Specifically, rescue experiments were carried out by co-transfecting the miR-139-5p mimic/inhibitor with or without GABRA1 overexpression (OE) and silencing (Si). The miR-139-5p inhibitor resulted in a significant increase in the migration and invasion, while miR-139-5p mimic resulted in decreased migration and invasion, of glioma cells. OE-GABRA1 remarkably rescued the suppressive effects of the miR-139-5p mimic on the invasive abilities of U251 and U87 cells (Figs. 9 and 10). In addition, cotransfection of the miR-139-5p inhibitor and OE-GABRA1 notably promoted migration and invasion compared with transfection of the miR-139-5p inhibitor alone. Cotransfection of the miR-139-5p inhibitor and Si-GABRA1 significantly inhibited the migration and invasion of U251 and U87 cells (Figs. 9 and 10). Taken together, our data indicated that miR-139-5p modulates glioma cell migration and invasion by targeting GABRA1.

MiR-139-5p regulates glioma cell proliferation and apoptosis through targeting of GABRA1
Cell proliferation and apoptosis experiments were then performed in U251 and U87 cell lines. The results of the quantitative measurement of EdU-positive cells indicated that downregulation of miR-139-5p promoted the proliferation, while upregulation of miR-139-5p inhibited the proliferation, of U251 and U87 cells (Fig. 11). Flow cytometric analysis showed that the cell apoptosis rate was notably higher in the miR-139-5p mimic group while the  19:213 rate was lower in the miR-139-5p inhibitor group compared with the NC group, indicating that miR-139-5p induced the apoptosis of U251 and U87 cells (Fig. 12). Since GABRA1 is directly targeted by miR-139-5p, we speculated that GABRA1 mediates the function of miR-139-5p. To test this hypothesis, we overexpressed GABRA1 in the miR-139-5p mimic group and found that GABRA1 overexpression rescued the miR-139-5pinduced inhibition of proliferation and increase in apoptosis (Figs. 11 and 12). Cotransfection of the miR-139-5p mimic and Si-GABRA1 significantly suppressed proliferation and promoted apoptosis, compared with transfection of miR-139-5p mimic alone (Figs. 11 and 12). The results revealed that the restoration of GABRA1 markedly reversed miR-139-5p-mediated glioma inhibitory effects.

Discussion
MiRNAs are a group of short and non-coding RNA molecules that influence the biological characteristics of GBM cells [8,9]. Mounting evidence has shown that miRNAs can act as cancer-suppressing factors as well as oncogenes [10][11][12]. MiRNAs play important roles in the tumorigenesis and progression, invasion, and metastasis of GBM by upregulating or downregulating cancer-related gene expression levels, and, therefore, have become promising biomarkers in the diagnosis and prognosis glioma [13].
In this work, bioinformatics technology was used to identify candidate miRNA biomarkers of glioblastoma via analysis of microarray data. As a result, 59 DEMs and 419 DEGs were identified in the GBM samples compared with healthy brain tissues.  The prediction and GO enrichment analyses performed on the top ten DEMs showed that the majority of DEMs were enriched in signal transduction (BP) and transcription factor activity (MF), and were primarily located in the nucleus (CC). To explore the interactions between DEMs and their corresponding target genes, we performed GO enrichment and KEGG pathway analyses using the list of targeted genes. GO term analysis indicated that the target genes were mainly involved in synaptic vesicle cycle and regulation of membrane potential. KEGG analysis of the target genes identified that they were involved in synaptic relevant pathways, including the glutamatergic synapse, dopaminergic synapse, and synaptic vesicle cycle signaling pathways. The results also indicated that the synaptic related genes aide in the prognosis of patients with glioma, which can improve objectivity in clinical judgment.
Subsequently, the integration analysis of the miRNAtarget gene regulatory pairs and DEGs revealed that 129 overlapping genes were regulated by nine DEMs. Further analysis revealed that miR-137, miR-139-5p, and miR-338-3p were associated with prognosis in patients with all grades of glioma.
MiRNA-137 is an miRNA that is widely expressed in the central nervous system, and is particularly specific to hippocampal tissue [14]. It has been reported that the biological function of miRNAs may be related to synaptic plasticity and transmission [15,16]. Previous studies have indicated lower expression of miRNA-137 in tumor tissues, as compared with healthy tissues, in gastric cancer, colon cancer, and oral squamous cell carcinoma, suggesting that miRNA-137 acts as a suppressor [17,18].
It has been reported that miR-338-3p functions as a tumor suppressor in several types of cancers [19][20][21]. Overexpression of miR-338-3p attenuated malignant biological behaviors of cells in gastric cancer and non-small cell lung cancer [19][20][21]. In the central nervous system, expression of miR-338-3p increases significantly as the dentate gyrus matures, and peaks in mature neurons [22]. Clinical data has also shown that low expression of miR-338-3p corresponds to a decrease in overall survival and disease-free survival [22].
Here, we found differing conclusions from previous research based on data from the CGGA. Lower expression of miRNA-137 and miR-338-3p was found in glioma tissue than in peritumoral tissue, suggesting that these two miRNAs may be tumor suppressor genes, consistent with previously reported in the literature. However, in the CGGA data, high expression of miRNA-137 and miR-338-3p was also associated with poor overall survival in glioma patients. The expression of miRNA-137 and miR-338-3p were positively correlated with glioma grade, suggesting that both miRNAs might be cancerpromoting factors.