High miR-3609 expression is associated with better prognosis in TNBC based on mining using systematic integrated public sequencing data.

MicroRNAs (miRNAs/miRs) are small endogenous RNAs that regulate gene expression post-transcriptionally. Abnormal miR-3609 expression is associated with the occurrence of pancreatic cancer, glioma and other diseases, such as polycystic ovary syndrome. However, the prognostic potential of miR-3609 has been reported in breast cancer. Thus, the present study aimed to investigate the differential expression and prognostic value of miR-3609 in patients with breast cancer from the UALCAN, cBioportal and Kaplan-Meier Plotter databases, respectively. Furthermore, the co-expression genes of miR-3609 in breast cancer were investigated using data from the LinkedOmics database, and functional enrichment analysis was performed using the LinkInterpreter module in LinkedOmics. The co-expression gene network was constructed using the Search Tool for the Retrieval of Interacting Genes/Proteins database, and the cytoHubba plug-in was used to identify the hub genes, which were visualized using Cytoscape software. The prognoses of the hub genes were performed using the Kaplan-Meier Plotter database. The Cell Counting Kit-8 and cell cycle assays were performed to confirm the functions of miR-3609 mimics transfection in MDA-MB-231 cells. Survival analysis using the Kaplan-Meier Plotter database demonstrated that high miR-3609 expression in triple-negative breast cancer (TNBC) was associated with a better prognosis. Furthermore, the experimental results indicated that high miR-3609 expression inhibited the proliferation of TNBC cells and induced cell cycle arrest of TNBC cells in the G0/G1 phase. Taken together, the results of the present study suggest that miR-3609 plays a vital role in mediating cell cycle arrest and inhibiting the proliferation of TNBC cells. Copyright: © Ding et al.

Introduction

Breast cancer is the most common cancer, accounting for 11% of all diagnosed cancer sites, and is the leading cause of cancer-associated mortality worldwide, accounting for 6.6% of all diagnosed cancer sites (1). Triple-negative breast cancer (TNBC) is associated with high metastatic risk and low overall survival (2). Treatment options are particularly limited for patients with TNBC, which is associated with earlier time to recurrence, higher risk of distant metastasis and worse prognosis after recurrence (3).

MicroRNAs (miRNAs/miRs) are small non-coding RNAs that are ~22 nucleotides in length, which play an important role in the post-transcriptional regulation of mRNA (4). miRNAs are vital to several aspects of cell biology, including apoptosis, proliferation, differentiation and cell cycle (5). Over the past decade, increasing evidence suggest that deregulation of miRNAs is a crucial part of tumor formation, maintenance, metastasis and immune tolerance, including breast cancer (6). Specific miRNA expression profiles are associated with distinct breast cancer subtypes (6). For example, miR-218, miR-342, miR-135b, miR-217, miR-299 and miR-190 are closely associated with estrogen receptor (ER) positive breast cancer, while the expression levels of miR-520f-520C, miR-377, miR-527-518a and miR-520g are closely associated with progesterone receptor (PR) positive breast cancer (7). It has been reported that abnormally expressed or mutated miRNAs play an important role in the pathogenesis of breast cancer (6). Our previous study demonstrated that miR-3609, which is abnormally expressed in endometrial dysplasia, oral squamous cell carcinoma and pancreatic cancer (8-10), was significantly downregulated in breast cancer (11). The results also demonstrated that restoration of miR-3609 expression sensitized breast cancer to doxorubicin by blocking the PD-L1 immune checkpoint in vitro (11). The present study extended the research in breast cancer based on public databases to perform a variety of bioinformatics analyses. This study aimed to determine the potential functions and distinct prognostic values of miR-3609 in breast cancer by analyzing the expression and mutation of the miR-3609 in patients with breast cancer.

Materials and methods

UALCAN

The UALCAN (http://ualcan.path.uab.edu) is a comprehensive and interactive web resource for analyzing genomics data from The Cancer Genome Atlas (TCGA) database (https://www.cancer.gov/), which was used to determine the expression of miR-3609 in breast cancer. It allows researchers to analyze the relative expression of potential genes of interest in tumors and normal samples, and in various tumor subgroups (12).

cBioPortal

The cBioPortal (https://www.cbioportal.org) is an open-access resource for interactive exploration of multidimensional cancer genomics datasets (13). Genetic alterations of miR-3609 was obtained from cBioPortal, based on TCGA database.

Kaplan-Meier plotter

The Kaplan Meier Plotter database (https://kmplot.com/analysis) is able to assess the effect of 54 k genes (mRNAs, miRNAs and proteins) on the survival of patients with 21 different types of cancer. Sources for the databases include Gene Expression Omnibus (GEO; https://www.ncbi.nlm.nih.gov/) and TCGA (14). The overall survival of patients with breast cancer was analyzed using the Kaplan-Meier Plotter database.

LinkedOmics

LinkedOmics (http://www.linkedomics.org/login.php) is a publicly available portal that contains multi-omics and clinical data from 32 types of cancer, and 11,158 patients from TCGA project (15). It also includes mass spectrometry-based proteomics data generated by the Clinical Proteomics Tumor Analysis Consortium for TCGA breast, colorectal and ovarian tumors.

The differentially expressed genes associated with miR-3609 were screened from TCGA BRCA cohort (n=755), using the LinkFinder module in the LinkedOmics database. Spearman's correlation coefficient analysis was performed to assess the correlation between the miR-3609 expression and genes differentially expressed in breast cancer.

Gene Ontology (GO; http://geneontology.org/) and Kyoto Encyclopedia of Genes and Genomes (KEGG; https://www.kegg.jp/) pathway enrichment analyses of the differentially expressed genes were performed using the LinkInterpreter module, the results of which were signed and ranked using the gene set enrichment analysis (GSEA) tool in the LinkedOmics.

Search Tool for the Retrieval of Interacting Genes/Proteins (STRING)

STRING (https://string-db.org) is a public database that predicts protein-protein interactions (PPIs), and aims to achieve a comprehensive and objective global network and present them with a unique set of computational predictions (16).

To establish a PPI network, the present study screened out co-expressed genes with interaction scores >0.4. Cytoscape software (version 3.8.0; https://cytoscape.org) which is an open source software platform for visualizing molecular interaction networks and biological pathways and integrating these networks with annotations, gene expression profiles and other state data, was used to visualize the PPI network, and hub genes were identified using cytoHubba plug-in (version 0.1; https://cytoscape.org) and the degree topological algorithm.

Cell culture and microRNA mimics transfection

The TNBC cell line, MDA-MB-231 was purchased from the American Type Culture Collection and maintained in DMEM (Hyclone; Cytiva) supplemented with 10% FBS (Gibco; Thermo Fisher Scientific, Inc.), 100 U/ml penicillin G and 100 µg/ml streptomycin (Beijing Solarbio Science & Technology Co., Ltd.), at 37˚C in a humidified 5% CO2 atmosphere.

miR-3609 mimics (40, 80 and 160 nM; 5'-CAAA GUGAUGAGUAAUACUGGCUG-3' and 5'-GCCAGUA UUACUCAUCACUUUGUU-3'), nonspecific miRNA mimics (5'-UUCUCCGAACGUGUCACGUTT-3' and 5'-ACGUGACACGUUCGGAGAATT-3') were synthesized by Shanghai GeneChem Co., Ltd., and transfected into MDA-MB-231 cells using Lipofectamine® RNAiMax reagent (Invitrogen; Thermo Fisher Scientific, Inc.), according to the manufacturer's instructions. After transfection at 37˚C for 24 h, miR-3609 expression in MDA-MB-231 cells transfected with 80 nM miR-3609 mimics increased the most (Fig. S1), MDA-MB-231 cells transfected with 80 nM miR-3609 mimics were used for subsequent analyses. Non-specific miRNA mimics, synthesized by Shanghai GeneChem Co., Ltd., were used as the negative controls. Transfected cells were used for subsequent experimentation 48 or 72 h post-transfection.

Cell proliferation assay

MDA-MB-231 cells were transfected with either miR-3609 mimics or a non-specific miRNA mimic (both at 80 nM), using Lipofectamine RNAiMax reagent (Invitrogen; Thermo Fisher Scientific, Inc.), according to the manufacturer's instructions. Transfected MDA-MB-231 cells were seeded into 96-well plates at a density of 2x103 cells/well and incubated for 24, 48, 72, 96 or 120 h at 37˚C in DMEM (Hyclone; Cytiva). Cell viability was assessed via the Cell Counting Kit-8 (CCK-8) assay (Suzhou Yuheng Biotechnology Co., Ltd.) for 2 h.

Cell cycle assay

Following transfection with either miR-3609 mimics or a non-specific miRNA mimic for 48 h, MDA-MB-231 cells were collected after centrifugation with 172.2 x g for 5 min at room temperature, and washed twice with pre-cooled PBS. Subsequently, 1.0x104 cells were fixed with 500 µl of 70% cold ethanol overnight at 4˚C, followed by addition of 500 µl of propidium iodide (PI) solution (11.6 µg/ml PI, 2 mg/ml RNaseA in PBS; Beijing Solarbio Science & Technology Co., Ltd.), and the mixture was incubated at room temperature for 30-60 min, in the dark. DNA content analysis was performed using the FACS Calibur instrument (BD FACSCantoII™) and CellQuest software (modfit 5.1; BD Biosciences). The experiments were performed in triplicate.

Statistical analysis

Statistical analysis was performed using GraphPad Prism 8.0 software (GraphPad Software, Inc.). Data are presented as the mean ± SD of at least three independent experiments. Unpaired Student's t-test was used to compare differences between two groups (Figs. 1A, and 4C and D). One-way ANOVA followed by Tukey's post hoc test was used to compare differences between multiple groups (Fig. S1). P<0.05 was considered to indicate a statistically significant difference.

Figure 1

Expression and alteration of miR-3609 in breast cancer. (A) miR-3609 expression in normal and breast cancer samples. (B and C) Genetic alterations of miR-3609, using cBioPortal. miR, microRNA; TCGA, The Cancer Genome Atlas.

Figure 4

Results of the cell cycle and Cell Counting Kit-8 assays. (A) Flow cytometric analysis of MDA-MB-231 cells transfected with 80 nM non-specific miRNA mimic and (B) miR-3609 mimics. (C) Percentages of cells in the G0/G1, S and G2/M phases (n=3, mean ± SD). (D) Transfection with miR-3609 mimics inhibited the proliferation of MDA-MB-231 cells. *P<0.05. miRNA/miR, microRNA; OD, optical density.

Results

Expression and alteration of miR-3609 in breast cancer

The present study determined the expression of miR-3609 between normal and primary breast tumor samples in TCGA database (Fig. 1A). The results revealed no significant differences between the two groups. The genetic alterations of miR-3609 in 1,108 breast invasive carcinoma samples in TCGA database were detected, using cBioPortal. The results demonstrated that 15 cases had miR-3609 alterations (1.4%) (Fig. 1B and C).

miR-3609 expression is associated with overall survival of patients with breast cancer

After investigating the overall survival of miR-3609 in breast cancer using the Kaplan-Meier Plotter database, the results demonstrated no significant differences in overall survival for miR-3609 among all patients with breast cancer from TCGA database. However, high miR-3609 expression was associated with better overall survival in patients with TNBC (Fig. 2A and B).

Figure 2

Prognosis of miR-3609 in breast cancer and TNBC. (A) Kaplan-Meier curve for overall survival in all patients with breast cancer. (B) Kaplan-Meier curve for overall survival in patients with TNBC. miR, microRNA; TNBC, triple-negative breast cancer; HR, hazard ratio.

Co-expression genes correlated with miR-3609 in breast cancer

The present study analyzed the co-expressed genes of miR-3609 in 755 patients with breast cancer, using the LinkedOmics database. The volcano plot depicts genes positively and negatively correlated with miR-3609 (Fig. 3A). The top 50 significant gene sets positively or negatively correlated with miR-3609 are depicted in the heat maps, respectively (Fig. 3B and C).

Figure 3

Genes correlated with miR-3609 in breast cancer. (A) Spearman's correlation coefficient analysis was performed to assess the correlation between miR-3609 expression and genes differentially expressed in breast cancer. (B and C) Heat maps depicting top 50 genes positively and negatively correlated with miR-3609 expression in breast cancer, respectively. Enriched Gene Ontology terms of miR-3609 co-expression genes in breast cancer, based on (D) biological processes, (E) cellular components and (F) molecular functions, and (G) KEGG pathway analysis. miR, microRNA; KEGG, Kyoto Encyclopedia of Genes and Genomes; FDR, false discovery rate.

GO and KEGG analyses of genes associated with miR-3609 in breast cancer

GO analysis demonstrated that differentially expressed genes associated with miR-3609 were mainly located in ‘spliceosomal complex’, ‘ribosome’, ‘condensed chromosome’ and ‘exoribonuclease complex’, and these cellular components may have participated in ‘RNA splicing’, ‘rRNA metabolic process’, ‘spinal cord development’, ‘G0-G1 transition’ or other biological processes. They acted as structural constituents of the ribosome and mRNA binding (Fig. 3D-F). KEGG pathway analysis demonstrated that the functions of the associated genes were primarily enriched in ribosome and spliceosome (Fig. 3G).

To confirm G0-G1 transition for biological processes of miR-3609, the cell cycle assay was performed. The results demonstrated that the proportion of MDA-MB-231 cells transfected with miR-3609 mimics in the G0/G1 phase was significantly higher (P=0.032), and that in the G2/M phase was significantly lower (P=0.018) compared with the control group (Fig. 4A-C), suggesting that miR-3609 induces cell cycle arrest of TNBC cells in the G0/G1 phase.

The CCK-8 assay was performed to assess the viability of MDA-MB-231 cells. The results demonstrated that the viability of MDA-MB-231 cells transfected with miR-3609 mimics was significantly lower compared with MDA-MB-231 cells transfected with non-specific miRNA mimic (P=0.0285), suggesting that high miR-3609 expression decreases the viability of MDA-MB-231 cells (Fig. 4D).

Construction of co-expression gene PPI network

The 422 significantly co-expressed genes were selected to construct aPPI network using the STRING database (P<0.001), and Cytoscape (cytoHubba plug-in) was used to identify the hub genes, along with the degree topological algorithm. Based on the degree score, the genes with the highest scores [interleukin enhancer binding factor 3 (ILF3), ELAV like RNA binding protein 1 (ELAVL1), heterogeneous nuclear ribonucleoprotein (HNRNP)L, HNRNPK, HNRNPH1 and polypyrimidine tract binding protein 1 (PTBP1)] were identified as potential hub genes (Fig. 5A and B). The gene of ILF3 encodes a double-stranded RNA (dsRNA)-binding protein that complexes with other proteins, dsRNAs, small noncoding RNAs and mRNAs to regulate gene expression and stabilize mRNAs (17). The protein encoded by ELAVL1 is a member of the ELAVL family of RNA-binding proteins that contain several RNA recognition motifs, and selectively bind AU-rich elements found in the 3' untranslated regions of mRNAs (18). Heterogeneous nuclear RNAs which include mRNA precursors and mature mRNAs are associated with specific proteins to form heterogenous ribonucleoprotein (hnRNP) complexes (19). HNRNPL is stably associated with hnRNP complexes, and is likely to play a major role in the formation, packaging, processing, and function of mRNA, along with other hnRNP proteins (20). HNRNPH1 encodes a member of a subfamily of ubiquitously expressed hnRNPs (19). Moreover, both HNRNPK and PTBP1 belong to the subfamily of ubiquitously expressed hnRNPs (21,22).

Figure 5

PPI network of co-expression genes. (A) PPI network. (B) The hub genes (red nodes), including ILF3, ELAVL1, HNRNPL, HNRNPK, HNRNPH1 and PTBP1. PPI, protein-protein interaction; ILF3, interleukin enhancer binding factor 3; ELAVL1, ELAV like RNA binding protein 1; HNRNP, heterogeneous nuclear ribonucleoprotein; PTBP1, polypyrimidine tract binding protein 1.

Prognostic analysis of hub genes in breast cancer

The overall survival of the hub genes in breast cancer was assessed using the Kaplan-Meier Plotter database. The results demonstrated that low HNRNPL expression was associated with a longer overall survival time (P<0.05; Fig. 6).

Figure 6

Prognostic analysis of the hub genes using the Kaplan-Meier Plotter database. (A) ILF3. (B) HNRNPH1. (C) PTBP1. (D) ELAVL1. (E) HNRNPK. (F) HNRNPL. ILF3, interleukin enhancer binding factor 3; ELAVL1, ELAV like RNA binding protein 1; HNRNP, heterogeneous nuclear ribonucleoprotein; PTBP1, polypyrimidine tract binding protein 1.

Discussion

Breast cancer remains the most common malignancy and the main cause of cancer-associated mortality in women worldwide, with incidence and mortality rates of 24.2% and 15.0%, respectively (23-25). TNBC is a heterogeneous breast cancer subtype, where the ER, PR and human epidermal growth factor receptor 2 are negatively expressed (26-28). TNBC represents 12-17% of all breast cancer cases (29), and has a more aggressive clinical course, with greater metastatic potential and poorer prognosis as demonstrated by the higher relapse (33.9% vs. 20.4%) and lower survival rates (42.2% vs. 28%) compared with patients with other breast cancers (30,31). Thus, it is important to identify novel biomarkers to predict the prognosis of patients with breast cancer, particularly for TNBC.

Increasing evidence suggest that non-coding RNAs are active participants in multiple stages of tumor immunity (6,32-35). Non-coding RNAs, including miRNAs, long non-coding RNAs and circular RNAs, differentially regulate multiple cellular processes in development and diseases via a variety of gene-regulation mechanisms (36). miRNAs are small non-coding RNA molecules that are ~22 nucleotides in length, which play important regulatory roles in several aspects of cell activities, such as cell differentiation, apoptosis and metabolism (37). Previous studies have reported that miRNAs, such as miR-1296 and miR-133a, are closely associated with the occurrence, malignant metastasis and poor prognosis of breast cancer, and overexpression of these miRNAs can inhibit the proliferation of breast cancer cells (38-40). It has also been reported that transfection with miR-3609 mimics can decrease cyclin-dependent kinase 1 expression in MDA-MB-231 cells (41). Our previous study demonstrated that miR-3609 can target the 3'-untranslated region of programmed death-ligand 1 (PD-L1) and effectively inhibit PD-L1 expression in breast cancer cells (11). The present study confirmed that high miR-3609 expression is associated with better overall survival in patients with TNBC. Thus, miR-3609 may be used as a potential prognostic biomarker for TNBC. The results of the present study demonstrated no significant difference in overall survival for miR-3609 expression among all patients with breast cancer. This may have been due to varied miR-3609 expression among diverse subgroups and subtypes (42). The present study also investigated the alterations of miR-3609 in breast cancer using the cBioPortal database. A total of 15 cases of miR-3609 alterations were found in the cBioPortal database. The co-expression genes correlated with miR-3609 expression in breast cancer was also assessed, using the LinkedOmics database. GSEA suggested that co-expression genes participated in RNA splicing, rRNA metabolic process and G0-G1 transition via GO and KEGG pathway analyses. The results of the CCK-8 and cell cycle assays demonstrated that high miR-3609 expression inhibited the proliferation of TNBC cells and induced cell cycle arrest of TNBC cells in the G0/G1 phase. Bioinformatics analysis was performed to identify hub genes by gene analysis of tumor data from public databases. ILF3, ELAVL1, HNRNPL, HNRNPK, HNRNPH1 and PTBP1 were identified as potential hub genes. Subsequently, the prognostic values of these hub genes were determined using the Kaplan-Meier Plotter database. Taken together, the results of the present study suggest that miR-3609 may be used as a potential prognostic marker and a therapeutic target for patients with TNBC.

The present study is not without limitations. First, the effects of miR-3609 on the protein expression levels of ILF3, ELAVL1, HNRNPL, HNRNPK, HNRNPH1 and PTBP1 were not investigated due to the lack of antibodies at the time of the study. Secondly, only one TNBC cell line was used to assess the effect of miR-3609 on cell cycle arrest. Thus, prospective studies will use other TNBC cell lines and in vivo models to verify the results presented here. In addition, further studies are required to determine the specific mechanisms of miR-3609 involved in the cell cycle of TNBC cells.

In conclusion, the present study integrated public sequencing data to guide the research of miR-3609 in breast cancer. The results presented here confirm that high miR-3609 expression is associated with a longer overall survival time in patients with TNBC. In addition, high miR-3609 expression inhibited the proliferation of TNBC cells and induced cell cycle arrest. Thus, miR-3609 may be used as a potential biomarker for TNBC diagnosis, and miRNAs-based therapeutic approaches may be an attractive alternative option for patients with TNBC.