FKBP-related ncRNA-mRNA axis in breast cancer CURRENT STATUS: REVIEW

Background Breast cancer (BC) is a disease with morbidity ranking the first of women worldwidely. FK506-binding protein (FKBP) family has been demonstrated to possess various functions by interacting with different molecular targets in BC. However, a comprehensive ncRNA-mRNA regulatory axis of FKBP has not yet been reported. Methods FKBP related miRNAs were obtained from miRWalk database. Then, potential lncRNAs, transcription factors as well as mRNAs of screened differentially expressed miRNAs (DE-miRNAs) were analysed by using LncBase v.2, miRGen v3 and miRWalk database. Additionally, differential expression and prognostic analysis of lncRNAs were evaluated using TANRIC database. Next, GO annotation and KEGG pathway analysis were processed using DAVID database. Protein-Protein Interaction (PPI) network was established and hub genes were identified using STRING database. Finally, differential expression and prognostic analysis of hub genes were further conducted using UALCAN and bc-GenExMiner v4.2 database, respectively. Results Eleven DE-miRNAs, consisting of four FKBP4 related DE-miRNAs and seven FKBP5 related DE-miRNAs, were screened. 482 predicted lncRNAs were found for DE-miRNAs. Then, expression and prognostic results of nine of top twenty lncRNAs of BC were significantly identified. LINC00662 and LINC00963 expression were significantly associated with patients’ overall survival (OS). Then, nine potential upstream transcription factors were identified in motifs of DE-miRNAs. 320 target genes were identified for GO annotation and KEGG pathway analysis, which were mainly enriched in cysteine-type endopeptidase activity involved in apoptotic process. Construction and analysis in PPI network showed that RAB7A was selected as a hub gene with the toppest connectivity scores. Differential expression analysis of nine in top ten hub genes of BC were significantly identified. RAB7A and ARRB1 expression were significantly related with BC patients’ OS. Conclusions In current study, we firstly

established a predicted FKBP-related ncRNA-mRNA regulatory network, thus exploring a comprehensive interpretation of molecular mechanisms and providing potential clues in seeking novel therapeutics for BC. In the future, much more experiments should be conducted to verify our findings.

Background
Breast cancer (BC) is the most widespread and deadly non-cutaneous tumor in worldwide women [1]. Early detection and comprehensive treatments, which consist of surgery, radiation, chemotherapy, endocrine therapy and targeted therapy, have significantly improved the prognosis in BC patients. However, BC is a heterogeneous disease of various different genetical, pathological, and clinical subtypes [2]. Even though intensive efforts have been made in both basic researches and clinical studies, it's still necessary to find more reliable markers to further improve therapeutics for BC patients.
FK506-binding protein (FKBP) family in Homo sapiens genomes has been found to target on various pathways in embryology, stress reaction, heart function, tumorigenesis and neuronal function [3]. In breast cancer, FKBP4 and FKBP5 are most extensively studied proteins among identified human FKBPs, which are demonstrated to interact with Hsp90 to affect steroid hormone receptor function [4]. MicroRNAs (miRNAs) are a cluster of small noncoding RNAs consisting of 20 to 24 nucleotides, regulating targeted gene expression by binding to several selective messenger RNAs (mRNAs) [5]. MiRNAs are also found to exert pivotal roles in the genesis and development of BC. For instance, the miRNA let-7's ability to restrain the expression of metastatic genes could be compromised when long non-coding RNA (lncRNA) H19 expression is upregulated, leading to high expression of c-Myc and activating migration and invasion of BC cells [6]. Despite many researches on miRNA expression and function of BC have been conducted, a comprehensive analysis of FKBP-related ncRNA-mRNA regulatory network via clinical information of BC is still lacking.
Construction of predicted ncRNA-mRNA axis contributing to BC might unravel the potential molecular mechanisms underlying processes of miRNAs' impact on BC.
Herein, a total of four FKBP4 related differentially expressed miRNAs (DE-miRNAs) and seven FKBP5 related DE-miRNAs were screened. Subsequently, expression and prognostic analytic result of nine of the top twenty lncRNAs in BC were significantly identified.
LINC00662 and LINC00963 expression were significantly associated with patients' overall survival (OS). Then, nine potential upstream transcription factors (TFs) were identified in motifs of DE-miRNAs. 320 target genes were selected for GO annotation and KEGG pathway analysis. Construction of Protein-Protein Interaction (PPI) network showed that RAB7A was recognized as a hub gene with the toppest connectivity scores. Differential expression analysis of nine in top ten hub genes of BC were significantly identified. RAB7A and ARRB1 expression were significantly associated with patients' OS. Finally, a potential FKBP-related ncRNA-mRNA regulatory axis contributing to the onset and development in BC was successfully achieved.

Verification of FKBPs Expression Levels
The mRNA expressions of twelve FKBPs were further verified using GEPIA, which is a recently developed database for analyzing the RNA sequencing expression data of carcinoma and adjacent samples in the TCGA and the GTEx projects [7].

Screening of Potential miRNAs of FKBPs and Targeted Genes of miRNAs
Both screened miRNAs of FKBP4 and FKBP5 and screened gene targets of miRNAs were conducted using miRWalk database, mainly using for experimentally verified miRNA-target interactions [8], specifically the screened miRNAs and genes were generated by intersection of miRDB and miRTarBase.

Screening Of Potential Lncrnas And Transcription Factors Of Mirnas
Predicted lncRNAs of screened miRNAs were all generated by using LncBase v.2, which is a tool used mainly for discovering the connection between miRNAs and lncRNAs [9]. The upstream TFs of screened miRNAs were analyzed by using miRGen v3, which is mainly conducted for discovering the connection between miRNAs and TFs [10].

Validation Of Lncrna Differential Expressions And Prognostic Functions
Both expression levels of top twenty lncRNAs in different subtypes and their prognostic roles of overall survival of BC patients were further validated by using The Atlas of ncRNA in Cancer (TANRIC), an open-access database for interactive exploration of lncRNAs of various cancer [11]. lncRNAs with |log2FC|>2 and P < 0.05 were regarded as statistically significant.

Go Annotation And Kegg Pathway Analysis
To better understand those screened candidate targeted genes, DAVID database was used to perform GO annotation and KEGG pathway analysis. The top ten enriched GO items were all listed in Fig. S4A-S4C. BP analysis revealed that candidate targeted genes of screened DE-miRNAs were significantly enriched in negative regulation of cysteine-type endopeptidase activity involved in apoptotic process (Fig. S4A). As for CC analysis, genes were significantly enriched in protein complex, postsynaptic density and cytoplasmic vesicle membrane (Fig. S4B). MF analysis for these genes consisted of protein binding, zinc ion binding and cadherin binding involved in cell-cell adhesion (Fig. S4C). KEGG pathway analysis was further utilized for candidate targeted genes of screened DE-miRNAs. As shown in Table 7, candidate targeted genes of screened DE-miRNAs were markedly enriched in Axon guidance. To better understanding the relationship among targeted genes of miRNAs, the PPI network was established by using the STRING database [13]. PPI node pairs with the score > 0.4 were chosen for further analysis. The top ten hub genes were verified based on the node number in the PPI network.

Verification Of Hub Gene Differential Expressions
The mRNA expressions of hub genes in BC were further verified by using UALCAN (http://ualcan.path.uab.edu), which is an interactive online database to perform in-depth analyses of gene expression data from TCGA [14].

Verification Of Hub Gene Prognostic Roles
The prognostic results of potential hub genes in BC were analyzed by using bc-GenExMiner v4.2 (bcgenex.centregauducheau.fr), a statistical mining tool of published BC transcriptomic data from TCGA and GEO [15].

Statistical Analysis
Majority of the statistical analysis was done through the above-mentioned bioinformatic tools, and only lncRNAs or miRNAs or genes with |log2FC|>2 and P < 0.05 were regarded as statistically significant.

Identification Of Upstream Transcription Factors Of De-mirnas
In current study, prediction of upstream TFs of screened DE-miRNAs was used by miRGen v3 database. Nine TFs for DE-miRNAs and corresponding motifs were presented in Fig. 2A-2I, which were NRF1, ELK4, E2F3, NR2F1, ZEB1, ZNF263, ZFX, POU2F2, and IRF1. As shown in Fig. 2J, NRF1 and ELK4 were the two most frequent predicted TFs of DE-miRNAs. The frequency of predicted TFs was also listed in Table 4.

Validation Of Hub Gene Expressions And Prognostic Roles
Using UALCAN database, we discovered that six of ten screened DE-miRNAs related hub genes were markedly upregulated in BC tissues than normal tissues. Three of ten screened DE-miRNAs related hub genes were significantly downregulated in BC tissues than normal tissues, whereas expression analysis of NTN1 showed no significant difference (Fig. 4A-4J).
To further identify potential hub genes, the prognostic functions of these hub genes in BC were conducted using bc-GenExMiner v4.2 database. As shown in Fig. 4K and 4L, a higher expression of RAB7A significantly indicated a worse prognosis while a higher expression of ARRB1 indicated a better prognosis of BC patients.
According to the predicted above-mentioned interactions, FKBP4 and FKBP5 related lncRNA-miRNA-mRNA regulatory axis related with development of BC were finally realized as presented in Fig. S6.

Discussion
It is widely acknowledged that there exist significant links between miRNA-mRNA regulatory axis and BC [16]. Recent studies have also suggested that lncRNAs could interact with other RNA transcripts via miRNA response element (MRE), which are proposed as the letters of a newfound RNA language [17]. For instance, lncRNA H19, transcriptional factor LIN28 as well as miRNA let-7 have been reported to form a doublenegative regulatory network, which palys a pivotal role during the maintenance breast cancer stem cells [18]. FKBPs have long been regarded as important regulators of the response to immunosuppressants FK506 and as molecular chaperones binding to different cellular receptors or targets [19]. More lately, various evidence has suggested that this complicated protein family might also play their roles in carcinogenesis, progression and chemoresistance of cancers [20][21][22].
Nevertheless, to our knowledge, a comprehensive FKBP-related ncRNA-mRNA regulatory axis in BC has not been established so far. In current study, we performed a differential expression analysis by using FKBPs mRNA data of GEPIA database. contributed to malignant proliferation of acute myeloid leukemia cells through upregulating ROCK1 [25]. Moreover, LINC00963 was found to facilitate osteosarcoma growth and progression via inhibiting miR-204-3p/FN1 axis [26].
Previous researches have suggested that the expression of miRNA could be modulated by TFs [27]. Therefore, we predicted nine TFs potentially regulating above-mentioned DE-miRNAs. Nuclear factor erythroid 2-like 1 (NRF1, including a short form Nrf1β/LCR-F1 and another long form TCF11) [28], was predicted as a TF potentially regulating expression of a relatively large proportion of screened DE-miRNAs. It has been demonstrated to act as an important player in regulating the expression and function of miRNAs. For example, a recent research has reported that NRF1 was participated in miR-219 signaling pathway, thereby inhibiting metastasis of BC cells [29]. Additionally, ETS-domain protein 4 (ELK4) was well elucidated to interact with miR-3188 in the development of atherosclerosis [30].
More researches on the functions of predicted TFs in BC are necessary to be further investigated.
Next, by integrating DE-mRNAs and targeted genes of DE-miRNAs, 320 candidate genes were identified. Subsequent GO and KEGG pathway analysis revealed that targeted genes were significantly enriched in cysteine-type endopeptidase activity involved in apoptotic process. A study performed by Siewiński et al. indicated that positive expression of high molecular weight cysteine proteinase inhibitor was observed on the tumor cell surface in serous and endometrioid metastatic ovarian cancer [31]. A plenty of investigations also suggested that apoptotic process correlated with BC [32][33][34], which further supported our current predicted findings.
Finally, PPI network was performed and top ten hub genes were verified. Moreover, differential expression analysis of these hub genes of BC were further conducted by using UALCAN database, including publicly available cancer OMICS data (TCGA and MET500).
Inspiringly, most of these genes have been demonstrated to act as key regulators of BC.
For instance, upregulated RAB7A was found correlated to poor prognosis of BC patients in this study, which is in accordance with the results of knockdown of RAB7A suppressing the proliferation and migration of BC cells [35]. In addition, analysis of hub genes' prognostic functions also implied significant tumor suppressive effect of ARRB1 in BC, which is in accordance with research results of Son et al [36]. Based on above-mentioned findings, we established a predicted FKBP-related ncRNA-mRNA regulatory network, which could be very important for probing novel mechanisms and possible therapeutic targets of BC.

Conclusion
In current study, we firstly established a predicted FKBP-related ncRNA-mRNA regulatory network, thus exploring a comprehensive interpretation of molecular mechanisms and providing potential clues in seeking novel therapeutics for BC. In the future, much more experiments should be conducted to verify our findings.

Supplementary Files
This is a list of supplementary files associated with the primary manuscript. Click to download.