Survey of miRNA-miRNA cooperative regulation principles across cancer types

Survey of miRNA-miRNA cooperative regulation principles across cancer types Abstract Cooperative regulation among multiple microRNAs (miRNAs) is a complex type of posttranscriptional regulation in human; however, the global view of the system-level regulatory principles across cancers is still unclear. Here, we investigated miRNA-miRNA cooperative regulatory landscape across 18 cancer types and summarized the regulatory principles of miRNAs. The miRNA-miRNA cooperative pan-cancer network exhibited a scale-free and modular architecture. Cancer types with similar tissue origins had high similarity in cooperative network structure and expression of cooperative miRNA pairs. In addition, cooperative miRNAs showed divergent properties, including higher expression, greater expression variation and a stronger regulatory strength towards targets and were likely to regulate cancer hallmark-related functions. We found a marked rewiring of miRNA-miRNA cooperation between various cancers and revealed conserved and rewired network miRNA hubs. We further identified the common hubs, cancer-specific hubs and other hubs, which tend to target known anticancer drug targets. Finally, miRNA cooperative modules were found to be associated with patient survival in several cancer types. Our study highlights the potential of pan-cancer miRNA-miRNA cooperative regulation as a novel paradigm that may aid in the discovery of tumorigenesis mechanisms and development of anticancer drugs. pan-cancer, cancer context miRNA–target regulation, miRNA-miRNA cooperative regulatory network, hub miRNAs, miRNA cooperative modules Introduction MicroRNAs (miRNAs) are small non-coding RNAs of approximately 22 nucleotides in length that play critical roles in various types of biological processes and complex diseases, including cancer [1–5]. A couple of studies have demonstrated that miRNAs perform their functions by forming miRNA-induced silencing complexes, thereby hindering mRNA translation or accelerating RNA degradation [6, 7]. With the development of high-throughput sequencing technology, an increasing number of miRNAs have been identified [8]. However, our knowledge of their functions in cancer is still limited. In addition, miRNA regulation exhibits great complexity in cancer. One miRNA can target multiple genes, and one gene can be regulated by more than one miRNA, indicating cooperative regulation among miRNAs [9]. Cooperative regulation among miRNAs has been identified in multiple species and various types of cancer. Wu et al. [10] first found that miRNAs from different miRNA clusters cooperatively regulate target gene expression, which signalled a new era of miRNA regulation. Bandi et al. [11] demonstrated that hsa-miR-34a and hsa-miR-15a/16 act synergistically to induce cell cycle arrest in non-small cell lung cancer. Kandi et al. [12] found that hsa-miR-125b and hsa-miR-99a encoded on chromosome 21 co-regulate vincristine resistance in childhood acute megakaryoblastic leukaemia. In addition, it has been demonstrated that hsa-miR-21 and hsa-miR-20a exhibited the largest synergistic effects in brain glioma cells [13]. Borzi et al. [14] investigated the interaction between hsa-miR-486-5p and hsa-miR-660-5p and assess their possible role and synergistic effect in lung cancer treatment. In another study, Lai et al. [15] validated that hsa-miR-205-5p and hsa-miR-342-3p achieve an enhanced repressive effect on E2F1 expression in a cooperative manner in an aggressive human non-small cell lung carcinoma cell line. At the same time, several miRNAs have also been identified as cancer biomarkers, such as hsa-miR-21-5p, hsa-miR-23a-3p and hsa-miR-27a-3p in pancreatic tumours and hsa-miR-142-3p and has-miR-494-3p in lymphoma [16, 17]. These observations indicate the presence of cooperative miRNA-mediated regulation, and studying the potential functional effects of miRNA cooperative regulation is an interesting and challenging area in miRNA research. However, it is a challenge to identify cooperative miRNA regulation using experimental methods because of the high number of miRNA combinations. Moreover, miRNA regulation is often reprogrammed in different tissues or different cancer types. To address this issue, several bioinformatics analysis pipelines have been proposed for inferring cooperative regulation among miRNAs, ranging from genomics to phenomics-based methods [18]. Several studies first constructed miRNA-miRNA interaction network using predicted miRNA–target data and then focused on hub miRNAs. Ying et al. [19] constructed a miRNA-miRNA co-regulating network in ovarian cancer, and hubs in the network were selected and identified as risk miRNAs, such as hsa-miR-579, hsa-miR-942, hsa-miR-105 and hsa-miR-150. Alshalalfa et al. [20] investigated mRNA-mediated miRNA-miRNA interactions in prostate cancer, and 11 miRNAs highly connected to other miRNAs were identified as diagnostic and prognostic biomarkers. In another study, Hua et al. [21] identified the miRNA-miRNA network for each breast cancer subtype and found that luminal A and basal-like subtype-specific networks exhibited changes in hub miRNAs. Yang et al. [22] constructed a miRNA-miRNA functional synergistic network and miRNA-miRNA co-regulatory network in breast cancer and proposed hub miRNAs (such as, hsa-miR-139-5p, hsa-miR-9, hsa-miR-375, hsa-miR-592 and hsa-miR-135a) that may be crucial in the initiation and progression of primary breast cancer. In addition, several studies first constructed miRNA-miRNA interaction networks using predicted miRNA–target data and then focused on modules. Zhang et al. [23] analysed a miRNA synergistic regulatory network in small cell lung cancer and identified a network module that included three miRNAs (hsa-let-7c, hsa-let-7b and hsa-let-7d). In addition, Li et al. [24] detected synergistic miRNA regulatory modules using ‘Mirsynery’ in ovarian, breast and thyroid cancer and identified several promising cancer-specific modules as prognostic biomarkers. Zhao et al. [25] constructed a miRNA-miRNA synergistic network in colorectal cancer and identified eight functional modules containing miRNAs that could implement specific functions as a whole in colorectal cancer (CRC). Especially, we proposed a functional module-based method for identification of cooperative miRNAs, and this method has been widely used in various types of cancer [9]. It is might expected that miRNA cooperative regulation is reprogrammed in different biological contexts. Context-specific miRNA-miRNA regulation, which means that the miRNA-miRNA regulation only occurred in a specific context, will provide a better approach for inferring the function of miRNAs in cancer. Thus, paired miRNA/mRNA expression levels across cancer samples were integrated to identify context-specific cooperative regulation among miRNAs [21, 24, 26]. Increasing cancer-specific miRNA-miRNA networks were analysed in various types of cancer, such as breast cancer, glioma and prostate cancer [20, 22, 27]. However, many questions regarding common or specific cooperative regulatory mechanisms in different cancer types have not been fully addressed, such as which collaborative miRNA regulation networks occur across all cancer types or in a specific cancer type? What about the functions, regulatory roles and biological insights of the pan-cancer miRNA-miRNA regulation, which include the whole common or specific cooperation cross all the cancer types? Moreover, increasing evidence has indicated that miRNA–gene regulation is critical for inferring cooperative miRNA regulation in cancer. Zhu et al. [28] applied the databases miRSel and ExprTargetDB jointly to obtain reliable miRNA–gene interaction data and determine miRNA synergy. Meng et al. [29] proposed that identification of context-dependent miRNA-target interactions is important, and they used common miRNA target-prediction programs, two databases and miRNA and mRNA expression data to obtain miRNA–target interaction information and identify functionally synergistic miRNA pairs. Song et al. [30] obtained miRNA target data from the Targetscan database to construct a miRNA-miRNA co-regulation network. Yang et al. [22] used the databases miRecords and miRWalk to obtain miRNA-gene interaction information and identified co-regulatory miRNA-miRNA networks. Integration of miRNA/mRNA gene expression to characterize the negative correlation between miRNA and target gene expression is widely used to identify context-specific regulation, thus contributing to the identification of cooperative miRNA regulation [31]. In addition to miRNAs, gene expression is also regulated by genetic and epigenetic regulators, such as DNA methylation [16]. Thus, to comprehensively understand miRNA–gene regulation as well as the cooperative regulation among miRNAs, we need to integrate multiple omics data sets across cancer types. The efforts of cancer genome projects, such as TCGA and ICGC, provide a great opportunity for us to investigate cooperative regulation among miRNAs across various cancer types [31, 32]. To address these issues, in this study, we focused on pan-cancer miRNA-miRNA cooperative regulation by assembling 18 cancer context cooperative regulation networks to investigate the biological insights they might present. The global principles of miRNA cooperative regulation were summarized. Based on the structure of the miRNA cooperative regulatory network and expression of cooperative miRNA pairs, we attempted to reveal the cancer tissue origins at different levels. In addition, we found that miRNAs in the pan-cancer miRNA-miRNA cooperative regulatory network show certain meaningful, expressional and functional properties, indicating that they might play important roles in tumorigenesis. Through comparative analyses of the structure of the pan-cancer and cancer-specific miRNA-miRNA cooperative regulatory networks, we identified the common hubs, cancer-specific hubs and other hubs and subsequently explored the application of hub miRNAs in clinical diagnosis and anticancer drug development. Finally, several miRNA modules were identified that might serve as novel prognostic biomarkers. These analyses and validations demonstrate how cancer-associated miRNA-miRNA cooperative regulatory networks can be used to accelerate discovery of miRNA-based biomarkers and potential therapeutics. Systematic investigation of miRNA-miRNA cooperative regulatory networks across multiple cancers can help us elucidate the commonalities and differences in the mechanisms driving specific cancer types. Results The cancer context-specific miRNA–target regulation To systematically characterize miRNA–target regulation across multiple cancer types, we used a three-step pipeline to identify cancer context-specific miRNA–target regulation by integration of genome-wide miRNA–target regulation and paired miRNA-mRNA expression profiles as well as copy number and DNA methylation data across 18 cancer types (Supplementary Figure S1A and Methods). Each cancer type is with paired mRNA and miRNA expression profiles as well as copy number and DNA methylation data measured for the same samples, forming a pan-cancer data compendium from 6348 tumour samples for 20 472 mRNAs and 867 miRNAs (Table 1). These cancer types were primarily involved in eight classes, namely, brain (brain lower-grade glioma, LGG; glioblastoma multiforme, GBM), genitourinary (kidney renal clear cell carcinoma, KIRC; kidney renal papillary cell carcinoma, KIRP; bladder urothelial carcinoma, BLCA; prostate adenocarcinoma, PRAD), thoracic (lung squamous cell carcinoma, LUSC; lung adenocarcinoma, LUAD), gastrointestinal (colon adenocarcinoma, COAD; stomach adenocarcinoma, STAD; liver hepatocellular carcinoma, LIHC), gynaecologic system (cervical squamous cell carcinoma and endocervical adenocarcinoma, CESC; ovarian serous cystadenoma carcinoma, OV; uterine corpus endometrial carcinoma, UCEC), breast carcinoma (BRCA), skin (skin cutaneous melanoma, SKCM) and head and neck cancer types (head and neck squamous cell carcinoma, HNSC; thyroid carcinoma, THCA). In addition to miRNA regulation, alterations in promoter DNA methylation and gene DNA copy number will likely alter the expression level of target genes and may introduce noise into the evaluation of actual miRNA regulation. Thus, the multivariate linear model can more accurately evaluate miRNA-mRNA expression associations in the presence of DNA copy number and promoter methylation aberrations that extensively influence mRNA expression under a given condition [33]. Under false discovery rate (FDR) < 0.05 for the significance of miRNA coefficients, candidate cancer context-specific miRNA-mRNA interactions were identified with a multivariate linear model (Supplementary Figure S1A and Methods). Considering miRNAs normally act as negative regulators, the expression profiles of genuinely regulating pairs were expected to be anti-correlated. As a result, we found that the majority of candidate cancer context miRNA-mRNA interactions are negatively correlated (Supplementary Figure S1B). Thus, we filtered the candidate miRNA-mRNA pairs and retained negatively correlated regulation pairs. On the other hand, several studies have recently reported that use of cross-linking and Argonaute (Ago) immunoprecipitation coupled with high-throughput sequencing (CLIP-Seq) could identify endogenous targets for miRNAs [34–37]. Indeed, we found that the cancer context miRNA-mRNA pairs have markedly more Ago CLIP-supported interactions than weakly associated pairs in all individual cancer types (Figure 1A). We further evaluated the precision of filtered cancer context miRNA-mRNA pairs with CLIP-supported sites based on low-throughput, experimentally validated miRNA–target interactions. We found that the precision of miRNA-mRNA pairs with CLIP-supported sites is average 7.5 times more than those without CLIP-supported sites (Figure 1B). For example, the proportion is increased to as high as 12.6 times in LIHC, and the minimum fold is to 4.7 in PRAD (Figure 1B). Thus, only candidate cancer context miRNA regulation with both negative correlations and CLIP-supported target sites were used for subsequent analyses (Supplementary Table S1), and all of these cancer context miRNA-mRNA interactions were assembled into the pan-cancer miRNA regulatory network. Table 1. Summary of analysed TCGA cancer data sets and cooperative miRNAs Cancer type Description Samples miRNAs mRNAs Cooperative miRNA pairs Cooperative miRNAs BLCA Bladder carcinoma 401 566 15 459 41 35 BRCA Breast carcinoma 603 506 15 455 51 45 CESC Cervical squamous cell carcinoma 292 517 15 425 3 6 COAD Colon adenocarcinoma 258 488 14 983 52 37 GBM Glioblastoma multiforme 78 474 15 090 50 50 HNSC Head and neck squamous cell carcinoma 474 520 15 447 24 25 KIRC Kidney renal clear cell carcinoma 179 429 15 180 56 45 KIRP Kidney renal papillary cell carcinoma 272 462 15 096 60 46 LGG Lower-grade glioma 509 518 15 147 222 99 LIHC Liver hepatocellular carcinoma 274 574 14 988 8 12 LUAD Lung adenocarcinoma 444 531 15 552 43 34 LUSC Lung squamous cell carcinoma 338 532 15 780 23 26 OV Ovarian cancer 293 715 15 742 53 51 PRAD Prostate adenocarcinoma 487 417 15 093 67 51 SKCM Skin cutaneous melanoma 351 597 15 239 10 15 STAD Stomach adenocarcinoma 367 499 15 820 73 46 THCA Thyroid carcinoma 498 485 14 920 37 35 UCEC Uterine corpus endometrial carcinoma 230 564 15 350 12 17 Total 6348 867a 20 472a 767 149 Cancer type Description Samples miRNAs mRNAs Cooperative miRNA pairs Cooperative miRNAs BLCA Bladder carcinoma 401 566 15 459 41 35 BRCA Breast carcinoma 603 506 15 455 51 45 CESC Cervical squamous cell carcinoma 292 517 15 425 3 6 COAD Colon adenocarcinoma 258 488 14 983 52 37 GBM Glioblastoma multiforme 78 474 15 090 50 50 HNSC Head and neck squamous cell carcinoma 474 520 15 447 24 25 KIRC Kidney renal clear cell carcinoma 179 429 15 180 56 45 KIRP Kidney renal papillary cell carcinoma 272 462 15 096 60 46 LGG Lower-grade glioma 509 518 15 147 222 99 LIHC Liver hepatocellular carcinoma 274 574 14 988 8 12 LUAD Lung adenocarcinoma 444 531 15 552 43 34 LUSC Lung squamous cell carcinoma 338 532 15 780 23 26 OV Ovarian cancer 293 715 15 742 53 51 PRAD Prostate adenocarcinoma 487 417 15 093 67 51 SKCM Skin cutaneous melanoma 351 597 15 239 10 15 STAD Stomach adenocarcinoma 367 499 15 820 73 46 THCA Thyroid carcinoma 498 485 14 920 37 35 UCEC Uterine corpus endometrial carcinoma 230 564 15 350 12 17 Total 6348 867a 20 472a 767 149 a miRNAs and mRNAs expressed in at least one cancer type. Table 1. Summary of analysed TCGA cancer data sets and cooperative miRNAs Cancer type Description Samples miRNAs mRNAs Cooperative miRNA pairs Cooperative miRNAs BLCA Bladder carcinoma 401 566 15 459 41 35 BRCA Breast carcinoma 603 506 15 455 51 45 CESC Cervical squamous cell carcinoma 292 517 15 425 3 6 COAD Colon adenocarcinoma 258 488 14 983 52 37 GBM Glioblastoma multiforme 78 474 15 090 50 50 HNSC Head and neck squamous cell carcinoma 474 520 15 447 24 25 KIRC Kidney renal clear cell carcinoma 179 429 15 180 56 45 KIRP Kidney renal papillary cell carcinoma 272 462 15 096 60 46 LGG Lower-grade glioma 509 518 15 147 222 99 LIHC Liver hepatocellular carcinoma 274 574 14 988 8 12 LUAD Lung adenocarcinoma 444 531 15 552 43 34 LUSC Lung squamous cell carcinoma 338 532 15 780 23 26 OV Ovarian cancer 293 715 15 742 53 51 PRAD Prostate adenocarcinoma 487 417 15 093 67 51 SKCM Skin cutaneous melanoma 351 597 15 239 10 15 STAD Stomach adenocarcinoma 367 499 15 820 73 46 THCA Thyroid carcinoma 498 485 14 920 37 35 UCEC Uterine corpus endometrial carcinoma 230 564 15 350 12 17 Total 6348 867a 20 472a 767 149 Cancer type Description Samples miRNAs mRNAs Cooperative miRNA pairs Cooperative miRNAs BLCA Bladder carcinoma 401 566 15 459 41 35 BRCA Breast carcinoma 603 506 15 455 51 45 CESC Cervical squamous cell carcinoma 292 517 15 425 3 6 COAD Colon adenocarcinoma 258 488 14 983 52 37 GBM Glioblastoma multiforme 78 474 15 090 50 50 HNSC Head and neck squamous cell carcinoma 474 520 15 447 24 25 KIRC Kidney renal clear cell carcinoma 179 429 15 180 56 45 KIRP Kidney renal papillary cell carcinoma 272 462 15 096 60 46 LGG Lower-grade glioma 509 518 15 147 222 99 LIHC Liver hepatocellular carcinoma 274 574 14 988 8 12 LUAD Lung adenocarcinoma 444 531 15 552 43 34 LUSC Lung squamous cell carcinoma 338 532 15 780 23 26 OV Ovarian cancer 293 715 15 742 53 51 PRAD Prostate adenocarcinoma 487 417 15 093 67 51 SKCM Skin cutaneous melanoma 351 597 15 239 10 15 STAD Stomach adenocarcinoma 367 499 15 820 73 46 THCA Thyroid carcinoma 498 485 14 920 37 35 UCEC Uterine corpus endometrial carcinoma 230 564 15 350 12 17 Total 6348 867a 20 472a 767 149 a miRNAs and mRNAs expressed in at least one cancer type. Figure 1. View largeDownload slide Cancer context miRNA–target regulation across cancer types. (A) Enrichment of Ago CLIP-supported miRNA–target interactions as a function of miRNA-mRNA expression association across cancer types. The miRNA–target interactions identified by a multivariate linear model with negative correlations were equally divided into 60 bins according to the P-values of regression. The y-axis represents the proportion of Ago CLIP-supported miRNA–target interactions occurring in each bin. The right x-axis represents the P-value of regression. (B) The left y-axis represents the number of experimentally validated miRNA–target interactions in the cancer context miRNA-mRNA interactions finally identified by us in each cancer, corresponding to the bars. The right y-axis represents the ratio between the precision of cancer context miRNA regulation with both negative correlations and CLIP-supported target sites and the precision of cancer context miRNA regulation only with negative correlations in each cancer, corresponding to the line. (C) The miRNA-mRNA interactions were divided into 11 categories according to the number of cancer types in which these interactions occur. The y-axis represents the proportion of each category of miRNA-mRNA regulation across diverse cancer types. (D) Proportion of miRNA-mRNA regulation occurring in different cancer types in different miRNA-mRNA regulation categories. (E) A pan-cancer miRNA-mRNA interaction subnetwork. A triangle node marks mRNA, and a circle node marks miRNA. An edge marks the miRNA-mRNA interaction. The miRNAs are shared by all the cancers, mRNAs occur in at least 16 cancers and miRNA-mRNA interactions occur in at least 15 cancers. Figure 1. View largeDownload slide Cancer context miRNA–target regulation across cancer types. (A) Enrichment of Ago CLIP-supported miRNA–target interactions as a function of miRNA-mRNA expression association across cancer types. The miRNA–target interactions identified by a multivariate linear model with negative correlations were equally divided into 60 bins according to the P-values of regression. The y-axis represents the proportion of Ago CLIP-supported miRNA–target interactions occurring in each bin. The right x-axis represents the P-value of regression. (B) The left y-axis represents the number of experimentally validated miRNA–target interactions in the cancer context miRNA-mRNA interactions finally identified by us in each cancer, corresponding to the bars. The right y-axis represents the ratio between the precision of cancer context miRNA regulation with both negative correlations and CLIP-supported target sites and the precision of cancer context miRNA regulation only with negative correlations in each cancer, corresponding to the line. (C) The miRNA-mRNA interactions were divided into 11 categories according to the number of cancer types in which these interactions occur. The y-axis represents the proportion of each category of miRNA-mRNA regulation across diverse cancer types. (D) Proportion of miRNA-mRNA regulation occurring in different cancer types in different miRNA-mRNA regulation categories. (E) A pan-cancer miRNA-mRNA interaction subnetwork. A triangle node marks mRNA, and a circle node marks miRNA. An edge marks the miRNA-mRNA interaction. The miRNAs are shared by all the cancers, mRNAs occur in at least 16 cancers and miRNA-mRNA interactions occur in at least 15 cancers. To have an overview of cancer context miRNA-mRNA pairs across 18 cancer types, we calculated the number of cancer types in which each miRNA-mRNA interaction occur. We found that the majority of cancer context miRNA-mRNA interactions occur in more than one cancer type (Figure 1C and Supplementary Figure S1C). Furthermore, approximately 24.9–71.0% of miRNA-mRNA interaction pairs were shared by at least five cancer types across each cancer context miRNA regulatory network (Figure 1C) and represented approximately 45.8% of pan-cancer miRNA-mRNA interactions. Conversely, the fraction of miRNA-mRNA interactions active in one cancer ranged from 2.7 to 30.5%, and most of them were present in GBM, OV, KIRC and UCEC, representing only approximately 12.5% of the pan-cancer miRNA-mRNA interaction network. These results suggest that most miRNA-mRNA interactions likely play roles in multiple cancer types, while a limited number of miRNA-mRNA interactions perform biological functions in only one specific cancer type. Moreover, miRNA-mRNA interactions were decreased, accompanied by a greater number of cancer types in which miRNA-mRNA interactions can occur. These results indicate that few core miRNA-mRNA interactions occur in most cancer types, but these miRNAs might cooperatively target genes across different cancers to perform a wide range of biological functions. We also calculated the contribution rate of each cancer in each type of miRNA regulation. As a result, we found that some cancers tend to have miRNA regulation in limited cancer types, such as GBM, LGG and OV, and there were also several cancers with the opposite tendency, for example, CESC, COAD and LUAD (Figure 1D), indicating that miRNA regulation has different important roles across cancers. Additionally, we found that most mRNAs are targeted by only a few miRNAs in cancer, as shown in Supplementary Figure S2. The number of mRNAs regulated by multiple miRNAs is small, but the proportion of these is relatively larger in LGG. Finally, we focused on miRNAs and mRNAs shared by multiple cancer types, and a pan-cancer miRNA-mRNA interaction subnetwork was constructed in which the miRNAs were shared by all the cancers, the mRNAs were present in at least 16 cancers and the miRNA-mRNA interactions occurred in at least 15 cancers (Figure 1E). MiRNA-miRNA cooperative interactions across cancer types reveal cancer clusters with similar tissue origin Next, we systematically identified miRNA-miRNA cooperation across 18 cancer types by using a multistage method [9]. Then, we assembled all identified cooperative miRNA pairs into a cooperative network in each cancer, where nodes represent miRNAs and edges represent their functional cooperative interactions (Table 1 and Supplementary Figure S3). Among the identified cooperative miRNA pairs, several have been reported in cancer, such as hsa-miR-105 in ovarian cancer, hsa-miR-205, hsa-miR-133a and hsa-miR-145 in prostate cancer, hsa-miR-9 in breast cancer [19, 20, 22]. Moreover, we constructed a miRNA-miRNA cooperative regulatory landscape in pan-cancer by integrating 18 cancer context cooperative networks (Figure 2). In total, the integrated network of pan-cancer contained 762 miRNA-miRNA cooperative interactions among 149 miRNAs (Figure 3A). Approximately 41.6% (62 of 149) of miRNAs exhibit cooperation in at least five cancers, and only 18.0% (27 of 149) of miRNAs exhibit cooperation in one cancer (Figure 3B upper panel), indicating that most miRNAs tend to cooperate with other miRNAs in several cancers. However, we discovered that for miRNA cooperation, 88.0% (674 of 762) are cancer-specific, and only 0.1% (1 of 762) of cooperative interactions are conserved in more than five cancers (Figure 3B, down panel). These results indicate that although miRNAs generally play cooperative roles in pan-cancer, they might cooperatively regulate with different miRNAs in a given cancer. The low conservation of miRNA cooperation may be explained in part by the cancer-specific expression of both miRNAs and genes. Especially, the cooperation between hsa-miR-200a-3p and hsa-miR-429 is shared by six cancers. Both of these miRNAs belong to miR-200 family, the role of which in tumour progression has been linked to several cancers [38, 39]. To obtain an overview of topological structure of this cooperative network in pan-cancer, we first examined degree distribution of this network, which revealed a power-law distribution. Therefore, similar to many large-scale networks, the cooperative network displayed scale-free characteristics (Figure 3C), indicating that the integrated cooperative network is not random but is characterized by a core set of organizing principles in its structure that distinguish it from randomly linked networks. The miRNA-miRNA synergistic network in ovarian cancer constructed by Ying et al. [19] was also a small-world network, with a few miRNAs interacting with a relatively large number of miRNA partners, and many miRNAs had few miRNA partners. Last, we found that the integrated cooperative network had higher clustering coefficients than randomly linked networks, as expected for module characteristics (Figure 3D). The scale-free and modular characteristics of the cooperative network in pan-cancer indicate that miRNA-miRNA cooperative interactions influence each other and effectively exchange regulated information both at a global and at a local scale. Figure 2. View largeDownload slide Cooperative miRNA-miRNA regulatory landscape in pan-cancer. Integrated Circos plot showing the landscape of miRNA-miRNA cooperative interactions in 18 cancer types. The outmost circle indicates the degree of cooperative miRNAs in the cooperative network in pan-cancer (the first track). The next outermost circle represents the number of cancer types in which the cooperative miRNAs occur (the second track). The next outermost circle represents the name of cooperative miRNAs (the third track). The 4th track to the 21st track represents 18 cancer types: if the miRNAs occur in this cancer type, a bar is presented that is labelled with the corresponding colour. In the centre of the figure, each arc indicates cooperative regulation between two miRNAs. The coloured arcs represent cooperative regulation occurring in a specific cancer type, whereas the black arcs represent regulation occurring in at least two cancer types. Figure 2. View largeDownload slide Cooperative miRNA-miRNA regulatory landscape in pan-cancer. Integrated Circos plot showing the landscape of miRNA-miRNA cooperative interactions in 18 cancer types. The outmost circle indicates the degree of cooperative miRNAs in the cooperative network in pan-cancer (the first track). The next outermost circle represents the number of cancer types in which the cooperative miRNAs occur (the second track). The next outermost circle represents the name of cooperative miRNAs (the third track). The 4th track to the 21st track represents 18 cancer types: if the miRNAs occur in this cancer type, a bar is presented that is labelled with the corresponding colour. In the centre of the figure, each arc indicates cooperative regulation between two miRNAs. The coloured arcs represent cooperative regulation occurring in a specific cancer type, whereas the black arcs represent regulation occurring in at least two cancer types. Figure 3. View largeDownload slide Integrated view of cooperative miRNA-miRNA regulation in pan-cancer. (A) The pan-caner miRNA-miRNA cooperative interaction network. This network consists of 762 cooperative interactions between 149 miRNAs. A circle node represents miRNA. The node size is proportional to the degree of miRNAs. The number of cancer types in which the miRNA occurs is mapped to the green gradient bar. An edge represents cooperative regulation between miRNAs. The width of the edge represents the number of cancer types in which the miRNA-miRNA pair occurs. (B) Distribution of cooperative miRNA and miRNA pairs. (Upper panel) The pie chart shows the proportion of cooperative miRNAs present in different numbers of cancer types. (Down panel) The pie chart shows the proportion of cooperative miRNA pairs present in different numbers of cancer types. (C) Distribution of degree for the observed pan-cancer miRNA-miRNA cooperative network (red circles and red line) and permuted networks (box plots). (D) The cluster coefficient of the pan-cancer miRNA-miRNA cooperative network is higher than those of the random networks (E) and (F). The Jaccard similarity coefficient matrix shows the similarity between each pair of miRNA-miRNA cooperative networks from the structure in (E). The Pearson's correlation of the cooperative miRNA pair matrix shows the similarity between each pair of miRNA-miRNA cooperative network from the expression in (F). Cancers with the same origin are labelled with grey shading and dashed boxes. (G) and (H) The subnetwork shared by COAD and STAD. A circle node represents miRNA. The black solid line presents miRNA-miRNA cooperative interactions shared by two cancers. The red line and light-blue line, respectively, represent COAD- and STAD-specific miRNA-miRNA cooperative interactions. Figure 3. View largeDownload slide Integrated view of cooperative miRNA-miRNA regulation in pan-cancer. (A) The pan-caner miRNA-miRNA cooperative interaction network. This network consists of 762 cooperative interactions between 149 miRNAs. A circle node represents miRNA. The node size is proportional to the degree of miRNAs. The number of cancer types in which the miRNA occurs is mapped to the green gradient bar. An edge represents cooperative regulation between miRNAs. The width of the edge represents the number of cancer types in which the miRNA-miRNA pair occurs. (B) Distribution of cooperative miRNA and miRNA pairs. (Upper panel) The pie chart shows the proportion of cooperative miRNAs present in different numbers of cancer types. (Down panel) The pie chart shows the proportion of cooperative miRNA pairs present in different numbers of cancer types. (C) Distribution of degree for the observed pan-cancer miRNA-miRNA cooperative network (red circles and red line) and permuted networks (box plots). (D) The cluster coefficient of the pan-cancer miRNA-miRNA cooperative network is higher than those of the random networks (E) and (F). The Jaccard similarity coefficient matrix shows the similarity between each pair of miRNA-miRNA cooperative networks from the structure in (E). The Pearson's correlation of the cooperative miRNA pair matrix shows the similarity between each pair of miRNA-miRNA cooperative network from the expression in (F). Cancers with the same origin are labelled with grey shading and dashed boxes. (G) and (H) The subnetwork shared by COAD and STAD. A circle node represents miRNA. The black solid line presents miRNA-miRNA cooperative interactions shared by two cancers. The red line and light-blue line, respectively, represent COAD- and STAD-specific miRNA-miRNA cooperative interactions. Several lines of evidence have indicated that cancer types with similar tissue origins share multiple molecular features, such as gene expression, at the level of both protein-coding genes and miRNAs [40–42]. However, whether cancer types with similar tissue origins exhibit similar miRNA cooperation patterns is still unknown. To address this question, we calculated a paired similarity score (Jaccard coefficient) based on miRNA cooperation in each cancer. We found that adenocarcinomas, such as COAD, STAD, LUAD and PRAD, showed greater similarity to each other (Figure 3E). This analysis indicated that cancer types with similar tissue origins show similar miRNA cooperation patterns. Especially, COAD and STAD shared 24 miRNAs, and eight of these cooperation patterns were common for the two cancers (Figure 3G and H). We further discovered that these miRNA cooperation patterns are closely related to carcinogenesis by regulating kappa B kinase NF-kappa B cascade, cytoskeleton organization and biogenesis. Previous studies have shown that hsa-miR-106a-5p, hsa-miR-30c-5p and hsa-miR-142-3p are involved in both COAD and STAD, and here, we found they cooperatively regulate cytoplasm organization and biogenesis processes and intercellular junction assembly and maintenance processes. These observations suggest that related miRNA-miRNA cooperative regulation mechanisms might operate in cancer types with similar tissue origins. Moreover, we also measured the similarity between diverse cancers based on the expression correlation coefficient of cooperative miRNAs and found that the expression pattern of pan-cancer cooperative miRNAs also revealed tissue-of-origin cancer. After clustering all cancers, we found that three squamous cell carcinomas, HNSC, LUSC and CESC, are closely clustered together (Figure 3F and Supplementary Figure S4). Another example is the cluster comprising KIRC and KICH, which are two types of kidney cancer. We found that their similarity in cooperative miRNA expression was higher than with other cancers. Thus, we concluded that both the structure of miRNA cooperative networks and the expression of cooperative miRNA pairs could reveal similar tissue origin of known and novel cancer clusters at different levels. Common regulatory principles of cooperative miRNA-miRNA pairs Through systematic analysis of the expression and function of cooperative regulatory miRNA-miRNA networks across cancers, some common regulatory principles of the cooperative regulatory miRNA-miRNA networks were summarized. First, we compared the expression levels of cooperative miRNAs with those of other miRNAs in each cancer, and the results show that the expression of cooperative miRNAs is higher than that of other miRNAs (Figure 4A). Moreover, the expression of cooperative miRNAs in multiple cancers tended to be higher than that of those involved in few cancers (Supplementary Figure S5). At the same time, the expression fluctuation of these cooperative miRNAs is significantly greater than that of non-cooperative miRNAs in each cancer (Figure 4B). Especially, cooperative miRNAs have the most expression fluctuation in LGG, and we proposed that these cooperative miRNAs might be important miRNAs in LGG [43]. We found that cooperative miRNA pairs exhibit stronger expression correlation than random miRNA pairs in each cancer (Figure 4C). Moreover, we also calculated the average regulatory strength of cooperative miRNAs towards their target genes and found significantly stronger regulation than non-cooperative miRNAs in 14 of the cancer types (Figure 4D). These results suggest that similar expression patterns might help miRNAs to perform cooperative functions. Figure 4. View largeDownload slide Common regulatory principles of cooperative regulatory miRNA-miRNA pairs. (A) The expression of cooperative miRNAs across diverse cancer types. The grey box represents the expression of miRNAs randomly selected. (B) Radar plot representing the expression fluctuation of cooperative miRNAs. The axes report the ratio of the number of cooperative miRNAs with smaller (or larger) variance to the number of total miRNAs with smaller (or larger) variance. The blue line represents cooperative miRNAs with larger variance. The red line represents cooperative miRNAs with smaller variance. (C) The correlation of expression of cooperative miRNA pairs across diverse cancer types. The grey box represents the expression correlation between non-cooperative miRNA pairs. The Pearson’s correlation coefficient was used. (D) The regulatory strength of cooperative miRNAs towards target genes across diverse cancer types. The grey box represents non-cooperative miRNA. The y-axis represents the mean of the regression coefficient (absolute value) of miRNAs in the multivariate linear model. (E) The cancer hallmark-related miRNA-miRNA cooperative subnetwork, which was assembled with the miRNA pairs that co-regulated at least one cancer hallmark-related GO term. (F) A heat-map containing the number of cooperative miRNAs that regulate GO terms related to cancer hallmarks across diverse cancer types. Figure 4. View largeDownload slide Common regulatory principles of cooperative regulatory miRNA-miRNA pairs. (A) The expression of cooperative miRNAs across diverse cancer types. The grey box represents the expression of miRNAs randomly selected. (B) Radar plot representing the expression fluctuation of cooperative miRNAs. The axes report the ratio of the number of cooperative miRNAs with smaller (or larger) variance to the number of total miRNAs with smaller (or larger) variance. The blue line represents cooperative miRNAs with larger variance. The red line represents cooperative miRNAs with smaller variance. (C) The correlation of expression of cooperative miRNA pairs across diverse cancer types. The grey box represents the expression correlation between non-cooperative miRNA pairs. The Pearson’s correlation coefficient was used. (D) The regulatory strength of cooperative miRNAs towards target genes across diverse cancer types. The grey box represents non-cooperative miRNA. The y-axis represents the mean of the regression coefficient (absolute value) of miRNAs in the multivariate linear model. (E) The cancer hallmark-related miRNA-miRNA cooperative subnetwork, which was assembled with the miRNA pairs that co-regulated at least one cancer hallmark-related GO term. (F) A heat-map containing the number of cooperative miRNAs that regulate GO terms related to cancer hallmarks across diverse cancer types. The miRNA gene family is an important component of small RNAs, and we further investigated the cooperation between homologous miRNAs. We found that the proportion of miRNAs that could cooperate with their family members is different cross-cancer types (Supplementary Figure S6). The proportion of miRNAs that could cooperate with their family members is up to 60% in BLCA. For example, hsa-miR-200a-3p and hsa-miR-429 belong to the mir-8 family and cooperatively regulate the insulin receptor signalling pathway. In addition, it is well known that cancer miRNAs can be onco-miRNAs or suppressive miRNAs. We further explored cooperative interactions among oncogenic miRNAs and suppressive miRNAs [44]. We found that miRNAs belonging to the same category tend to cooperatively interact in several cancer types, such as BRCA, COAD, GBM and PRAD (Supplementary Figure S7). For example, hsa-miR-21-5p and hsa-miR-27a-3p, which are both oncogenic miRNAs, cooperatively interact in BRCA. The complexity of cancer can be reduced and represented by a few cancer hallmarks that enable tumour growth and metastasis [45]. These hallmarks provide a framework for understanding the remarkable diversity of cancers. To further investigate the functional roles of cooperative miRNA pairs in cancer hallmarks, we selected the miRNA pairs that co-regulated at least one cancer hallmark-related Gene Ontology (GO) term. These miRNA pairs were assembled into a cancer hallmark-related subnetwork (Figure 4E and Supplementary Table S2). We found that several cancer hallmarks are deeply regulated by the cooperative miRNA subnetwork. For example, 22 cooperative miRNAs in STAD are involved in regulating cancer hallmarks, and eight of these regulated tissue invasion and metastasis, while in LGG, 12 of 16 miRNAs regulated insensitivity to antigrowth signals (Figure 4F). Moreover, miRNA cooperation regulated the same cancer hallmarks in different cancers; for instance, 26 cooperative miRNAs negatively regulated both apoptosis and programmed cell death, evading apoptosis in eight types of cancer. In addition to cancer hallmarks reflecting the common characteristics of cancer, cooperative miRNAs can also control specific signalling pathways (Supplementary Figure S8). Two pairs of miRNAs, hsa-miR-301a-3p and hsa-miR-130a-3p pair, and hsa-miR-200b-3p and hsa-miR-429 pair, regulate the epidermal growth factor receptor signalling pathway in OV. These results indicate cancer specificities for cooperative miRNA-mediated regulation of cancer hallmarks and signalling pathways. Moreover, same miRNA pairs might regulate different hallmarks in different cancers. For example, hsa-miR-106a-5p and hsa-miR-340-5p cooperatively regulate cytoplasm organization and biogenesis in STAD but positively regulate transcription signal transduction in KIRP. Thus, we conclude that cooperative miRNA pairs may contribute to carcinogenesis by regulating cancer hallmark functions in multiple cancers as well as in cancer-specific patterns. Conserved and rewired cooperative network miRNA hubs Hub nodes whose connectivity is extremely high always play important roles in biological network. For example, in the protein–protein interaction networks of various organisms, hubs tend to be essential proteins [46]. To systematically investigate the key nodes in the miRNA-miRNA cooperative network, we first identified the top 30.0% of nodes with the highest connectivity within the pan-cancer cooperative network as hub miRNAs. Next, we assessed the extent to which these hub miRNAs are shared among the different cancer-specific miRNA-miRNA cooperative networks. These hubs were grouped into three categories: common hubs, cancer-specific hubs and other hubs. The common hubs exhibit central roles of miRNAs across more than one cancer; the specific hubs identify miRNAs with specific central roles in a given cancer. For example, hsa-miR-301a-3p is involved in 14 cancer types, is a hub in the pan-cancer network and is also a hub in HNSC, OV and KIRP. Thus, hsa-miR-301a-3p was defined as a common hub based on our method. hsa-miR-186-5p is also a common hub that participates in 17 cancer cooperative networks and is a hub miRNA in 2 (COAD, LUAD) cancer cooperative networks (Figure 5A). Common hub hsa-miR-186-5p has 29 cooperative miRNA partners in the pan-cancer cooperative network. In addition, 36.7% (18 of 49) of hubs are specific hubs, and the remaining 14.3% are other hubs (Figure 5A). For example, hsa-miR-107 is a specific hub in OV. To further investigate the roles of these three types of hubs, we associated each miRNA hub with biological function based on the function modules they cooperatively regulate. The results showed that the three categories of hubs regulated many common functions, such as regulating apoptosis, regulation of programmed cell death and some signalling pathways (Figure 5B and Supplementary Table S3). At the same time, TAM was also used to perform functional annotation for these hub miRNAs directly [47]. The annotation results showed that hub miRNA regulate many common cancer-related functions, such as cell proliferation, DNA repair, apoptosis and immune response. These results indicated that these miRNAs might play important roles in cancer. On the other hand, some functions are specifically controlled by different categories of hubs, such as angiogenesis and vasculature development, which are specifically regulated by the hubs hsa-miR-193b-3p and hsa-miR-29b-3p (Figure 5B). Thus, cancer-specific miRNA hubs may contribute to carcinogenesis by regulating cancer-specific functions. Figure 5. View largeDownload slide The conserved and rewired network hubs in each cancer type. (A) A summary bubble-bar plot showing the distribution of hub miRNAs in each cancer-specific miRNA-miRNA cooperative network. The top bars show the number of hub miRNAs (blue bars) and experimentally validated miRNAs related to cancer (grey bars) occurring in each cancer miRNA-miRNA cooperative network. In addition, the bars on the right show the degree of hub miRNAs in the pan-cancer miRNA-miRNA cooperative network. The bubble size indicates the degree of hub miRNAs standardized by the total number of miRNAs in each cancer miRNA-miRNA cooperative network, and different colours correspond to different hub categories. (B) Selected GO functional categories regulated by common hubs, cancer-specific hubs and other hubs. The bar on the right shows GO functional categories regulated by at least two types of hubs, and the bar on the left shows GO functional categories regulated by only one type of hub. The length of the bars represents the number of hub miRNAs. (C) Comparison of the number of anticancer drug targets in three types of hub miRNA targets. (D) An example of a common hub miRNA, hsa-miR-186-5p. (E) Clustering trees based on the expression of hub miRNAs that are differently expressed between cancer and normal tissues in BLCA, LUAD and UCEC. Dark cyan represents the normal samples, and red represents the tumour samples. (F) The ROCs of the hub miRNAs differentially expressed between cancer and normal tissues in BLCA, LUAD and UCEC. Figure 5. View largeDownload slide The conserved and rewired network hubs in each cancer type. (A) A summary bubble-bar plot showing the distribution of hub miRNAs in each cancer-specific miRNA-miRNA cooperative network. The top bars show the number of hub miRNAs (blue bars) and experimentally validated miRNAs related to cancer (grey bars) occurring in each cancer miRNA-miRNA cooperative network. In addition, the bars on the right show the degree of hub miRNAs in the pan-cancer miRNA-miRNA cooperative network. The bubble size indicates the degree of hub miRNAs standardized by the total number of miRNAs in each cancer miRNA-miRNA cooperative network, and different colours correspond to different hub categories. (B) Selected GO functional categories regulated by common hubs, cancer-specific hubs and other hubs. The bar on the right shows GO functional categories regulated by at least two types of hubs, and the bar on the left shows GO functional categories regulated by only one type of hub. The length of the bars represents the number of hub miRNAs. (C) Comparison of the number of anticancer drug targets in three types of hub miRNA targets. (D) An example of a common hub miRNA, hsa-miR-186-5p. (E) Clustering trees based on the expression of hub miRNAs that are differently expressed between cancer and normal tissues in BLCA, LUAD and UCEC. Dark cyan represents the normal samples, and red represents the tumour samples. (F) The ROCs of the hub miRNAs differentially expressed between cancer and normal tissues in BLCA, LUAD and UCEC. We further systematically investigated the association between miRNA cooperation and anticancer drugs. We found that 72 of 105 FDA-approved anticancer drugs could target hub miRNA–target genes. The common hub miRNAs contained the largest number of anticancer drug target genes compared with specific hubs and other hubs (Figure 5C). The results suggest that targets of common hub miRNAs are more likely to be drug-able. For example, some target genes of the common miRNA hsa-miR-186-5p, including FGFR1 and GART, have been used as anticancer drugs targets of Sorafenib in KIRC/KIRP and Pemetrexed in LUAD/LUSC (Figure 5D). According to the pan-cancer characteristics of hsa-miR-186-5p, target genes of this miRNA are favoured to be drugable and may be potential targets of anticancer drugs in the other six cancer types. In addition, we also found that 31 of 49 hub miRNAs could regulate cancer hallmarks. The common hub miRNAs regulate more hallmarks than the specific hub and other hub miRNAs (Supplementary Table S4 and Supplementary Figure S9A). For example, hsa-miR-106a-5p regulates self-sufficiency in growth signals and tissue invasion and metastasis. We found that 49.0% (24 of 49) of hub miRNAs are experimentally validated cancer-related miRNAs, including 12 common hubs, 9 specific hubs and 3 other hubs (Supplementary Table S4 and Supplementary Figure S9B). For example, the common hub hsa-miR-137 has been be associated with multiple cancer types, including GBM, BLCA, COAD, STAD, LIHC, OV, BRCA and SKCM. Finally, we further investigated the expression of hub miRNAs. We also found that the average expression level of other hub miRNAs is higher than that of specific hubs and common hubs (Supplementary Table S4 and Supplementary Figure S9C). These results suggested that hub miRNAs might play fundamental roles in multiple cancer types. It is widely accepted that the occurrence of cancer is related to the concerted activity of many miRNAs [48]. Functionally coherent miRNA cooperative modules can be used as biomarkers to distinguish cancer from normal tissues. Here, we evaluated the efficiency of the hubs using the area under the receiver operating characteristic (ROC) curve (AUC). ROC curves of the hubs exhibit an AUC value >0.9 in BLCA, LUAD and UCEC (Figure 5E and F). MiRNA cooperative modules can be used as prognostic biomarkers in multiple cancer types The notion that it may be possible to reduce cancer mortality by identifying and monitoring survival-related biomarkers rests on the idea that a module biomarker is a better predictor of survival than an individual gene [49]. Given that integrated cooperative network is associated with cancer-related biological processes and has the characteristics of pan-cancer, as shown in the above analysis results, we were interested in identifying prognostic modules from the pan-cancer cooperative network. Modules in the miRNA cooperative network represent groups of functionally related miRNAs dedicated to specific biological processes. We therefore used CFinder software with suggested parameters to identify modules in the pan-cancer cooperative network [50]. In total, with K = 6, two miRNA cooperative modules were identified (Figure 6A and Supplementary Figure S10). One module comprised seven common hubs, and the other module included five common hubs and two specific hubs. For instance, Figure 6A shows a module that includes 20 cooperative pairs among seven common hub miRNAs, each of which are involved in nine types of cancer miRNA cooperative networks on average. Functional analysis of the module revealed that these miRNAs are involved in multiple cancer-associated hallmarks, such as evading apoptosis, self-sufficiency in growth signals and insensitivity to antigrowth signals (Figure 6A). We further found that all the miRNAs in this module are well-known cancer-related miRNAs. The results demonstrate that this module is closely related to cancer development and represents pan-cancer characteristics. Figure 6. View largeDownload slide The miRNA-miRNA cooperative module as a biomarker. (A) An example of a miRNA-miRNA cooperative module. The Circos plot shows the characteristic of the miRNAs in the module. In the centre of the figure, the miRNA-miRNA cooperative module is shown. In this module, a circle node represents miRNA; the node size is proportional to the degree of miRNAs; the number of cancer types in which the miRNA occurs is mapped to the green gradient bar (details in Figure 3A); an edge represents the cooperative regulation between miRNAs. The circle is equally divided into seven parts by the seven miRNAs. Each sector corresponds to one miRNA. The innermost circle indicates cancer hallmarks regulated by the miRNAs (the first track). The next innermost circle represents the association between miRNAs and cancer. If it has been indicated that the miRNA is related to the cancer, a bar is presented that is labelled with the corresponding colour (the second track). The next circle represents the association between miRNAs and cancer-specific cooperative miRNA-miRNA networks. If the miRNA occurred in the cooperative miRNA-miRNA network of this cancer type, a bar is presented that is labelled with the corresponding colour (the third track). The next circle represents the degree of miRNAs in the corresponding cancer-specific cooperative miRNA-miRNA networks, which is presented as a bar plot (the fourth track). (B) The expression pattern of miRNAs from the module that were significant in single-variable Cox regression and their performance in survival analysis based on five types of cancer patients in training set and test set. (Upper panel) Colour-gram of miRNA expression profiles of patients ordered by risk score (details in the ‘Materials and methods’ section). The median of the risk score cut-off was used to divide patients into low-risk (grey line) and high-risk groups (brick-red line). (Down panel) Kaplan–Meier estimates of overall survival of the two patient subgroups identified by the miRNA signature in training set and test set. Figure 6. View largeDownload slide The miRNA-miRNA cooperative module as a biomarker. (A) An example of a miRNA-miRNA cooperative module. The Circos plot shows the characteristic of the miRNAs in the module. In the centre of the figure, the miRNA-miRNA cooperative module is shown. In this module, a circle node represents miRNA; the node size is proportional to the degree of miRNAs; the number of cancer types in which the miRNA occurs is mapped to the green gradient bar (details in Figure 3A); an edge represents the cooperative regulation between miRNAs. The circle is equally divided into seven parts by the seven miRNAs. Each sector corresponds to one miRNA. The innermost circle indicates cancer hallmarks regulated by the miRNAs (the first track). The next innermost circle represents the association between miRNAs and cancer. If it has been indicated that the miRNA is related to the cancer, a bar is presented that is labelled with the corresponding colour (the second track). The next circle represents the association between miRNAs and cancer-specific cooperative miRNA-miRNA networks. If the miRNA occurred in the cooperative miRNA-miRNA network of this cancer type, a bar is presented that is labelled with the corresponding colour (the third track). The next circle represents the degree of miRNAs in the corresponding cancer-specific cooperative miRNA-miRNA networks, which is presented as a bar plot (the fourth track). (B) The expression pattern of miRNAs from the module that were significant in single-variable Cox regression and their performance in survival analysis based on five types of cancer patients in training set and test set. (Upper panel) Colour-gram of miRNA expression profiles of patients ordered by risk score (details in the ‘Materials and methods’ section). The median of the risk score cut-off was used to divide patients into low-risk (grey line) and high-risk groups (brick-red line). (Down panel) Kaplan–Meier estimates of overall survival of the two patient subgroups identified by the miRNA signature in training set and test set. Subsequently, we examined whether the modules are associated with the survival of cancer patients. In each cancer, we randomly assigned 70% of the specimens to a training set and 30% of the specimens to a test set; we only used the training data set for detecting the miRNA signature of cancer. We found that miRNAs in this module can be used to classify cancer samples into two groups with significantly different overall survival rates in five cancers (log-rank test, P < 0.05) (Figure 6B). In addition, we also validated the miRNA signature in the test set where samples were classified into the high-risk group or low-risk group using the same miRNAs used in the training set. As an example, four miRNAs with six cooperative interactions in this module were identified as being associated with LGG patient survival. All the three miRNAs, hsa-miR-142-3p, hsa-miR-200a-3p and hsa-miR-30e-5p, are LGG-risk miRNAs for which increased expression is associated with an increased risk of survival. Using the median risk value scores of these miRNAs as the threshold, the LGG patients were divided into a high-risk group and a low-risk group, and these two groups exhibited significantly different survival times in training set and test set (Figure 6B, P = 0.002 and P = 0.009). Three of these miRNAs, hsa-miR-142-3p, hsa-miR-30e-5p and hsa-miR-200a-3p, occur in the LGG miRNA-miRNA cooperation network, the degree of which, respectively, is 9, 5 and 5. The three miRNAs are all known cancer-associated miRNAs. Thus, this suggests that a cooperative network-based approach could effectively identify prognostic biomarkers. Discussion In the current study, we identified miRNA-miRNA cooperation in each cancer via the multiple-step method. Systematic construction and analysis of miRNA-miRNA cooperation networks across multiple cancers can help to elucidate the commonalities and differences in cancer mechanisms. As general biological networks, the pan-cancer miRNA-miRNA cooperation network is scale-free and has a small-world property. Watts and Strogatz [51] have analysed how fast disturbances spread through small-world networks and revealed that the time wasted for spreading of a disturbance in a small-world network is close to the theoretically possible minimum for any graph with the same number of nodes and edges. Therefore, small-world may allow the cooperation of miRNAs to respond quickly to disturbances. The increased number of molecular tumour maps presents a great opportunity for us to investigate molecular characteristic across various cancer types. Different molecular characteristics have shown different conservation levels and specificity among cancers. In our study, through comparison of miRNA cooperation networks across cancers, we showed that only a small proportion of the cooperative activities are shared by multiple cancers. miRNA cooperation varied greatly from one cancer type to another. Park et al. [52] found that >90% of the interactions between cancer driver alterations were only detected in a single cancer type. Additionally, Xu et al. [53] found that ∼34.56% of competitive endogenous RNA (ceRNA) regulation occurred in only one cancer and only 3.78% of ceRNA-ceRNA interactions were conserved in >10 cancers. Together, these studies illustrate how the interactions between various molecules change across different types of cancer. However, some conservation of molecules among different tumour types has been found. Hamilton et al. [54] found that 45% of broadly conserved pan-cancer oncomiRs share strong homology in their seed motifs. Davoli et al. [55] identified distinct sets of drivers in each tumour type, but the majority of these were also identified in the pan-cancer analysis as lower confidence candidates, and there is still significant overlap among different tumour types. This suggested that molecules may be conserved, and these molecule interactions might be different in different types of cancer. In 2017, the US Food and Drug Administration has approved for the first time a cancer treatment based on common molecular biomarkers rather than the location in the body where the tumour originated. This is important news for the cancer community. Molecular markers may be able to more accurately reflect the curative effect of drugs than traditional tissue origins. Heim et al. [56] found that tumours at a particular anatomic site are genetically more similar to tumours from different organs and tissues than to tumours of the same origin. Snyder and colleagues [57] discovered a neoantigen signature that was specifically present in tumours with a sustained response to CTLA-4 blockade and that examining the exomes of melanoma patients could inform the response to the anti–CTLA-4 blockade therapy. However, Hoadley et al. performed an integrative analysis using molecular profiles of 12 cancer types, revealing a unified classification into 11 major subtypes. Five subtypes were nearly identical to their tissue-of-origin counterparts, but several distinct cancer types were found to converge into common subtypes [40]. In our study, we found that some cancer types with similar tissue origins exhibited higher cooperative network structural similarity and expression of cooperative miRNA pairs and that some cancer types with the same tissue origin are different in miRNA cooperative regulation. Cancer types with similar tissue origins do not necessarily have the same molecular characteristics. Our findings provide independent evidence for development of drug and treatment strategies. Moreover, most cooperative miRNA pairs also exhibited the same expression tendency, which allows for a rapid response to disturbances, and we discovered that cooperative miRNAs have more strong regulation to the target genes in function modules. Because miRNA regulation mostly reduces target mRNA expression, miRNA–target genes with low mRNA expression levels would more likely be under strong regulation by miRNAs [58]. In our multivariate linear model, when the absolute value of the regression coefficient of miRNA was larger, the gene expression targeted by this miRNA was lower. Therefore, when the absolute value of the regression coefficient of miRNA is larger, we think that the miRNA exhibits a stronger regulation effect on the target genes. We discovered that cooperative miRNAs exhibit stronger regulation of the target genes in function modules by comparing regression coefficients. In addition, Mendell [59] has proposed that the functions of a given miRNA could be attributed to strong regulation of a select few dominant targets or, alternatively, more subtle regulation of many targets simultaneously. The target genes in function modules may be as dominant targets. Our analysis also provides a hub-based view to elucidate the common and specific miRNA hubs across cancers. Hubs are topologically centred in the miRNA cooperative network and have maximal informational links with other miRNAs. Through in-depth analyses of the structure of the pan-cancer miRNA cooperative network, we identified the cancer-common and specific hubs. Hub miRNAs can be used as potential biomarkers to distinguish cancer from normal tissues. It has been shown that serum miRNA levels are stable, reproducible and consistent in humans [60]. Meanwhile, miRNAs can be released from cancer cells to body fluids through secretion of exosomes. Owing to the surprising stability of miRNAs in tissues and body fluids, miRNAs have emerged as advantageous non-invasive cancer biomarkers with immeasurable clinical potential. For example, in early gastric cancer, almost all of the serum-based traditional biomarkers, including carcinoembryonic antigen, are not widely recognized in the clinical screening or diagnosis because of their limited specificity and sensitivity [61]. The stability of miRNAs in serum and gastric juice has been confirmed, and the expression levels of certain miRNA diagnostic biomarkers can assist in screening for gastric cancer [62]. Several miRNAs circulating in the blood of gastric cancer patients can be applied as diagnostic biomarkers, including hsa-miR-21, hsa-miR-106a/b, hsa-miR-421, hsa-miR-129, hsa-miR-199a-3p, hsa-miR-218 and hsa-miR-223 [63–67]. Some common hub miRNAs have been identified as diagnostic biomarkers in several cancers. hsa-miR-200a is a common hub miRNA, and Yun et al. [68] found that hsa-miR-200a is an independent predictor of non-muscle-invasive bladder cancer recurrence. Park et al. [69] found that the detection of hsa-miR-125a and hsa-miR-200a in saliva can be used as non-invasive and rapid diagnostic biomarkers of oral cancer. The hub miRNA hsa-miR-34a has been identified as a novel biomarker of relapse in surgically resected non-small cell lung cancer [70]. Additionally, both of the common hub miRNAs hsa-miR-128 and hsa-miR-106a-5p have been used as diagnostic markers in GBM and lung squamous cell carcinoma. Wang et al. [71] demonstrated that cell-free hsa-miR-21, hsa-miR-128 and hsa-miR-342-3p in plasma are specific and sensitive markers for GBM diagnosis, suggesting that these miRNAs may be used as non-invasive biomarkers of GBM. Zhang et al. [72] identified a three-miRNA panel comprising hsa-miR-106a-5p, hsa-miR-20a-5p and hsa-miR-93-5p as a potential diagnostic marker for male lung squamous cell carcinoma (SCC) patients. In summary, we presented the miRNA-miRNA cooperation landscape across human major cancers and showed the importance of various aspects. A bird’s-eye view of the functional miRNA cooperative networks of large sample sets encompassing multiple tumour lineages are important for understanding the key roles of miRNAs in disease. Our study opens new avenues for leveraging publicly available genomic data to investigate the functions and mechanisms of miRNAs across human cancers. Materials and methods MiRNA–target interactions The CLIP-supported miRNA–target pairs were downloaded from starBase V2.0 [73]. The target sites from at least one program (TargetScan [74] (Release 6.2), miRanda [75], Pictar2.0 [76], PITA [77] or RNA22 [78]) that resided within any entry of the Ago CLIP peaks were considered CLIP-supported sites. A miRNA–target pair was considered as the CLIP-supported pair if and only if the miRNA had at least one CLIP-supported site in the target. After converting the miRNA ID, 417 469 miRNA–target interactions between 383 miRNAs and 13 797 genes were finally reserved for subsequent analysis. High-quality experimentally validated miRNA–target interactions were obtained with R package multiMiR [79] from miRecords [80] miRTarbase [81] and TarBase [82], including 49 505 miRNA–target interactions between 13 523 genes and 459 miRNAs. Omics data across cancer types To identify miRNA–target regulation in each cancer type, we obtained paired miRNA-mRNA expression profiles as well as copy number and DNA methylation data from the TCGA Data Portal (2015). Based on the availability of these four types of molecular measurements and the sample size, 18 cancer types were analysed. Overall, 95.0% of the cancer types had >100 samples with the four abovementioned molecular measurements obtained at the same time, and thus, these relative sufficient data were used for further identification of miRNA–target regulation in the context of each cancer (Table 1). All miRNA expression data sets were obtained from TCGA Data Portal. MiRNA expression was profiled using microarrays in GBM and OVA and by small RNA sequencing in the remaining cancer types. For microarray data sets, TCGA Level 3 expression data were directly used. For miRNA sequencing data sets, expression values of mature miRNAs were used, which for a given miRNA was calculated as the ratio of reads mapping to it relative to the total number of reads mapping to annotated mature miRNAs in the given sample. We further removed miRNAs with a read count <10 that were detected in >95.0% of samples in each cancer type data set. The microarray and sequencing data expression values were log2 transformed for subsequent analysis. All mRNA expression data sets were obtained from TCGA Data Portal and profiled via RNA sequencing. For mRNA expression, mapped and gene-level-summarized (Level 3, TPM) RNA-seq data sets were used. We further removed mRNAs with a read count <20 that were detected in >95.0% of samples in each cancer type data set. The expression values were log2 transformed for subsequent analysis. DNA methylation data sets were also obtained from TCGA Data Portal (Level 3), and the methylation values for each protein-coding gene were defined as the average beta-values of probes mapping to the corresponding gene promoter (±2 kb of annotated transcription start sites). In addition, DNA copy number data sets were obtained from Firehose (http://gdac.broadinstitute.org/runs/analyses__2015_04_02/). We used Level 4 non-discretized gene-summarized log2-transformed aCGH copy number calls (tumour/normal ratio) computed by the Gistic2 algorithm [83]. Functional annotation of target genes The Biological Process (BP) terms for GO were downloaded from the MSigDB (v5.1) database [84]. As in previous studies, process categories from GO were restricted to BP terms such that the number of genes annotated to a term was at least 5 and no >500. Ultimately, 792 filtered GO BP terms were used for further analysis. Moreover, a list of GO terms determined to be related to the hallmarks of cancer was obtained from a previous study [85]. Protein–protein interaction network of target genes We assembled protein–protein interaction data from HPRD [86] and further removed self-loop interactions. Then, the gene symbols for each interaction were mapped to their corresponding Entrez gene identifiers. Finally, the maximum component of the whole protein–protein interaction network was extracted, which contained 35 865 interactions among 9028 genes. Collection of cancer-associated miRNAs Several database systems have chosen to provide a comprehensive resource indicating miRNA dysregulation in various human diseases. We manually collected the cancer-associated miRNAs from four databases, including HMDD [87], miR2Disease [88], miRCancer [89] and OncomiRDB [44] (Supplementary Table S5). Overview of the identification of miRNA cooperation across cancer types A multiple-step model was proposed to construct miRNA-miRNA functional cooperative networks in each cancer type. First, for each miRNA, a multivariate linear regression model was used to identify cancer context miRNA–target regulation (Supplementary Figure S1A). Then, the functionally cooperative miRNA pairs were identified as follows: for each miRNA pair, we initially identified miRNA pairs that significantly shared target genes; the shared target genes of the miRNA pair were used to identify candidate functional modules; and then, the candidate module sets were further filtered using two topological features in the protein interaction network (Supplementary Figure S1D). Here, we defined a pair of miRNAs as cooperative if they significantly co-regulated at least one functional module. Finally, all cooperative miRNA pairs were assembled into a functional miRNA cooperative network in each cancer type, where nodes represent miRNAs and edges represent their functional cooperative interactions. Identification of miRNA–target regulation in each cancer type For identifying cancer-specific miRNA–target regulation, we first measured the association between miRNA and mRNA expression for each miRNA-mRNA pair across a set of tumours in individual cancer types using a multivariate linear model that also took into account variation in mRNA expression induced by changes in DNA copy number and promoter methylation at the mRNA gene (Supplementary Figure S1A). For each predicted miRNA-mRNA regulation of miRNA u and mRNA i in a specific cancer, mRNA i expression y changes as a linear function of miRNA u expression x, DNA copy number x and DNA methylation x across tumour samples of a given cancer type: yi=βuxu,i+βcnxcn,i+βmexme,i+ β0+εi, i=1,…, n. In this linear model, β, βcn and β are the corresponding regression coefficients for miRNA expression, copy number and methylation variables, respectively, and β is the intercept. All P-values of the miRNA coefficient β in each cancer were subject to FDR correction. Under FDR < 0.05 for the signification of miRNA coefficients, the candidate cancer context miRNA-mRNA interactions were identified. Last, only candidate miRNA–target pairs with both negative correlation and CLIP-supported target sites were considered as cancer context miRNA regulation. Identification of functionally cooperative miRNA pairs across 18 cancer types Next, the functionally cooperative miRNA pairs were identified according to the protocol described in our previous study [9], which briefly included three steps: first, we initially identified miRNA pairs that share target mRNAs. Second, the shared target mRNAs of the miRNA pair were used to identify candidate functional modules. Third, the candidate module sets were further filtered using two topological features in the protein interaction network. A pair of miRNAs was defined as cooperative only if they significantly co-regulated at least one functional module. For a given pair of miRNAs, A and B, we first obtained their co-regulated target subset (A∩B), which was required to have at least three genes. Then, target subset functional enrichment was performed using the hypergeometric test across each selected GO term. At a given significance level, we not only obtained the enriched GO terms but also captured the subset in A∩B annotated to each term, GAB, the set of candidate function modules. We next further filtered these candidate function modules according to two topological features in the protein interaction network: (i) the minimum distance from every gene to others in the subset is no >2; (ii) the characteristic path length (CPL) is significantly shorter than random. The P-value indicating significance was calculated using the edge-switching method and was defined as the fraction of the CPL for the same subset in random networks that was shorter than that in the real network. We generated 1000 random networks using the software Mfinder (available at http://www.weizmann.ac.il/mcb/UriAlon/). After performing function enrichment and two topological restrictions in the network, miRNAs A and B were considered cooperative if they co-regulated at least one functional module. Topological measurements of the functional miRNA cooperative network In this study, we investigated several topological features of the miRNA-miRNA cooperative network using the R package ‘igraph’. First, we evaluated whether the degree distribution of the cooperative network satisfied a power-law model. For each miRNA in a miRNA-miRNA cooperative network, degree was defined as the number of edges linked to it. Clustering coefficient was defined as the ratio between the number of edges linking adjacent nodes and the total number of possible edges among them [90]. Classification of miRNA hubs Furthermore, we selected the top 30.0% of miRNAs with the highest degrees in the pan-cancer cooperative network as the hub miRNAs. To systematically analyse the hubs across cancers, we split the hubs into three groups: (i) cancer-specific hubs, in which the miRNAs were only hubs in the cooperative network of only one cancer type; (ii) common miRNA hubs, in which the miRNAs were hubs in more than one cancer type; and (iii) others. To define the first category, miRNAs belonging to hubs in the pan-cancer network that were ranked in the top 15.0% in one network but were not ranked in the top 15.0% in any other cancers were selected. Conversely, common hubs were required to be hubs in at least one cancer. The remaining miRNA hubs belonged to others. Identification of hub miRNA biomarkers For each cancer, we randomly selected tumour samples with the same number of normal samples ten times. Each time, t-test was used to estimate differential expression. The miRNAs with an FDR < 0.05 all 10 times were considered to be differentially expressed. To evaluate the potential of hub miRNAs with differential expression to be biomarkers, a scoring classifier was constructed. For each cancer, we first performed Z-score transformation on the expression levels across the samples for each hub miRNA with differential expression and then summarized the Z-scores as the integrated expression signature (score). Then, the samples could be divided into two classes (normal and tumour) by choosing a cut-off. In addition, an ROC curve that was drawn by plotting sensitivity against specificity was used for classifier evaluation. This procedure was performed using the R package pROC [91]. In addition, cluster analysis was performed using the R package pheatmap with default parameters. Identification of miRNA cooperative modules Next, miRNA functional cooperative modules, which were defined as cliques, were identified in the pan-cancer miRNA cooperative network using the clique percolation clustering method [92]. Cliques are all of the complete subgraphs, and all cliques in the network were identified using CFinder (2.0.6) [50]. Identification of survival-related miRNA cooperative modules To evaluate the prognostic performance of our miRNA cooperative modules identified above in each cancer, we first used univariate Cox regression analysis to evaluate the association between survival and the expression level of each miRNA in the module. A regression coefficient with a plus sign indicated that increased expression is associated with an increased risk of survival (risk miRNAs), and conversely, a minus sign indicated that increased expression is associated with a decreased risk of survival (protective miRNAs). After selecting the miRNAs that were significantly correlated with survival (P < 0.05) in the module, a mathematical formula for survival prediction was constructed, taking into account both the strength and positive or negative association of each miRNA with survival. More specifically, we assigned a risk score to each patient according to a linear combination of the miRNA expression values ExpmiRNAi weighted by the regression coefficients from the aforementioned univariate Cox regression analysis. The risk score for each patient was calculated as follows: Risk_score=∑i=1nβi*ExpmiRNAi, where βi is the Cox regression coefficient of miRNA i, and n is the number of miRNAs that were significantly associated with survival in the module. All patients were thus dichotomized into high-risk and low-risk groups using the median risk score as the cut-off point. Patients with higher risk scores were expected to have poor survival outcomes. The Kaplan–Meier method was further used to estimate the overall survival time for the two subgroups. The differences in the survival times were analysed using a log-rank test. Randomization test To determine if the miRNA-miRNA cooperative network is a small-world network, we used the duplication model to construct random graphs with the function of sample_gnm() in R. This is a well-known model with power-law degree distributions that provides small-world networks. We generated 1000 instances and computed the average clustering coefficients. The significant P-value was the fraction of the average clustering coefficients in random conditions, which is greater than the value in the real condition. Statistical analysis Because the regression coefficients do not follow a normal distribution, we used a Wilcoxon rank-sum test to evaluate differences in the mean of regression coefficients of cooperative miRNAs from that of non-cooperative miRNAs in each network. The miRNAs were descending sorted by expression variance, and the top 20.0% of the miRNAs were defined as miRNAs with greater variance. Fisher's exact test was used to determine if there were non-random associations between two cooperative properties and variance of expression. Key Points The miRNA-miRNA cooperative landscape was surveyed across cancer types. Pan-cancer miRNA-miRNA cooperative network exhibited a scale-free and modular architecture. Cooperative miRNAs exhibit higher expression, greater expression variation and a stronger regulatory strength towards targets as well as regulated cancer hallmark functions. Hub miRNAs regulate many cancer-associated functions and tend to target known anticancer drug targets. miRNA cooperative modules were found to serve as novel prognostic biomarkers. Supplementary Data Supplementary data are available online at https://academic.oup.com/bib. Funding The National High Technology Research and Development Program of China (863 Program, grant number 2014AA021102), the National Program on Key Basic Research Project (973 Program, grant number 2014CB910504), the funds for Creative Research Groups of the National Natural Science Foundation of China (grant number 81421063), the National Natural Science Foundation of China (grant numbers 91439117, 61473106, 31571331, 61502126, 31601065); the China Postdoctoral Science Foundation (grant numbers 2017M611385); Heilongjiang Province Postdoctoral Fund (grant number LBH-Z16152), Weihan Yu Youth Science Fund Project of Harbin Medical University, Harbin Special Funds of Innovative Talents on Science and Technology Research Project (grant number RC2015QN003080), the Harbin Medical University Innovation Science Research Fund (grant number 2016JCZX45) and the Funds for the Graduate Innovation Fund of Heilongjiang Province (grant number YJSCX2015-7HYD). Tingting Shao is an instructor in the College of Bioinformatics Science and Technology at Harbin Medical University, China. Her research activity focused on ncRNA regulation in complex diseases. Guangjuan Wang is an MS student in the College of Bioinformatics Science and Technology at Harbin Medical University, China. Her research interests focus on bioinformatics methods. Hong Chen is a PhD in the College of Bioinformatics Science and Technology at Harbin Medical University, China. Her research focuses on bioinformatics. Yunjin Xie is an undergraduate in the College of Bioinformatics Science and Technology at Harbin Medical University, China. Xiyun Jin is a master student in the College of Bioinformatics Science and Technology at Harbin Medical University, China. Her research focuses on bioinformatics. Jing Bai is an instructor in the College of Bioinformatics Science and Technology at Harbin Medical University, China. Her research interests focus on method development and ncRNA regulation. Juan Xu is an associate professor in the College of Bioinformatics Science and Technology at Harbin Medical University, China. Her researches focus on computational biology in complex diseases. Xia Li is a professor in the College of Bioinformatics Science and Technology at Harbin Medical University, China. Her research interests focus on ncRNA regulation in complex diseases. Jian Huang is a professor in the Center for Informational Biology at University of Electronic Science and Technology of China, China Yan Jin is a professor in the Department of Medical Genetics at Harbin Medical University, China. Yongsheng Li is an associate professor in the College of Bioinformatics Science and Technology at Harbin Medical University, China. His research interests focus on ncRNA regulation and bioinformatics methods development. References 1 Fujii YR. Quantum Language of MicroRNA: application for new cancer therapeutic targets . Methods Mol Biol 2018 ; 1733 : 145 – 57 . Google Scholar CrossRef Search ADS PubMed 2 He Z , Wang Y , Huang G , et al. The lncRNA UCA1 interacts with miR-182 to modulate glioma proliferation and migration by targeting iASPP . Arch Biochem Biophysics 2017 ; 623–624 : 1 – 8 . Google Scholar CrossRef Search ADS 3 Feng J , Xue S , Pang Q , et al. miR-141-3p inhibits fibroblast proliferation and migration by targeting GAB1 in keloids . Biochem Biophys Res Commun 2017 ; 490 ( 2 ): 302 – 8 . Google Scholar CrossRef Search ADS PubMed 4 Ziebarth JD , Bhattacharya A , Chen A , et al. PolymiRTS database 2.0: linking polymorphisms in microRNA target sites with human diseases and complex traits . Nucleic Acids Res 2012 ; 40 ( D1 ): D216 – 21 . Google Scholar CrossRef Search ADS PubMed 5 Deng Y , Bai H , Hu H. rs11671784 G/A variation in miR-27a decreases chemo-sensitivity of bladder cancer by decreasing miR-27a and increasing the target RUNX-1 expression . Biochem Biophys Res Commun 2015 ; 458 ( 2 ): 321 – 7 . Google Scholar CrossRef Search ADS PubMed 6 Naveed A , Ur-Rahman S , Abdullah S , et al. A concise review of MicroRNA exploring the insights of MicroRNA regulations in bacterial, viral and metabolic diseases . Mol Biotechnol 2017 ; 59 ( 11–12 ): 518 – 29 . Google Scholar CrossRef Search ADS PubMed 7 Kim D , Chang HR , Baek D. Rules for functional microRNA targeting . BMB Rep 2017 ; 50 ( 11 ): 554 – 9 . Google Scholar CrossRef Search ADS PubMed 8 Kozomara A , Griffiths-Jones S. miRBase: annotating high confidence microRNAs using deep sequencing data . Nucleic Acids Res 2014 ; 42 ( D1 ): D68 – 73 . Google Scholar CrossRef Search ADS PubMed 9 Xu J , Li CX , Li YS , et al. MiRNA-miRNA synergistic network: construction via co-regulating functional modules and disease miRNA topological features . Nucleic Acids Res 2011 ; 39 ( 3 ): 825 – 36 . Google Scholar CrossRef Search ADS PubMed 10 Wu S , Huang S , Ding J , et al. Multiple microRNAs modulate p21Cip1/Waf1 expression by directly targeting its 3' untranslated region . Oncogene 2010 ; 29 ( 15 ): 2302 – 8 . Google Scholar CrossRef Search ADS PubMed 11 Bandi N , Vassella E. miR-34a and miR-15a/16 are co-regulated in non-small cell lung cancer and control cell cycle progression in a synergistic and Rb-dependent manner . Mol Cancer 2011 ; 10 ( 1 ): 55 . Google Scholar CrossRef Search ADS PubMed 12 Kandi R , Gutti U , Saladi RGV , et al. MiR-125b and miR-99a encoded on chromosome 21 co-regulate vincristine resistance in childhood acute megakaryoblastic leukemia . Hematol Oncol Stem Cell Ther 2015 ; 8 ( 2 ): 95 – 7 . Google Scholar CrossRef Search ADS PubMed 13 Zhao Y , Cui X , Zhu W , et al. Synergistic regulatory effects of microRNAs on brain glioma cells . Mol Med Rep 2017 ; 16 ( 2 ): 1409 – 16 . Google Scholar CrossRef Search ADS PubMed 14 Borzi C , Calzolari L , Centonze G , et al. mir-660-p53-mir-486 network: a new key regulatory pathway in lung tumorigenesis . Int J Mol Sci 2017 ; 18 ( 1 ): 222 . Google Scholar CrossRef Search ADS 15 Lai X , Gupta SK , Schmitz U , et al. MiR-205-5p and miR-342-3p cooperate in the repression of the E2F1 transcription factor in the context of anticancer chemotherapy resistance . Theranostics 2018 ; 8 ( 4 ): 1106 – 20 . Google Scholar CrossRef Search ADS PubMed 16 Frampton AE , Castellano L , Colombo T , et al. MicroRNAs cooperatively inhibit a network of tumor suppressor genes to promote pancreatic tumor growth and progression . Gastroenterology 2014 ; 146 ( 1 ): 268 – 77.e18 . Google Scholar CrossRef Search ADS PubMed 17 Chen HH , Huang WT , Yang LW , et al. The PTEN-AKT-mTOR/RICTOR pathway in nasal natural killer cell lymphoma is activated by miR-494-3p via PTEN but inhibited by miR-142-3p via RICTOR . Am J Pathol 2015 ; 185 ( 5 ): 1487 – 99 . Google Scholar CrossRef Search ADS PubMed 18 Xu J , Shao T , Ding N , et al. miRNA–miRNA crosstalk: from genomics to phenomics . Brief Bioinform 2017 ; 18 ( 6 ): 1002 – 11 . Google Scholar PubMed 19 Ying H , Lyu J , Ying T , et al. RETRACTED ARTICLE: risk miRNA screening of ovarian cancer based on miRNA functional synergistic network . J Ovarian Res 2014 ; 7 ( 1 ): 226 . Google Scholar CrossRef Search ADS 20 Alshalalfa M. MicroRNA response elements-mediated miRNA-miRNA interactions in prostate cancer . Adv Bioinform 2012 ; 2012 : 839837 . 21 Hua L , Zhou P , Li L , et al. Prioritizing breast cancer subtype related miRNAs using miRNA–mRNA dysregulated relationships extracted from their dual expression profiling . J Theor Biol 2013 ; 331 : 1 – 11 . Google Scholar CrossRef Search ADS PubMed 22 Yang Y , Xing Y , Liang C , et al. Crucial microRNAs and genes of human primary breast cancer explored by microRNA-mRNA integrated analysis . Tumor Biol 2015 ; 36 ( 7 ): 5571 – 9 . Google Scholar CrossRef Search ADS 23 Zhang TF , Cheng KW , Shi WY , et al. MiRNA synergistic network construction and enrichment analysis for common target genes in small-cell lung cancer . Asian Pac J Cancer Prev 2012 ; 13 ( 12 ): 6375 – 8 . Google Scholar CrossRef Search ADS PubMed 24 Li Y , Liang C , Wong K-C , et al. Mirsynergy: detecting synergistic miRNA regulatory modules by overlapping neighbourhood expansion . Bioinformatics 2014 ; 30 ( 18 ): 2627 – 35 . Google Scholar CrossRef Search ADS PubMed 25 Zhao X , Song H , Zuo Z , et al. Identification of miRNA-miRNA synergistic relationships in colorectal cancer . Int J Biol Macromol 2013 ; 55 : 98 – 103 . Google Scholar CrossRef Search ADS PubMed 26 Na YJ , Kim JH. Understanding cooperativity of microRNAs via microRNA association networks . BMC Genomics 2013 ; 14 ( Suppl 5 ): S17 . Google Scholar CrossRef Search ADS PubMed 27 Xiao Y , Ping Y , Fan H , et al. Identifying dysfunctional miRNA-mRNA regulatory modules by inverse activation, cofunction, and high interconnection of target genes: a case study of glioblastoma . Neuro Oncol 2013 ; 15 ( 7 ): 818 – 28 . Google Scholar CrossRef Search ADS PubMed 28 Zhu W , Zhao Y , Xu Y , et al. Dissection of protein interactomics highlights microRNA synergy . PLoS One 2013 ; 8 ( 5 ): e63342 . Google Scholar CrossRef Search ADS PubMed 29 Meng X , Wang J , Yuan C , et al. CancerNet: a database for decoding multilevel molecular interactions across diverse cancer types . Oncogenesis 2015 ; 4 : e177. Google Scholar CrossRef Search ADS PubMed 30 Song R , Catchpoole DR , Kennedy PJ , et al. Identification of lung cancer miRNA-miRNA co-regulation networks through a progressive data refining approach . J Theor Biol 2015 ; 380 : 271 – 9 . Google Scholar CrossRef Search ADS PubMed 31 The Cancer Genome Atlas Research Network , Weinstein JN , Collisson EA , et al. The Cancer Genome Atlas Pan-Cancer analysis project . Nat Genet 2013 ; 45 : 1113 – 20 . Google Scholar CrossRef Search ADS PubMed 32 The International Cancer Genome Consortium , Hudson TJ , Anderson W , et al. International network of cancer genome projects . Nature 2010 ; 464 : 993 – 8 . Google Scholar CrossRef Search ADS PubMed 33 Jacobsen A , Silber J , Harinath G , et al. Analysis of microRNA-target interactions across diverse cancer types . Nat Struct Mol Biol 2013 ; 20 ( 11 ): 1325 – 32 . Google Scholar CrossRef Search ADS PubMed 34 Yang JH , Li JH , Shao P , et al. starBase: a database for exploring microRNA-mRNA interaction maps from Argonaute CLIP-Seq and Degradome-Seq data . Nucleic Acids Res 2011 ; 39 : D202 – 9 . Google Scholar CrossRef Search ADS PubMed 35 Chi SW , Zang JB , Mele A , et al. Argonaute HITS-CLIP decodes microRNA-mRNA interaction maps . Nature 2009 ; 460 ( 7254 ): 479 – 86 . Google Scholar CrossRef Search ADS PubMed 36 Luna JM , Barajas JM , Teng KY , et al. Argonaute CLIP defines a deregulated miR-122-bound transcriptome that correlates with patient survival in human liver cancer . Mol Cell 2017 ; 67 ( 3 ): 400 – 10.e7 . Google Scholar CrossRef Search ADS PubMed 37 Zhang XQ , Yang JH. Discovering circRNA-microRNA interactions from CLIP-Seq Data . Methods Mol Biol 2018 ; 1724 : 193 – 207 . Google Scholar CrossRef Search ADS PubMed 38 Korpal M , Lee ES , Hu G , et al. The miR-200 family inhibits epithelial-mesenchymal transition and cancer cell migration by direct targeting of E-cadherin transcriptional repressors ZEB1 and ZEB2 . J Biol Chem 2008 ; 283 ( 22 ): 14910 – 4 . Google Scholar CrossRef Search ADS PubMed 39 Yoneyama K , Ishibashi O , Kawase R , et al. miR-200a, miR-200b and miR-429 are onco-miRs that target the PTEN gene in endometrioid endometrial carcinoma . Anticancer Res 2015 ; 35 : 1401 – 10 . Google Scholar PubMed 40 Hoadley KA , Yau C , Wolf DM , et al. Multiplatform analysis of 12 cancer types reveals molecular classification within and across tissues of origin . Cell 2014 ; 158 : 929 – 44 . Google Scholar CrossRef Search ADS PubMed 41 Rosenfeld N , Aharonov R , Meiri E , et al. MicroRNAs accurately identify cancer tissue origin . Nat Biotechnol 2008 ; 26 ( 4 ): 462 – 9 . Google Scholar CrossRef Search ADS PubMed 42 Xu Q , Chen J , Ni S , et al. Pan-cancer transcriptome analysis reveals a gene expression signature for the identification of tumor tissue origin . Mod Pathol 2016 ; 29 ( 6 ): 546 – 56 . Google Scholar CrossRef Search ADS PubMed 43 Mar JC , Matigian NA , Mackay-Sim A , et al. Variance of gene expression identifies altered network constraints in neurological disease . PLoS Genetics 2011 ; 7 ( 8 ): e1002207 . Google Scholar CrossRef Search ADS PubMed 44 Wang D , Gu J , Wang T , et al. OncomiRDB: a database for the experimentally verified oncogenic and tumor-suppressive microRNAs . Bioinformatics 2014 ; 30 ( 15 ): 2237 – 8 . Google Scholar CrossRef Search ADS PubMed 45 Wang E , Zaman N , McGee S , et al. Predictive genomics: a cancer hallmark network framework for predicting tumor clinical phenotypes using genome sequencing data . Semin Cancer Biol 2015 ; 30 : 4 – 12 . Google Scholar CrossRef Search ADS PubMed 46 Zotenko E , Mestre J , O'Leary DP , et al. Why do hubs in the yeast protein interaction network tend to be essential: reexamining the connection between the network topology and essentiality . PLoS Comput Biol 2008 ; 4 ( 8 ): e1000140 . Google Scholar CrossRef Search ADS PubMed 47 Lu M , Shi B , Wang J , et al. TAM: a method for enrichment and depletion analysis of a microRNA category in a list of microRNAs . BMC Bioinformatics 2010 ; 11 ( 1 ): 419 . Google Scholar CrossRef Search ADS PubMed 48 Li J , Lenferink AEG , Deng Y , et al. Corrigendum: identification of high-quality cancer prognostic markers and metastasis network modules . Nat Commun 2012 ; 3 ( 1 ): 655 . Google Scholar CrossRef Search ADS 49 Wu G , Stein L. A network module-based method for identifying cancer prognostic signatures . Genome Biol 2012 ; 13 ( 12 ): R112. Google Scholar CrossRef Search ADS PubMed 50 Niklas N , Hafenscher J , Barna A , et al. cFinder: definition and quantification of multiple haplotypes in a mixed sample . BMC Res Notes 2015 ; 8 ( 1 ): 422 . Google Scholar CrossRef Search ADS PubMed 51 Watts DJ , Strogatz SH. Collective dynamics of ‘small-world’ networks . Nature 1998 ; 393 ( 6684 ): 440 – 2 . Google Scholar CrossRef Search ADS PubMed 52 Park S , Lehner B. Cancer type‐dependent genetic interactions between cancer driver alterations indicate plasticity of epistasis across cell types . Mol Syst Biol 2015 ; 11 ( 7 ): 824. Google Scholar CrossRef Search ADS PubMed 53 Xu J , Li Y , Lu J , et al. The mRNA related ceRNA-ceRNA landscape and significance across 20 major cancer types . Nucleic Acids Res 2015 ; 43 ( 17 ): 8169 – 82 . Google Scholar CrossRef Search ADS PubMed 54 Hamilton MP , Rajapakshe K , Hartig SM , et al. Identification of a pan-cancer oncogenic microRNA superfamily anchored by a central core seed motif . Nat Commun 2013 ; 4 : 2730. Google Scholar CrossRef Search ADS PubMed 55 Davoli T , Xu Andrew W , Mengwasser KE , et al. Cumulative haploinsufficiency and triplosensitivity drive aneuploidy patterns and shape the cancer genome . Cell 2013 ; 155 ( 4 ): 948 – 62 . Google Scholar CrossRef Search ADS PubMed 56 Heim D , Budczies J , Stenzinger A , et al. Cancer beyond organ and tissue specificity: next‐generation‐sequencing gene mutation data reveal complex genetic similarities across major cancers . Int J Cancer 2014 ; 135 ( 10 ): 2362 – 9 . Google Scholar CrossRef Search ADS PubMed 57 Zhao X , Modur V , Carayannopoulos LN , et al. Biomarkers in pharmaceutical research . Clin Chem 2015 ; 61 ( 11 ): 1343 – 53 . Google Scholar CrossRef Search ADS PubMed 58 Saito T , Sætrom P. Target gene expression levels and competition between transfected and endogenous microRNAs are strong confounding factors in microRNA high-throughput experiments . Silence 2012 ; 3 ( 1 ): 3. Google Scholar CrossRef Search ADS PubMed 59 Mendell JT. miRiad roles for the miR-17-92 cluster in development and disease . Cell 2008 ; 133 ( 2 ): 217 – 22 . Google Scholar CrossRef Search ADS PubMed 60 Di Lena M , Travaglio E , Altomare DF. New strategies for colorectal cancer screening . World J Gastroenterol 2013 ; 19 ( 12 ): 1855 – 60 . Google Scholar CrossRef Search ADS PubMed 61 Wu HH , Lin WC , Tsai KW. Advances in molecular biomarkers for gastric cancer: miRNAs as emerging novel cancer markers . Expert Rev Mol Med 2014 ; 16 : e1 . Google Scholar CrossRef Search ADS PubMed 62 Chen X , Ba Y , Ma L , et al. Characterization of microRNAs in serum: a novel class of biomarkers for diagnosis of cancer and other diseases . Cell Res 2008 ; 18 ( 10 ): 997 – 1006 . Google Scholar CrossRef Search ADS PubMed 63 Cui L , Zhang X , Ye G , et al. Gastric juice MicroRNAs as potential biomarkers for the screening of gastric cancer . Cancer 2013 ; 119 ( 9 ): 1618 – 26 . Google Scholar CrossRef Search ADS PubMed 64 Zhang X , Cui L , Ye G , et al. Gastric juice microRNA-421 is a new biomarker for screening gastric cancer, Tumour biology: the journal of the International Society for Oncodevelopmental . Biol Med 2012 ; 33 ( 6 ): 2349 – 55 . 65 Yu X , Luo L , Wu Y , et al. Gastric juice miR-129 as a potential biomarker for screening gastric cancer . Med Oncol 2013 ; 30 ( 1 ): 365 . Google Scholar CrossRef Search ADS PubMed 66 Li C , Li JF , Cai Q , et al. MiRNA-199a-3p in plasma as a potential diagnostic biomarker for gastric cancer . Ann Surg Oncol 2013 ; 20 ( S3 ): 397 – 405 . Google Scholar CrossRef Search ADS 67 Li BS , Zhao Yl , Guo G , et al. Plasma microRNAs, miR-223, miR-21 and miR-218, as novel potential biomarkers for gastric cancer detection . PLoS One 2012 ; 7 ( 7 ): e41629 . Google Scholar CrossRef Search ADS PubMed 68 Yun SJ , Jeong P , Kim WT , et al. Cell-free microRNAs in urine as diagnostic and prognostic biomarkers of bladder cancer . Int J Oncol 2012 ; 41 ( 5 ): 1871 – 8 . Google Scholar CrossRef Search ADS PubMed 69 Park NJ , Zhou H , Elashoff D , et al. Salivary microRNA: discovery, characterization, and clinical utility for oral cancer detection . Clin Cancer Res 2009 ; 15 ( 17 ): 5473 . Google Scholar CrossRef Search ADS PubMed 70 Gallardo E , Navarro A , Viñolas N , et al. miR-34a as a prognostic marker of relapse in surgically resected non-small-cell lung cancer . Carcinogenesis 2009 ; 30 ( 11 ): 1903 – 9 . Google Scholar CrossRef Search ADS PubMed 71 Wang Q , Li P , Li A , et al. Plasma specific miRNAs as predictive biomarkers for diagnosis and prognosis of glioma . J Exp Clin Cancer Res 2012 ; 31 ( 1 ): 97 . Google Scholar CrossRef Search ADS PubMed 72 Zhang L , Shan X , Wang J , et al. A three-microRNA signature for lung squamous cell carcinoma diagnosis in Chinese male patients . Oncotarget 2017 ; 8 ( 49 ): 86897 – 907 . Google Scholar PubMed 73 Li JH , Liu S , Zhou H , et al. starBase v2.0: decoding miRNA-ceRNA, miRNA-ncRNA and protein-RNA interaction networks from large-scale CLIP-Seq data . Nucleic Acids Res 2014 ; 42 ( D1 ): D92 – 7 . Google Scholar CrossRef Search ADS PubMed 74 Grimson A , Farh KK , Johnston WK , et al. MicroRNA targeting specificity in mammals: determinants beyond seed pairing . Mol Cell 2007 ; 27 ( 1 ): 91 – 105 . Google Scholar CrossRef Search ADS PubMed 75 Betel D , Koppal A , Agius P , et al. Comprehensive modeling of microRNA targets predicts functional non-conserved and non-canonical sites . Genome Biol 2010 ; 11 ( 8 ): R90 . Google Scholar CrossRef Search ADS PubMed 76 Anders G , Mackowiak SD , Jens M , et al. doRiNA: a database of RNA interactions in post-transcriptional regulation . Nucleic Acids Res 2012 ; 40 ( D1 ): D180 – 6 . Google Scholar CrossRef Search ADS PubMed 77 Kertesz M , Iovino N , Unnerstall U , et al. The role of site accessibility in microRNA target recognition . Nat Genet 2007 ; 39 ( 10 ): 1278 – 84 . Google Scholar CrossRef Search ADS PubMed 78 Miranda KC , Huynh T , Tay Y , et al. A pattern-based method for the identification of MicroRNA binding sites and their corresponding heteroduplexes . Cell 2006 ; 126 ( 6 ): 1203 – 17 . Google Scholar CrossRef Search ADS PubMed 79 Ru Y , Kechris KJ , Tabakoff B , et al. The multiMiR R package and database: integration of microRNA–target interactions along with their disease and drug associations . Nucleic Acids Res 2014 ; 42 ( 17 ): e133 . Google Scholar CrossRef Search ADS PubMed 80 Xiao F , Zuo Z , Cai G , et al. miRecords: an integrated resource for microRNA-target interactions . Nucleic Acids Res 2009 ; 37 : D105 – 10 . Google Scholar CrossRef Search ADS PubMed 81 Hsu SD , Lin FM , Wu WY , et al. miRTarBase: a database curates experimentally validated microRNA–target interactions . Nucleic Acids Res 2011 ; 39(Suppl 1) : D163 – 9 . Google Scholar CrossRef Search ADS 82 Sethupathy P , Corda B , Hatzigeorgiou AG. TarBase: a comprehensive database of experimentally supported animal microRNA targets . RNA 2005 ; 12 ( 2 ): 192 – 7 . Google Scholar CrossRef Search ADS PubMed 83 Mermel CH , Schumacher SE , Hill B , et al. GISTIC2.0 facilitates sensitive and confident localization of the targets of focal somatic copy-number alteration in human cancers . Genome Biol 2011 ; 12 ( 4 ): R41 . Google Scholar CrossRef Search ADS PubMed 84 Liberzon A , Subramanian A , Pinchback R , et al. Molecular signatures database (MSigDB) 3.0 . Bioinformatics 2011 ; 27 ( 12 ): 1739 – 40 . Google Scholar CrossRef Search ADS PubMed 85 Plaisier CL , Pan M , Baliga NS. A miRNA-regulatory network explains how dysregulated miRNAs perturb oncogenic processes across diverse cancers . Genome Res 2012 ; 22 ( 11 ): 2302 – 14 . Google Scholar CrossRef Search ADS PubMed 86 Keshava Prasad TS , Goel R , Kandasamy K , et al. Human protein reference database–2009 update . Nucleic Acids Res 2009 ; 37 : D767 – 72 . Google Scholar CrossRef Search ADS PubMed 87 Li Y , Qiu C , Tu J , et al. HMDD v2.0: a database for experimentally supported human microRNA and disease associations . Nucleic Acids Res 2014 ; 42 : D1070 – 4 . Google Scholar CrossRef Search ADS PubMed 88 Jiang Q , Wang Y , Hao Y , et al. miR2Disease: a manually curated database for microRNA deregulation in human disease . Nucleic Acids Res 2009 ; 37 : D98 – 104 . Google Scholar CrossRef Search ADS PubMed 89 Xie B , Ding Q , Han H , et al. miRCancer: a microRNA-cancer association database constructed by text mining on literature . Bioinformatics 2013 ; 29 ( 5 ): 638 – 44 . Google Scholar CrossRef Search ADS PubMed 90 Hsu CW , Juan HF , Huang HC. Characterization of microRNA-regulated protein-protein interaction network . Proteomics 2008 ; 8 ( 10 ): 1975 – 9 . Google Scholar CrossRef Search ADS PubMed 91 Robin X , Turck N , Hainard A , et al. pROC: an open-source package for R and S+ to analyze and compare ROC curves . BMC Bioinformatics 2011 ; 12 ( 1 ): 77 . Google Scholar CrossRef Search ADS PubMed 92 Palla G , Derenyi I , Farkas I , et al. Uncovering the overlapping community structure of complex networks in nature and society . Nature 2005 ; 435 ( 7043 ): 814 – 18 . Google Scholar CrossRef Search ADS PubMed © The Author(s) 2018. Published by Oxford University Press. All rights reserved. For permissions, please email: journals.permissions@oup.com This article is published and distributed under the terms of the Oxford University Press, Standard Journals Publication Model (https://academic.oup.com/journals/pages/about_us/legal/notices) http://www.deepdyve.com/assets/images/DeepDyve-Logo-lg.png Briefings in Bioinformatics Oxford University Press

Survey of miRNA-miRNA cooperative regulation principles across cancer types

Loading next page...
 
/lp/ou_press/survey-of-mirna-mirna-cooperative-regulation-principles-across-cancer-6T03GU7T0Y
Publisher
Oxford University Press
Copyright
© The Author(s) 2018. Published by Oxford University Press. All rights reserved. For permissions, please email: journals.permissions@oup.com
ISSN
1467-5463
eISSN
1477-4054
D.O.I.
10.1093/bib/bby038
Publisher site
See Article on Publisher Site

Abstract

Abstract Cooperative regulation among multiple microRNAs (miRNAs) is a complex type of posttranscriptional regulation in human; however, the global view of the system-level regulatory principles across cancers is still unclear. Here, we investigated miRNA-miRNA cooperative regulatory landscape across 18 cancer types and summarized the regulatory principles of miRNAs. The miRNA-miRNA cooperative pan-cancer network exhibited a scale-free and modular architecture. Cancer types with similar tissue origins had high similarity in cooperative network structure and expression of cooperative miRNA pairs. In addition, cooperative miRNAs showed divergent properties, including higher expression, greater expression variation and a stronger regulatory strength towards targets and were likely to regulate cancer hallmark-related functions. We found a marked rewiring of miRNA-miRNA cooperation between various cancers and revealed conserved and rewired network miRNA hubs. We further identified the common hubs, cancer-specific hubs and other hubs, which tend to target known anticancer drug targets. Finally, miRNA cooperative modules were found to be associated with patient survival in several cancer types. Our study highlights the potential of pan-cancer miRNA-miRNA cooperative regulation as a novel paradigm that may aid in the discovery of tumorigenesis mechanisms and development of anticancer drugs. pan-cancer, cancer context miRNA–target regulation, miRNA-miRNA cooperative regulatory network, hub miRNAs, miRNA cooperative modules Introduction MicroRNAs (miRNAs) are small non-coding RNAs of approximately 22 nucleotides in length that play critical roles in various types of biological processes and complex diseases, including cancer [1–5]. A couple of studies have demonstrated that miRNAs perform their functions by forming miRNA-induced silencing complexes, thereby hindering mRNA translation or accelerating RNA degradation [6, 7]. With the development of high-throughput sequencing technology, an increasing number of miRNAs have been identified [8]. However, our knowledge of their functions in cancer is still limited. In addition, miRNA regulation exhibits great complexity in cancer. One miRNA can target multiple genes, and one gene can be regulated by more than one miRNA, indicating cooperative regulation among miRNAs [9]. Cooperative regulation among miRNAs has been identified in multiple species and various types of cancer. Wu et al. [10] first found that miRNAs from different miRNA clusters cooperatively regulate target gene expression, which signalled a new era of miRNA regulation. Bandi et al. [11] demonstrated that hsa-miR-34a and hsa-miR-15a/16 act synergistically to induce cell cycle arrest in non-small cell lung cancer. Kandi et al. [12] found that hsa-miR-125b and hsa-miR-99a encoded on chromosome 21 co-regulate vincristine resistance in childhood acute megakaryoblastic leukaemia. In addition, it has been demonstrated that hsa-miR-21 and hsa-miR-20a exhibited the largest synergistic effects in brain glioma cells [13]. Borzi et al. [14] investigated the interaction between hsa-miR-486-5p and hsa-miR-660-5p and assess their possible role and synergistic effect in lung cancer treatment. In another study, Lai et al. [15] validated that hsa-miR-205-5p and hsa-miR-342-3p achieve an enhanced repressive effect on E2F1 expression in a cooperative manner in an aggressive human non-small cell lung carcinoma cell line. At the same time, several miRNAs have also been identified as cancer biomarkers, such as hsa-miR-21-5p, hsa-miR-23a-3p and hsa-miR-27a-3p in pancreatic tumours and hsa-miR-142-3p and has-miR-494-3p in lymphoma [16, 17]. These observations indicate the presence of cooperative miRNA-mediated regulation, and studying the potential functional effects of miRNA cooperative regulation is an interesting and challenging area in miRNA research. However, it is a challenge to identify cooperative miRNA regulation using experimental methods because of the high number of miRNA combinations. Moreover, miRNA regulation is often reprogrammed in different tissues or different cancer types. To address this issue, several bioinformatics analysis pipelines have been proposed for inferring cooperative regulation among miRNAs, ranging from genomics to phenomics-based methods [18]. Several studies first constructed miRNA-miRNA interaction network using predicted miRNA–target data and then focused on hub miRNAs. Ying et al. [19] constructed a miRNA-miRNA co-regulating network in ovarian cancer, and hubs in the network were selected and identified as risk miRNAs, such as hsa-miR-579, hsa-miR-942, hsa-miR-105 and hsa-miR-150. Alshalalfa et al. [20] investigated mRNA-mediated miRNA-miRNA interactions in prostate cancer, and 11 miRNAs highly connected to other miRNAs were identified as diagnostic and prognostic biomarkers. In another study, Hua et al. [21] identified the miRNA-miRNA network for each breast cancer subtype and found that luminal A and basal-like subtype-specific networks exhibited changes in hub miRNAs. Yang et al. [22] constructed a miRNA-miRNA functional synergistic network and miRNA-miRNA co-regulatory network in breast cancer and proposed hub miRNAs (such as, hsa-miR-139-5p, hsa-miR-9, hsa-miR-375, hsa-miR-592 and hsa-miR-135a) that may be crucial in the initiation and progression of primary breast cancer. In addition, several studies first constructed miRNA-miRNA interaction networks using predicted miRNA–target data and then focused on modules. Zhang et al. [23] analysed a miRNA synergistic regulatory network in small cell lung cancer and identified a network module that included three miRNAs (hsa-let-7c, hsa-let-7b and hsa-let-7d). In addition, Li et al. [24] detected synergistic miRNA regulatory modules using ‘Mirsynery’ in ovarian, breast and thyroid cancer and identified several promising cancer-specific modules as prognostic biomarkers. Zhao et al. [25] constructed a miRNA-miRNA synergistic network in colorectal cancer and identified eight functional modules containing miRNAs that could implement specific functions as a whole in colorectal cancer (CRC). Especially, we proposed a functional module-based method for identification of cooperative miRNAs, and this method has been widely used in various types of cancer [9]. It is might expected that miRNA cooperative regulation is reprogrammed in different biological contexts. Context-specific miRNA-miRNA regulation, which means that the miRNA-miRNA regulation only occurred in a specific context, will provide a better approach for inferring the function of miRNAs in cancer. Thus, paired miRNA/mRNA expression levels across cancer samples were integrated to identify context-specific cooperative regulation among miRNAs [21, 24, 26]. Increasing cancer-specific miRNA-miRNA networks were analysed in various types of cancer, such as breast cancer, glioma and prostate cancer [20, 22, 27]. However, many questions regarding common or specific cooperative regulatory mechanisms in different cancer types have not been fully addressed, such as which collaborative miRNA regulation networks occur across all cancer types or in a specific cancer type? What about the functions, regulatory roles and biological insights of the pan-cancer miRNA-miRNA regulation, which include the whole common or specific cooperation cross all the cancer types? Moreover, increasing evidence has indicated that miRNA–gene regulation is critical for inferring cooperative miRNA regulation in cancer. Zhu et al. [28] applied the databases miRSel and ExprTargetDB jointly to obtain reliable miRNA–gene interaction data and determine miRNA synergy. Meng et al. [29] proposed that identification of context-dependent miRNA-target interactions is important, and they used common miRNA target-prediction programs, two databases and miRNA and mRNA expression data to obtain miRNA–target interaction information and identify functionally synergistic miRNA pairs. Song et al. [30] obtained miRNA target data from the Targetscan database to construct a miRNA-miRNA co-regulation network. Yang et al. [22] used the databases miRecords and miRWalk to obtain miRNA-gene interaction information and identified co-regulatory miRNA-miRNA networks. Integration of miRNA/mRNA gene expression to characterize the negative correlation between miRNA and target gene expression is widely used to identify context-specific regulation, thus contributing to the identification of cooperative miRNA regulation [31]. In addition to miRNAs, gene expression is also regulated by genetic and epigenetic regulators, such as DNA methylation [16]. Thus, to comprehensively understand miRNA–gene regulation as well as the cooperative regulation among miRNAs, we need to integrate multiple omics data sets across cancer types. The efforts of cancer genome projects, such as TCGA and ICGC, provide a great opportunity for us to investigate cooperative regulation among miRNAs across various cancer types [31, 32]. To address these issues, in this study, we focused on pan-cancer miRNA-miRNA cooperative regulation by assembling 18 cancer context cooperative regulation networks to investigate the biological insights they might present. The global principles of miRNA cooperative regulation were summarized. Based on the structure of the miRNA cooperative regulatory network and expression of cooperative miRNA pairs, we attempted to reveal the cancer tissue origins at different levels. In addition, we found that miRNAs in the pan-cancer miRNA-miRNA cooperative regulatory network show certain meaningful, expressional and functional properties, indicating that they might play important roles in tumorigenesis. Through comparative analyses of the structure of the pan-cancer and cancer-specific miRNA-miRNA cooperative regulatory networks, we identified the common hubs, cancer-specific hubs and other hubs and subsequently explored the application of hub miRNAs in clinical diagnosis and anticancer drug development. Finally, several miRNA modules were identified that might serve as novel prognostic biomarkers. These analyses and validations demonstrate how cancer-associated miRNA-miRNA cooperative regulatory networks can be used to accelerate discovery of miRNA-based biomarkers and potential therapeutics. Systematic investigation of miRNA-miRNA cooperative regulatory networks across multiple cancers can help us elucidate the commonalities and differences in the mechanisms driving specific cancer types. Results The cancer context-specific miRNA–target regulation To systematically characterize miRNA–target regulation across multiple cancer types, we used a three-step pipeline to identify cancer context-specific miRNA–target regulation by integration of genome-wide miRNA–target regulation and paired miRNA-mRNA expression profiles as well as copy number and DNA methylation data across 18 cancer types (Supplementary Figure S1A and Methods). Each cancer type is with paired mRNA and miRNA expression profiles as well as copy number and DNA methylation data measured for the same samples, forming a pan-cancer data compendium from 6348 tumour samples for 20 472 mRNAs and 867 miRNAs (Table 1). These cancer types were primarily involved in eight classes, namely, brain (brain lower-grade glioma, LGG; glioblastoma multiforme, GBM), genitourinary (kidney renal clear cell carcinoma, KIRC; kidney renal papillary cell carcinoma, KIRP; bladder urothelial carcinoma, BLCA; prostate adenocarcinoma, PRAD), thoracic (lung squamous cell carcinoma, LUSC; lung adenocarcinoma, LUAD), gastrointestinal (colon adenocarcinoma, COAD; stomach adenocarcinoma, STAD; liver hepatocellular carcinoma, LIHC), gynaecologic system (cervical squamous cell carcinoma and endocervical adenocarcinoma, CESC; ovarian serous cystadenoma carcinoma, OV; uterine corpus endometrial carcinoma, UCEC), breast carcinoma (BRCA), skin (skin cutaneous melanoma, SKCM) and head and neck cancer types (head and neck squamous cell carcinoma, HNSC; thyroid carcinoma, THCA). In addition to miRNA regulation, alterations in promoter DNA methylation and gene DNA copy number will likely alter the expression level of target genes and may introduce noise into the evaluation of actual miRNA regulation. Thus, the multivariate linear model can more accurately evaluate miRNA-mRNA expression associations in the presence of DNA copy number and promoter methylation aberrations that extensively influence mRNA expression under a given condition [33]. Under false discovery rate (FDR) < 0.05 for the significance of miRNA coefficients, candidate cancer context-specific miRNA-mRNA interactions were identified with a multivariate linear model (Supplementary Figure S1A and Methods). Considering miRNAs normally act as negative regulators, the expression profiles of genuinely regulating pairs were expected to be anti-correlated. As a result, we found that the majority of candidate cancer context miRNA-mRNA interactions are negatively correlated (Supplementary Figure S1B). Thus, we filtered the candidate miRNA-mRNA pairs and retained negatively correlated regulation pairs. On the other hand, several studies have recently reported that use of cross-linking and Argonaute (Ago) immunoprecipitation coupled with high-throughput sequencing (CLIP-Seq) could identify endogenous targets for miRNAs [34–37]. Indeed, we found that the cancer context miRNA-mRNA pairs have markedly more Ago CLIP-supported interactions than weakly associated pairs in all individual cancer types (Figure 1A). We further evaluated the precision of filtered cancer context miRNA-mRNA pairs with CLIP-supported sites based on low-throughput, experimentally validated miRNA–target interactions. We found that the precision of miRNA-mRNA pairs with CLIP-supported sites is average 7.5 times more than those without CLIP-supported sites (Figure 1B). For example, the proportion is increased to as high as 12.6 times in LIHC, and the minimum fold is to 4.7 in PRAD (Figure 1B). Thus, only candidate cancer context miRNA regulation with both negative correlations and CLIP-supported target sites were used for subsequent analyses (Supplementary Table S1), and all of these cancer context miRNA-mRNA interactions were assembled into the pan-cancer miRNA regulatory network. Table 1. Summary of analysed TCGA cancer data sets and cooperative miRNAs Cancer type Description Samples miRNAs mRNAs Cooperative miRNA pairs Cooperative miRNAs BLCA Bladder carcinoma 401 566 15 459 41 35 BRCA Breast carcinoma 603 506 15 455 51 45 CESC Cervical squamous cell carcinoma 292 517 15 425 3 6 COAD Colon adenocarcinoma 258 488 14 983 52 37 GBM Glioblastoma multiforme 78 474 15 090 50 50 HNSC Head and neck squamous cell carcinoma 474 520 15 447 24 25 KIRC Kidney renal clear cell carcinoma 179 429 15 180 56 45 KIRP Kidney renal papillary cell carcinoma 272 462 15 096 60 46 LGG Lower-grade glioma 509 518 15 147 222 99 LIHC Liver hepatocellular carcinoma 274 574 14 988 8 12 LUAD Lung adenocarcinoma 444 531 15 552 43 34 LUSC Lung squamous cell carcinoma 338 532 15 780 23 26 OV Ovarian cancer 293 715 15 742 53 51 PRAD Prostate adenocarcinoma 487 417 15 093 67 51 SKCM Skin cutaneous melanoma 351 597 15 239 10 15 STAD Stomach adenocarcinoma 367 499 15 820 73 46 THCA Thyroid carcinoma 498 485 14 920 37 35 UCEC Uterine corpus endometrial carcinoma 230 564 15 350 12 17 Total 6348 867a 20 472a 767 149 Cancer type Description Samples miRNAs mRNAs Cooperative miRNA pairs Cooperative miRNAs BLCA Bladder carcinoma 401 566 15 459 41 35 BRCA Breast carcinoma 603 506 15 455 51 45 CESC Cervical squamous cell carcinoma 292 517 15 425 3 6 COAD Colon adenocarcinoma 258 488 14 983 52 37 GBM Glioblastoma multiforme 78 474 15 090 50 50 HNSC Head and neck squamous cell carcinoma 474 520 15 447 24 25 KIRC Kidney renal clear cell carcinoma 179 429 15 180 56 45 KIRP Kidney renal papillary cell carcinoma 272 462 15 096 60 46 LGG Lower-grade glioma 509 518 15 147 222 99 LIHC Liver hepatocellular carcinoma 274 574 14 988 8 12 LUAD Lung adenocarcinoma 444 531 15 552 43 34 LUSC Lung squamous cell carcinoma 338 532 15 780 23 26 OV Ovarian cancer 293 715 15 742 53 51 PRAD Prostate adenocarcinoma 487 417 15 093 67 51 SKCM Skin cutaneous melanoma 351 597 15 239 10 15 STAD Stomach adenocarcinoma 367 499 15 820 73 46 THCA Thyroid carcinoma 498 485 14 920 37 35 UCEC Uterine corpus endometrial carcinoma 230 564 15 350 12 17 Total 6348 867a 20 472a 767 149 a miRNAs and mRNAs expressed in at least one cancer type. Table 1. Summary of analysed TCGA cancer data sets and cooperative miRNAs Cancer type Description Samples miRNAs mRNAs Cooperative miRNA pairs Cooperative miRNAs BLCA Bladder carcinoma 401 566 15 459 41 35 BRCA Breast carcinoma 603 506 15 455 51 45 CESC Cervical squamous cell carcinoma 292 517 15 425 3 6 COAD Colon adenocarcinoma 258 488 14 983 52 37 GBM Glioblastoma multiforme 78 474 15 090 50 50 HNSC Head and neck squamous cell carcinoma 474 520 15 447 24 25 KIRC Kidney renal clear cell carcinoma 179 429 15 180 56 45 KIRP Kidney renal papillary cell carcinoma 272 462 15 096 60 46 LGG Lower-grade glioma 509 518 15 147 222 99 LIHC Liver hepatocellular carcinoma 274 574 14 988 8 12 LUAD Lung adenocarcinoma 444 531 15 552 43 34 LUSC Lung squamous cell carcinoma 338 532 15 780 23 26 OV Ovarian cancer 293 715 15 742 53 51 PRAD Prostate adenocarcinoma 487 417 15 093 67 51 SKCM Skin cutaneous melanoma 351 597 15 239 10 15 STAD Stomach adenocarcinoma 367 499 15 820 73 46 THCA Thyroid carcinoma 498 485 14 920 37 35 UCEC Uterine corpus endometrial carcinoma 230 564 15 350 12 17 Total 6348 867a 20 472a 767 149 Cancer type Description Samples miRNAs mRNAs Cooperative miRNA pairs Cooperative miRNAs BLCA Bladder carcinoma 401 566 15 459 41 35 BRCA Breast carcinoma 603 506 15 455 51 45 CESC Cervical squamous cell carcinoma 292 517 15 425 3 6 COAD Colon adenocarcinoma 258 488 14 983 52 37 GBM Glioblastoma multiforme 78 474 15 090 50 50 HNSC Head and neck squamous cell carcinoma 474 520 15 447 24 25 KIRC Kidney renal clear cell carcinoma 179 429 15 180 56 45 KIRP Kidney renal papillary cell carcinoma 272 462 15 096 60 46 LGG Lower-grade glioma 509 518 15 147 222 99 LIHC Liver hepatocellular carcinoma 274 574 14 988 8 12 LUAD Lung adenocarcinoma 444 531 15 552 43 34 LUSC Lung squamous cell carcinoma 338 532 15 780 23 26 OV Ovarian cancer 293 715 15 742 53 51 PRAD Prostate adenocarcinoma 487 417 15 093 67 51 SKCM Skin cutaneous melanoma 351 597 15 239 10 15 STAD Stomach adenocarcinoma 367 499 15 820 73 46 THCA Thyroid carcinoma 498 485 14 920 37 35 UCEC Uterine corpus endometrial carcinoma 230 564 15 350 12 17 Total 6348 867a 20 472a 767 149 a miRNAs and mRNAs expressed in at least one cancer type. Figure 1. View largeDownload slide Cancer context miRNA–target regulation across cancer types. (A) Enrichment of Ago CLIP-supported miRNA–target interactions as a function of miRNA-mRNA expression association across cancer types. The miRNA–target interactions identified by a multivariate linear model with negative correlations were equally divided into 60 bins according to the P-values of regression. The y-axis represents the proportion of Ago CLIP-supported miRNA–target interactions occurring in each bin. The right x-axis represents the P-value of regression. (B) The left y-axis represents the number of experimentally validated miRNA–target interactions in the cancer context miRNA-mRNA interactions finally identified by us in each cancer, corresponding to the bars. The right y-axis represents the ratio between the precision of cancer context miRNA regulation with both negative correlations and CLIP-supported target sites and the precision of cancer context miRNA regulation only with negative correlations in each cancer, corresponding to the line. (C) The miRNA-mRNA interactions were divided into 11 categories according to the number of cancer types in which these interactions occur. The y-axis represents the proportion of each category of miRNA-mRNA regulation across diverse cancer types. (D) Proportion of miRNA-mRNA regulation occurring in different cancer types in different miRNA-mRNA regulation categories. (E) A pan-cancer miRNA-mRNA interaction subnetwork. A triangle node marks mRNA, and a circle node marks miRNA. An edge marks the miRNA-mRNA interaction. The miRNAs are shared by all the cancers, mRNAs occur in at least 16 cancers and miRNA-mRNA interactions occur in at least 15 cancers. Figure 1. View largeDownload slide Cancer context miRNA–target regulation across cancer types. (A) Enrichment of Ago CLIP-supported miRNA–target interactions as a function of miRNA-mRNA expression association across cancer types. The miRNA–target interactions identified by a multivariate linear model with negative correlations were equally divided into 60 bins according to the P-values of regression. The y-axis represents the proportion of Ago CLIP-supported miRNA–target interactions occurring in each bin. The right x-axis represents the P-value of regression. (B) The left y-axis represents the number of experimentally validated miRNA–target interactions in the cancer context miRNA-mRNA interactions finally identified by us in each cancer, corresponding to the bars. The right y-axis represents the ratio between the precision of cancer context miRNA regulation with both negative correlations and CLIP-supported target sites and the precision of cancer context miRNA regulation only with negative correlations in each cancer, corresponding to the line. (C) The miRNA-mRNA interactions were divided into 11 categories according to the number of cancer types in which these interactions occur. The y-axis represents the proportion of each category of miRNA-mRNA regulation across diverse cancer types. (D) Proportion of miRNA-mRNA regulation occurring in different cancer types in different miRNA-mRNA regulation categories. (E) A pan-cancer miRNA-mRNA interaction subnetwork. A triangle node marks mRNA, and a circle node marks miRNA. An edge marks the miRNA-mRNA interaction. The miRNAs are shared by all the cancers, mRNAs occur in at least 16 cancers and miRNA-mRNA interactions occur in at least 15 cancers. To have an overview of cancer context miRNA-mRNA pairs across 18 cancer types, we calculated the number of cancer types in which each miRNA-mRNA interaction occur. We found that the majority of cancer context miRNA-mRNA interactions occur in more than one cancer type (Figure 1C and Supplementary Figure S1C). Furthermore, approximately 24.9–71.0% of miRNA-mRNA interaction pairs were shared by at least five cancer types across each cancer context miRNA regulatory network (Figure 1C) and represented approximately 45.8% of pan-cancer miRNA-mRNA interactions. Conversely, the fraction of miRNA-mRNA interactions active in one cancer ranged from 2.7 to 30.5%, and most of them were present in GBM, OV, KIRC and UCEC, representing only approximately 12.5% of the pan-cancer miRNA-mRNA interaction network. These results suggest that most miRNA-mRNA interactions likely play roles in multiple cancer types, while a limited number of miRNA-mRNA interactions perform biological functions in only one specific cancer type. Moreover, miRNA-mRNA interactions were decreased, accompanied by a greater number of cancer types in which miRNA-mRNA interactions can occur. These results indicate that few core miRNA-mRNA interactions occur in most cancer types, but these miRNAs might cooperatively target genes across different cancers to perform a wide range of biological functions. We also calculated the contribution rate of each cancer in each type of miRNA regulation. As a result, we found that some cancers tend to have miRNA regulation in limited cancer types, such as GBM, LGG and OV, and there were also several cancers with the opposite tendency, for example, CESC, COAD and LUAD (Figure 1D), indicating that miRNA regulation has different important roles across cancers. Additionally, we found that most mRNAs are targeted by only a few miRNAs in cancer, as shown in Supplementary Figure S2. The number of mRNAs regulated by multiple miRNAs is small, but the proportion of these is relatively larger in LGG. Finally, we focused on miRNAs and mRNAs shared by multiple cancer types, and a pan-cancer miRNA-mRNA interaction subnetwork was constructed in which the miRNAs were shared by all the cancers, the mRNAs were present in at least 16 cancers and the miRNA-mRNA interactions occurred in at least 15 cancers (Figure 1E). MiRNA-miRNA cooperative interactions across cancer types reveal cancer clusters with similar tissue origin Next, we systematically identified miRNA-miRNA cooperation across 18 cancer types by using a multistage method [9]. Then, we assembled all identified cooperative miRNA pairs into a cooperative network in each cancer, where nodes represent miRNAs and edges represent their functional cooperative interactions (Table 1 and Supplementary Figure S3). Among the identified cooperative miRNA pairs, several have been reported in cancer, such as hsa-miR-105 in ovarian cancer, hsa-miR-205, hsa-miR-133a and hsa-miR-145 in prostate cancer, hsa-miR-9 in breast cancer [19, 20, 22]. Moreover, we constructed a miRNA-miRNA cooperative regulatory landscape in pan-cancer by integrating 18 cancer context cooperative networks (Figure 2). In total, the integrated network of pan-cancer contained 762 miRNA-miRNA cooperative interactions among 149 miRNAs (Figure 3A). Approximately 41.6% (62 of 149) of miRNAs exhibit cooperation in at least five cancers, and only 18.0% (27 of 149) of miRNAs exhibit cooperation in one cancer (Figure 3B upper panel), indicating that most miRNAs tend to cooperate with other miRNAs in several cancers. However, we discovered that for miRNA cooperation, 88.0% (674 of 762) are cancer-specific, and only 0.1% (1 of 762) of cooperative interactions are conserved in more than five cancers (Figure 3B, down panel). These results indicate that although miRNAs generally play cooperative roles in pan-cancer, they might cooperatively regulate with different miRNAs in a given cancer. The low conservation of miRNA cooperation may be explained in part by the cancer-specific expression of both miRNAs and genes. Especially, the cooperation between hsa-miR-200a-3p and hsa-miR-429 is shared by six cancers. Both of these miRNAs belong to miR-200 family, the role of which in tumour progression has been linked to several cancers [38, 39]. To obtain an overview of topological structure of this cooperative network in pan-cancer, we first examined degree distribution of this network, which revealed a power-law distribution. Therefore, similar to many large-scale networks, the cooperative network displayed scale-free characteristics (Figure 3C), indicating that the integrated cooperative network is not random but is characterized by a core set of organizing principles in its structure that distinguish it from randomly linked networks. The miRNA-miRNA synergistic network in ovarian cancer constructed by Ying et al. [19] was also a small-world network, with a few miRNAs interacting with a relatively large number of miRNA partners, and many miRNAs had few miRNA partners. Last, we found that the integrated cooperative network had higher clustering coefficients than randomly linked networks, as expected for module characteristics (Figure 3D). The scale-free and modular characteristics of the cooperative network in pan-cancer indicate that miRNA-miRNA cooperative interactions influence each other and effectively exchange regulated information both at a global and at a local scale. Figure 2. View largeDownload slide Cooperative miRNA-miRNA regulatory landscape in pan-cancer. Integrated Circos plot showing the landscape of miRNA-miRNA cooperative interactions in 18 cancer types. The outmost circle indicates the degree of cooperative miRNAs in the cooperative network in pan-cancer (the first track). The next outermost circle represents the number of cancer types in which the cooperative miRNAs occur (the second track). The next outermost circle represents the name of cooperative miRNAs (the third track). The 4th track to the 21st track represents 18 cancer types: if the miRNAs occur in this cancer type, a bar is presented that is labelled with the corresponding colour. In the centre of the figure, each arc indicates cooperative regulation between two miRNAs. The coloured arcs represent cooperative regulation occurring in a specific cancer type, whereas the black arcs represent regulation occurring in at least two cancer types. Figure 2. View largeDownload slide Cooperative miRNA-miRNA regulatory landscape in pan-cancer. Integrated Circos plot showing the landscape of miRNA-miRNA cooperative interactions in 18 cancer types. The outmost circle indicates the degree of cooperative miRNAs in the cooperative network in pan-cancer (the first track). The next outermost circle represents the number of cancer types in which the cooperative miRNAs occur (the second track). The next outermost circle represents the name of cooperative miRNAs (the third track). The 4th track to the 21st track represents 18 cancer types: if the miRNAs occur in this cancer type, a bar is presented that is labelled with the corresponding colour. In the centre of the figure, each arc indicates cooperative regulation between two miRNAs. The coloured arcs represent cooperative regulation occurring in a specific cancer type, whereas the black arcs represent regulation occurring in at least two cancer types. Figure 3. View largeDownload slide Integrated view of cooperative miRNA-miRNA regulation in pan-cancer. (A) The pan-caner miRNA-miRNA cooperative interaction network. This network consists of 762 cooperative interactions between 149 miRNAs. A circle node represents miRNA. The node size is proportional to the degree of miRNAs. The number of cancer types in which the miRNA occurs is mapped to the green gradient bar. An edge represents cooperative regulation between miRNAs. The width of the edge represents the number of cancer types in which the miRNA-miRNA pair occurs. (B) Distribution of cooperative miRNA and miRNA pairs. (Upper panel) The pie chart shows the proportion of cooperative miRNAs present in different numbers of cancer types. (Down panel) The pie chart shows the proportion of cooperative miRNA pairs present in different numbers of cancer types. (C) Distribution of degree for the observed pan-cancer miRNA-miRNA cooperative network (red circles and red line) and permuted networks (box plots). (D) The cluster coefficient of the pan-cancer miRNA-miRNA cooperative network is higher than those of the random networks (E) and (F). The Jaccard similarity coefficient matrix shows the similarity between each pair of miRNA-miRNA cooperative networks from the structure in (E). The Pearson's correlation of the cooperative miRNA pair matrix shows the similarity between each pair of miRNA-miRNA cooperative network from the expression in (F). Cancers with the same origin are labelled with grey shading and dashed boxes. (G) and (H) The subnetwork shared by COAD and STAD. A circle node represents miRNA. The black solid line presents miRNA-miRNA cooperative interactions shared by two cancers. The red line and light-blue line, respectively, represent COAD- and STAD-specific miRNA-miRNA cooperative interactions. Figure 3. View largeDownload slide Integrated view of cooperative miRNA-miRNA regulation in pan-cancer. (A) The pan-caner miRNA-miRNA cooperative interaction network. This network consists of 762 cooperative interactions between 149 miRNAs. A circle node represents miRNA. The node size is proportional to the degree of miRNAs. The number of cancer types in which the miRNA occurs is mapped to the green gradient bar. An edge represents cooperative regulation between miRNAs. The width of the edge represents the number of cancer types in which the miRNA-miRNA pair occurs. (B) Distribution of cooperative miRNA and miRNA pairs. (Upper panel) The pie chart shows the proportion of cooperative miRNAs present in different numbers of cancer types. (Down panel) The pie chart shows the proportion of cooperative miRNA pairs present in different numbers of cancer types. (C) Distribution of degree for the observed pan-cancer miRNA-miRNA cooperative network (red circles and red line) and permuted networks (box plots). (D) The cluster coefficient of the pan-cancer miRNA-miRNA cooperative network is higher than those of the random networks (E) and (F). The Jaccard similarity coefficient matrix shows the similarity between each pair of miRNA-miRNA cooperative networks from the structure in (E). The Pearson's correlation of the cooperative miRNA pair matrix shows the similarity between each pair of miRNA-miRNA cooperative network from the expression in (F). Cancers with the same origin are labelled with grey shading and dashed boxes. (G) and (H) The subnetwork shared by COAD and STAD. A circle node represents miRNA. The black solid line presents miRNA-miRNA cooperative interactions shared by two cancers. The red line and light-blue line, respectively, represent COAD- and STAD-specific miRNA-miRNA cooperative interactions. Several lines of evidence have indicated that cancer types with similar tissue origins share multiple molecular features, such as gene expression, at the level of both protein-coding genes and miRNAs [40–42]. However, whether cancer types with similar tissue origins exhibit similar miRNA cooperation patterns is still unknown. To address this question, we calculated a paired similarity score (Jaccard coefficient) based on miRNA cooperation in each cancer. We found that adenocarcinomas, such as COAD, STAD, LUAD and PRAD, showed greater similarity to each other (Figure 3E). This analysis indicated that cancer types with similar tissue origins show similar miRNA cooperation patterns. Especially, COAD and STAD shared 24 miRNAs, and eight of these cooperation patterns were common for the two cancers (Figure 3G and H). We further discovered that these miRNA cooperation patterns are closely related to carcinogenesis by regulating kappa B kinase NF-kappa B cascade, cytoskeleton organization and biogenesis. Previous studies have shown that hsa-miR-106a-5p, hsa-miR-30c-5p and hsa-miR-142-3p are involved in both COAD and STAD, and here, we found they cooperatively regulate cytoplasm organization and biogenesis processes and intercellular junction assembly and maintenance processes. These observations suggest that related miRNA-miRNA cooperative regulation mechanisms might operate in cancer types with similar tissue origins. Moreover, we also measured the similarity between diverse cancers based on the expression correlation coefficient of cooperative miRNAs and found that the expression pattern of pan-cancer cooperative miRNAs also revealed tissue-of-origin cancer. After clustering all cancers, we found that three squamous cell carcinomas, HNSC, LUSC and CESC, are closely clustered together (Figure 3F and Supplementary Figure S4). Another example is the cluster comprising KIRC and KICH, which are two types of kidney cancer. We found that their similarity in cooperative miRNA expression was higher than with other cancers. Thus, we concluded that both the structure of miRNA cooperative networks and the expression of cooperative miRNA pairs could reveal similar tissue origin of known and novel cancer clusters at different levels. Common regulatory principles of cooperative miRNA-miRNA pairs Through systematic analysis of the expression and function of cooperative regulatory miRNA-miRNA networks across cancers, some common regulatory principles of the cooperative regulatory miRNA-miRNA networks were summarized. First, we compared the expression levels of cooperative miRNAs with those of other miRNAs in each cancer, and the results show that the expression of cooperative miRNAs is higher than that of other miRNAs (Figure 4A). Moreover, the expression of cooperative miRNAs in multiple cancers tended to be higher than that of those involved in few cancers (Supplementary Figure S5). At the same time, the expression fluctuation of these cooperative miRNAs is significantly greater than that of non-cooperative miRNAs in each cancer (Figure 4B). Especially, cooperative miRNAs have the most expression fluctuation in LGG, and we proposed that these cooperative miRNAs might be important miRNAs in LGG [43]. We found that cooperative miRNA pairs exhibit stronger expression correlation than random miRNA pairs in each cancer (Figure 4C). Moreover, we also calculated the average regulatory strength of cooperative miRNAs towards their target genes and found significantly stronger regulation than non-cooperative miRNAs in 14 of the cancer types (Figure 4D). These results suggest that similar expression patterns might help miRNAs to perform cooperative functions. Figure 4. View largeDownload slide Common regulatory principles of cooperative regulatory miRNA-miRNA pairs. (A) The expression of cooperative miRNAs across diverse cancer types. The grey box represents the expression of miRNAs randomly selected. (B) Radar plot representing the expression fluctuation of cooperative miRNAs. The axes report the ratio of the number of cooperative miRNAs with smaller (or larger) variance to the number of total miRNAs with smaller (or larger) variance. The blue line represents cooperative miRNAs with larger variance. The red line represents cooperative miRNAs with smaller variance. (C) The correlation of expression of cooperative miRNA pairs across diverse cancer types. The grey box represents the expression correlation between non-cooperative miRNA pairs. The Pearson’s correlation coefficient was used. (D) The regulatory strength of cooperative miRNAs towards target genes across diverse cancer types. The grey box represents non-cooperative miRNA. The y-axis represents the mean of the regression coefficient (absolute value) of miRNAs in the multivariate linear model. (E) The cancer hallmark-related miRNA-miRNA cooperative subnetwork, which was assembled with the miRNA pairs that co-regulated at least one cancer hallmark-related GO term. (F) A heat-map containing the number of cooperative miRNAs that regulate GO terms related to cancer hallmarks across diverse cancer types. Figure 4. View largeDownload slide Common regulatory principles of cooperative regulatory miRNA-miRNA pairs. (A) The expression of cooperative miRNAs across diverse cancer types. The grey box represents the expression of miRNAs randomly selected. (B) Radar plot representing the expression fluctuation of cooperative miRNAs. The axes report the ratio of the number of cooperative miRNAs with smaller (or larger) variance to the number of total miRNAs with smaller (or larger) variance. The blue line represents cooperative miRNAs with larger variance. The red line represents cooperative miRNAs with smaller variance. (C) The correlation of expression of cooperative miRNA pairs across diverse cancer types. The grey box represents the expression correlation between non-cooperative miRNA pairs. The Pearson’s correlation coefficient was used. (D) The regulatory strength of cooperative miRNAs towards target genes across diverse cancer types. The grey box represents non-cooperative miRNA. The y-axis represents the mean of the regression coefficient (absolute value) of miRNAs in the multivariate linear model. (E) The cancer hallmark-related miRNA-miRNA cooperative subnetwork, which was assembled with the miRNA pairs that co-regulated at least one cancer hallmark-related GO term. (F) A heat-map containing the number of cooperative miRNAs that regulate GO terms related to cancer hallmarks across diverse cancer types. The miRNA gene family is an important component of small RNAs, and we further investigated the cooperation between homologous miRNAs. We found that the proportion of miRNAs that could cooperate with their family members is different cross-cancer types (Supplementary Figure S6). The proportion of miRNAs that could cooperate with their family members is up to 60% in BLCA. For example, hsa-miR-200a-3p and hsa-miR-429 belong to the mir-8 family and cooperatively regulate the insulin receptor signalling pathway. In addition, it is well known that cancer miRNAs can be onco-miRNAs or suppressive miRNAs. We further explored cooperative interactions among oncogenic miRNAs and suppressive miRNAs [44]. We found that miRNAs belonging to the same category tend to cooperatively interact in several cancer types, such as BRCA, COAD, GBM and PRAD (Supplementary Figure S7). For example, hsa-miR-21-5p and hsa-miR-27a-3p, which are both oncogenic miRNAs, cooperatively interact in BRCA. The complexity of cancer can be reduced and represented by a few cancer hallmarks that enable tumour growth and metastasis [45]. These hallmarks provide a framework for understanding the remarkable diversity of cancers. To further investigate the functional roles of cooperative miRNA pairs in cancer hallmarks, we selected the miRNA pairs that co-regulated at least one cancer hallmark-related Gene Ontology (GO) term. These miRNA pairs were assembled into a cancer hallmark-related subnetwork (Figure 4E and Supplementary Table S2). We found that several cancer hallmarks are deeply regulated by the cooperative miRNA subnetwork. For example, 22 cooperative miRNAs in STAD are involved in regulating cancer hallmarks, and eight of these regulated tissue invasion and metastasis, while in LGG, 12 of 16 miRNAs regulated insensitivity to antigrowth signals (Figure 4F). Moreover, miRNA cooperation regulated the same cancer hallmarks in different cancers; for instance, 26 cooperative miRNAs negatively regulated both apoptosis and programmed cell death, evading apoptosis in eight types of cancer. In addition to cancer hallmarks reflecting the common characteristics of cancer, cooperative miRNAs can also control specific signalling pathways (Supplementary Figure S8). Two pairs of miRNAs, hsa-miR-301a-3p and hsa-miR-130a-3p pair, and hsa-miR-200b-3p and hsa-miR-429 pair, regulate the epidermal growth factor receptor signalling pathway in OV. These results indicate cancer specificities for cooperative miRNA-mediated regulation of cancer hallmarks and signalling pathways. Moreover, same miRNA pairs might regulate different hallmarks in different cancers. For example, hsa-miR-106a-5p and hsa-miR-340-5p cooperatively regulate cytoplasm organization and biogenesis in STAD but positively regulate transcription signal transduction in KIRP. Thus, we conclude that cooperative miRNA pairs may contribute to carcinogenesis by regulating cancer hallmark functions in multiple cancers as well as in cancer-specific patterns. Conserved and rewired cooperative network miRNA hubs Hub nodes whose connectivity is extremely high always play important roles in biological network. For example, in the protein–protein interaction networks of various organisms, hubs tend to be essential proteins [46]. To systematically investigate the key nodes in the miRNA-miRNA cooperative network, we first identified the top 30.0% of nodes with the highest connectivity within the pan-cancer cooperative network as hub miRNAs. Next, we assessed the extent to which these hub miRNAs are shared among the different cancer-specific miRNA-miRNA cooperative networks. These hubs were grouped into three categories: common hubs, cancer-specific hubs and other hubs. The common hubs exhibit central roles of miRNAs across more than one cancer; the specific hubs identify miRNAs with specific central roles in a given cancer. For example, hsa-miR-301a-3p is involved in 14 cancer types, is a hub in the pan-cancer network and is also a hub in HNSC, OV and KIRP. Thus, hsa-miR-301a-3p was defined as a common hub based on our method. hsa-miR-186-5p is also a common hub that participates in 17 cancer cooperative networks and is a hub miRNA in 2 (COAD, LUAD) cancer cooperative networks (Figure 5A). Common hub hsa-miR-186-5p has 29 cooperative miRNA partners in the pan-cancer cooperative network. In addition, 36.7% (18 of 49) of hubs are specific hubs, and the remaining 14.3% are other hubs (Figure 5A). For example, hsa-miR-107 is a specific hub in OV. To further investigate the roles of these three types of hubs, we associated each miRNA hub with biological function based on the function modules they cooperatively regulate. The results showed that the three categories of hubs regulated many common functions, such as regulating apoptosis, regulation of programmed cell death and some signalling pathways (Figure 5B and Supplementary Table S3). At the same time, TAM was also used to perform functional annotation for these hub miRNAs directly [47]. The annotation results showed that hub miRNA regulate many common cancer-related functions, such as cell proliferation, DNA repair, apoptosis and immune response. These results indicated that these miRNAs might play important roles in cancer. On the other hand, some functions are specifically controlled by different categories of hubs, such as angiogenesis and vasculature development, which are specifically regulated by the hubs hsa-miR-193b-3p and hsa-miR-29b-3p (Figure 5B). Thus, cancer-specific miRNA hubs may contribute to carcinogenesis by regulating cancer-specific functions. Figure 5. View largeDownload slide The conserved and rewired network hubs in each cancer type. (A) A summary bubble-bar plot showing the distribution of hub miRNAs in each cancer-specific miRNA-miRNA cooperative network. The top bars show the number of hub miRNAs (blue bars) and experimentally validated miRNAs related to cancer (grey bars) occurring in each cancer miRNA-miRNA cooperative network. In addition, the bars on the right show the degree of hub miRNAs in the pan-cancer miRNA-miRNA cooperative network. The bubble size indicates the degree of hub miRNAs standardized by the total number of miRNAs in each cancer miRNA-miRNA cooperative network, and different colours correspond to different hub categories. (B) Selected GO functional categories regulated by common hubs, cancer-specific hubs and other hubs. The bar on the right shows GO functional categories regulated by at least two types of hubs, and the bar on the left shows GO functional categories regulated by only one type of hub. The length of the bars represents the number of hub miRNAs. (C) Comparison of the number of anticancer drug targets in three types of hub miRNA targets. (D) An example of a common hub miRNA, hsa-miR-186-5p. (E) Clustering trees based on the expression of hub miRNAs that are differently expressed between cancer and normal tissues in BLCA, LUAD and UCEC. Dark cyan represents the normal samples, and red represents the tumour samples. (F) The ROCs of the hub miRNAs differentially expressed between cancer and normal tissues in BLCA, LUAD and UCEC. Figure 5. View largeDownload slide The conserved and rewired network hubs in each cancer type. (A) A summary bubble-bar plot showing the distribution of hub miRNAs in each cancer-specific miRNA-miRNA cooperative network. The top bars show the number of hub miRNAs (blue bars) and experimentally validated miRNAs related to cancer (grey bars) occurring in each cancer miRNA-miRNA cooperative network. In addition, the bars on the right show the degree of hub miRNAs in the pan-cancer miRNA-miRNA cooperative network. The bubble size indicates the degree of hub miRNAs standardized by the total number of miRNAs in each cancer miRNA-miRNA cooperative network, and different colours correspond to different hub categories. (B) Selected GO functional categories regulated by common hubs, cancer-specific hubs and other hubs. The bar on the right shows GO functional categories regulated by at least two types of hubs, and the bar on the left shows GO functional categories regulated by only one type of hub. The length of the bars represents the number of hub miRNAs. (C) Comparison of the number of anticancer drug targets in three types of hub miRNA targets. (D) An example of a common hub miRNA, hsa-miR-186-5p. (E) Clustering trees based on the expression of hub miRNAs that are differently expressed between cancer and normal tissues in BLCA, LUAD and UCEC. Dark cyan represents the normal samples, and red represents the tumour samples. (F) The ROCs of the hub miRNAs differentially expressed between cancer and normal tissues in BLCA, LUAD and UCEC. We further systematically investigated the association between miRNA cooperation and anticancer drugs. We found that 72 of 105 FDA-approved anticancer drugs could target hub miRNA–target genes. The common hub miRNAs contained the largest number of anticancer drug target genes compared with specific hubs and other hubs (Figure 5C). The results suggest that targets of common hub miRNAs are more likely to be drug-able. For example, some target genes of the common miRNA hsa-miR-186-5p, including FGFR1 and GART, have been used as anticancer drugs targets of Sorafenib in KIRC/KIRP and Pemetrexed in LUAD/LUSC (Figure 5D). According to the pan-cancer characteristics of hsa-miR-186-5p, target genes of this miRNA are favoured to be drugable and may be potential targets of anticancer drugs in the other six cancer types. In addition, we also found that 31 of 49 hub miRNAs could regulate cancer hallmarks. The common hub miRNAs regulate more hallmarks than the specific hub and other hub miRNAs (Supplementary Table S4 and Supplementary Figure S9A). For example, hsa-miR-106a-5p regulates self-sufficiency in growth signals and tissue invasion and metastasis. We found that 49.0% (24 of 49) of hub miRNAs are experimentally validated cancer-related miRNAs, including 12 common hubs, 9 specific hubs and 3 other hubs (Supplementary Table S4 and Supplementary Figure S9B). For example, the common hub hsa-miR-137 has been be associated with multiple cancer types, including GBM, BLCA, COAD, STAD, LIHC, OV, BRCA and SKCM. Finally, we further investigated the expression of hub miRNAs. We also found that the average expression level of other hub miRNAs is higher than that of specific hubs and common hubs (Supplementary Table S4 and Supplementary Figure S9C). These results suggested that hub miRNAs might play fundamental roles in multiple cancer types. It is widely accepted that the occurrence of cancer is related to the concerted activity of many miRNAs [48]. Functionally coherent miRNA cooperative modules can be used as biomarkers to distinguish cancer from normal tissues. Here, we evaluated the efficiency of the hubs using the area under the receiver operating characteristic (ROC) curve (AUC). ROC curves of the hubs exhibit an AUC value >0.9 in BLCA, LUAD and UCEC (Figure 5E and F). MiRNA cooperative modules can be used as prognostic biomarkers in multiple cancer types The notion that it may be possible to reduce cancer mortality by identifying and monitoring survival-related biomarkers rests on the idea that a module biomarker is a better predictor of survival than an individual gene [49]. Given that integrated cooperative network is associated with cancer-related biological processes and has the characteristics of pan-cancer, as shown in the above analysis results, we were interested in identifying prognostic modules from the pan-cancer cooperative network. Modules in the miRNA cooperative network represent groups of functionally related miRNAs dedicated to specific biological processes. We therefore used CFinder software with suggested parameters to identify modules in the pan-cancer cooperative network [50]. In total, with K = 6, two miRNA cooperative modules were identified (Figure 6A and Supplementary Figure S10). One module comprised seven common hubs, and the other module included five common hubs and two specific hubs. For instance, Figure 6A shows a module that includes 20 cooperative pairs among seven common hub miRNAs, each of which are involved in nine types of cancer miRNA cooperative networks on average. Functional analysis of the module revealed that these miRNAs are involved in multiple cancer-associated hallmarks, such as evading apoptosis, self-sufficiency in growth signals and insensitivity to antigrowth signals (Figure 6A). We further found that all the miRNAs in this module are well-known cancer-related miRNAs. The results demonstrate that this module is closely related to cancer development and represents pan-cancer characteristics. Figure 6. View largeDownload slide The miRNA-miRNA cooperative module as a biomarker. (A) An example of a miRNA-miRNA cooperative module. The Circos plot shows the characteristic of the miRNAs in the module. In the centre of the figure, the miRNA-miRNA cooperative module is shown. In this module, a circle node represents miRNA; the node size is proportional to the degree of miRNAs; the number of cancer types in which the miRNA occurs is mapped to the green gradient bar (details in Figure 3A); an edge represents the cooperative regulation between miRNAs. The circle is equally divided into seven parts by the seven miRNAs. Each sector corresponds to one miRNA. The innermost circle indicates cancer hallmarks regulated by the miRNAs (the first track). The next innermost circle represents the association between miRNAs and cancer. If it has been indicated that the miRNA is related to the cancer, a bar is presented that is labelled with the corresponding colour (the second track). The next circle represents the association between miRNAs and cancer-specific cooperative miRNA-miRNA networks. If the miRNA occurred in the cooperative miRNA-miRNA network of this cancer type, a bar is presented that is labelled with the corresponding colour (the third track). The next circle represents the degree of miRNAs in the corresponding cancer-specific cooperative miRNA-miRNA networks, which is presented as a bar plot (the fourth track). (B) The expression pattern of miRNAs from the module that were significant in single-variable Cox regression and their performance in survival analysis based on five types of cancer patients in training set and test set. (Upper panel) Colour-gram of miRNA expression profiles of patients ordered by risk score (details in the ‘Materials and methods’ section). The median of the risk score cut-off was used to divide patients into low-risk (grey line) and high-risk groups (brick-red line). (Down panel) Kaplan–Meier estimates of overall survival of the two patient subgroups identified by the miRNA signature in training set and test set. Figure 6. View largeDownload slide The miRNA-miRNA cooperative module as a biomarker. (A) An example of a miRNA-miRNA cooperative module. The Circos plot shows the characteristic of the miRNAs in the module. In the centre of the figure, the miRNA-miRNA cooperative module is shown. In this module, a circle node represents miRNA; the node size is proportional to the degree of miRNAs; the number of cancer types in which the miRNA occurs is mapped to the green gradient bar (details in Figure 3A); an edge represents the cooperative regulation between miRNAs. The circle is equally divided into seven parts by the seven miRNAs. Each sector corresponds to one miRNA. The innermost circle indicates cancer hallmarks regulated by the miRNAs (the first track). The next innermost circle represents the association between miRNAs and cancer. If it has been indicated that the miRNA is related to the cancer, a bar is presented that is labelled with the corresponding colour (the second track). The next circle represents the association between miRNAs and cancer-specific cooperative miRNA-miRNA networks. If the miRNA occurred in the cooperative miRNA-miRNA network of this cancer type, a bar is presented that is labelled with the corresponding colour (the third track). The next circle represents the degree of miRNAs in the corresponding cancer-specific cooperative miRNA-miRNA networks, which is presented as a bar plot (the fourth track). (B) The expression pattern of miRNAs from the module that were significant in single-variable Cox regression and their performance in survival analysis based on five types of cancer patients in training set and test set. (Upper panel) Colour-gram of miRNA expression profiles of patients ordered by risk score (details in the ‘Materials and methods’ section). The median of the risk score cut-off was used to divide patients into low-risk (grey line) and high-risk groups (brick-red line). (Down panel) Kaplan–Meier estimates of overall survival of the two patient subgroups identified by the miRNA signature in training set and test set. Subsequently, we examined whether the modules are associated with the survival of cancer patients. In each cancer, we randomly assigned 70% of the specimens to a training set and 30% of the specimens to a test set; we only used the training data set for detecting the miRNA signature of cancer. We found that miRNAs in this module can be used to classify cancer samples into two groups with significantly different overall survival rates in five cancers (log-rank test, P < 0.05) (Figure 6B). In addition, we also validated the miRNA signature in the test set where samples were classified into the high-risk group or low-risk group using the same miRNAs used in the training set. As an example, four miRNAs with six cooperative interactions in this module were identified as being associated with LGG patient survival. All the three miRNAs, hsa-miR-142-3p, hsa-miR-200a-3p and hsa-miR-30e-5p, are LGG-risk miRNAs for which increased expression is associated with an increased risk of survival. Using the median risk value scores of these miRNAs as the threshold, the LGG patients were divided into a high-risk group and a low-risk group, and these two groups exhibited significantly different survival times in training set and test set (Figure 6B, P = 0.002 and P = 0.009). Three of these miRNAs, hsa-miR-142-3p, hsa-miR-30e-5p and hsa-miR-200a-3p, occur in the LGG miRNA-miRNA cooperation network, the degree of which, respectively, is 9, 5 and 5. The three miRNAs are all known cancer-associated miRNAs. Thus, this suggests that a cooperative network-based approach could effectively identify prognostic biomarkers. Discussion In the current study, we identified miRNA-miRNA cooperation in each cancer via the multiple-step method. Systematic construction and analysis of miRNA-miRNA cooperation networks across multiple cancers can help to elucidate the commonalities and differences in cancer mechanisms. As general biological networks, the pan-cancer miRNA-miRNA cooperation network is scale-free and has a small-world property. Watts and Strogatz [51] have analysed how fast disturbances spread through small-world networks and revealed that the time wasted for spreading of a disturbance in a small-world network is close to the theoretically possible minimum for any graph with the same number of nodes and edges. Therefore, small-world may allow the cooperation of miRNAs to respond quickly to disturbances. The increased number of molecular tumour maps presents a great opportunity for us to investigate molecular characteristic across various cancer types. Different molecular characteristics have shown different conservation levels and specificity among cancers. In our study, through comparison of miRNA cooperation networks across cancers, we showed that only a small proportion of the cooperative activities are shared by multiple cancers. miRNA cooperation varied greatly from one cancer type to another. Park et al. [52] found that >90% of the interactions between cancer driver alterations were only detected in a single cancer type. Additionally, Xu et al. [53] found that ∼34.56% of competitive endogenous RNA (ceRNA) regulation occurred in only one cancer and only 3.78% of ceRNA-ceRNA interactions were conserved in >10 cancers. Together, these studies illustrate how the interactions between various molecules change across different types of cancer. However, some conservation of molecules among different tumour types has been found. Hamilton et al. [54] found that 45% of broadly conserved pan-cancer oncomiRs share strong homology in their seed motifs. Davoli et al. [55] identified distinct sets of drivers in each tumour type, but the majority of these were also identified in the pan-cancer analysis as lower confidence candidates, and there is still significant overlap among different tumour types. This suggested that molecules may be conserved, and these molecule interactions might be different in different types of cancer. In 2017, the US Food and Drug Administration has approved for the first time a cancer treatment based on common molecular biomarkers rather than the location in the body where the tumour originated. This is important news for the cancer community. Molecular markers may be able to more accurately reflect the curative effect of drugs than traditional tissue origins. Heim et al. [56] found that tumours at a particular anatomic site are genetically more similar to tumours from different organs and tissues than to tumours of the same origin. Snyder and colleagues [57] discovered a neoantigen signature that was specifically present in tumours with a sustained response to CTLA-4 blockade and that examining the exomes of melanoma patients could inform the response to the anti–CTLA-4 blockade therapy. However, Hoadley et al. performed an integrative analysis using molecular profiles of 12 cancer types, revealing a unified classification into 11 major subtypes. Five subtypes were nearly identical to their tissue-of-origin counterparts, but several distinct cancer types were found to converge into common subtypes [40]. In our study, we found that some cancer types with similar tissue origins exhibited higher cooperative network structural similarity and expression of cooperative miRNA pairs and that some cancer types with the same tissue origin are different in miRNA cooperative regulation. Cancer types with similar tissue origins do not necessarily have the same molecular characteristics. Our findings provide independent evidence for development of drug and treatment strategies. Moreover, most cooperative miRNA pairs also exhibited the same expression tendency, which allows for a rapid response to disturbances, and we discovered that cooperative miRNAs have more strong regulation to the target genes in function modules. Because miRNA regulation mostly reduces target mRNA expression, miRNA–target genes with low mRNA expression levels would more likely be under strong regulation by miRNAs [58]. In our multivariate linear model, when the absolute value of the regression coefficient of miRNA was larger, the gene expression targeted by this miRNA was lower. Therefore, when the absolute value of the regression coefficient of miRNA is larger, we think that the miRNA exhibits a stronger regulation effect on the target genes. We discovered that cooperative miRNAs exhibit stronger regulation of the target genes in function modules by comparing regression coefficients. In addition, Mendell [59] has proposed that the functions of a given miRNA could be attributed to strong regulation of a select few dominant targets or, alternatively, more subtle regulation of many targets simultaneously. The target genes in function modules may be as dominant targets. Our analysis also provides a hub-based view to elucidate the common and specific miRNA hubs across cancers. Hubs are topologically centred in the miRNA cooperative network and have maximal informational links with other miRNAs. Through in-depth analyses of the structure of the pan-cancer miRNA cooperative network, we identified the cancer-common and specific hubs. Hub miRNAs can be used as potential biomarkers to distinguish cancer from normal tissues. It has been shown that serum miRNA levels are stable, reproducible and consistent in humans [60]. Meanwhile, miRNAs can be released from cancer cells to body fluids through secretion of exosomes. Owing to the surprising stability of miRNAs in tissues and body fluids, miRNAs have emerged as advantageous non-invasive cancer biomarkers with immeasurable clinical potential. For example, in early gastric cancer, almost all of the serum-based traditional biomarkers, including carcinoembryonic antigen, are not widely recognized in the clinical screening or diagnosis because of their limited specificity and sensitivity [61]. The stability of miRNAs in serum and gastric juice has been confirmed, and the expression levels of certain miRNA diagnostic biomarkers can assist in screening for gastric cancer [62]. Several miRNAs circulating in the blood of gastric cancer patients can be applied as diagnostic biomarkers, including hsa-miR-21, hsa-miR-106a/b, hsa-miR-421, hsa-miR-129, hsa-miR-199a-3p, hsa-miR-218 and hsa-miR-223 [63–67]. Some common hub miRNAs have been identified as diagnostic biomarkers in several cancers. hsa-miR-200a is a common hub miRNA, and Yun et al. [68] found that hsa-miR-200a is an independent predictor of non-muscle-invasive bladder cancer recurrence. Park et al. [69] found that the detection of hsa-miR-125a and hsa-miR-200a in saliva can be used as non-invasive and rapid diagnostic biomarkers of oral cancer. The hub miRNA hsa-miR-34a has been identified as a novel biomarker of relapse in surgically resected non-small cell lung cancer [70]. Additionally, both of the common hub miRNAs hsa-miR-128 and hsa-miR-106a-5p have been used as diagnostic markers in GBM and lung squamous cell carcinoma. Wang et al. [71] demonstrated that cell-free hsa-miR-21, hsa-miR-128 and hsa-miR-342-3p in plasma are specific and sensitive markers for GBM diagnosis, suggesting that these miRNAs may be used as non-invasive biomarkers of GBM. Zhang et al. [72] identified a three-miRNA panel comprising hsa-miR-106a-5p, hsa-miR-20a-5p and hsa-miR-93-5p as a potential diagnostic marker for male lung squamous cell carcinoma (SCC) patients. In summary, we presented the miRNA-miRNA cooperation landscape across human major cancers and showed the importance of various aspects. A bird’s-eye view of the functional miRNA cooperative networks of large sample sets encompassing multiple tumour lineages are important for understanding the key roles of miRNAs in disease. Our study opens new avenues for leveraging publicly available genomic data to investigate the functions and mechanisms of miRNAs across human cancers. Materials and methods MiRNA–target interactions The CLIP-supported miRNA–target pairs were downloaded from starBase V2.0 [73]. The target sites from at least one program (TargetScan [74] (Release 6.2), miRanda [75], Pictar2.0 [76], PITA [77] or RNA22 [78]) that resided within any entry of the Ago CLIP peaks were considered CLIP-supported sites. A miRNA–target pair was considered as the CLIP-supported pair if and only if the miRNA had at least one CLIP-supported site in the target. After converting the miRNA ID, 417 469 miRNA–target interactions between 383 miRNAs and 13 797 genes were finally reserved for subsequent analysis. High-quality experimentally validated miRNA–target interactions were obtained with R package multiMiR [79] from miRecords [80] miRTarbase [81] and TarBase [82], including 49 505 miRNA–target interactions between 13 523 genes and 459 miRNAs. Omics data across cancer types To identify miRNA–target regulation in each cancer type, we obtained paired miRNA-mRNA expression profiles as well as copy number and DNA methylation data from the TCGA Data Portal (2015). Based on the availability of these four types of molecular measurements and the sample size, 18 cancer types were analysed. Overall, 95.0% of the cancer types had >100 samples with the four abovementioned molecular measurements obtained at the same time, and thus, these relative sufficient data were used for further identification of miRNA–target regulation in the context of each cancer (Table 1). All miRNA expression data sets were obtained from TCGA Data Portal. MiRNA expression was profiled using microarrays in GBM and OVA and by small RNA sequencing in the remaining cancer types. For microarray data sets, TCGA Level 3 expression data were directly used. For miRNA sequencing data sets, expression values of mature miRNAs were used, which for a given miRNA was calculated as the ratio of reads mapping to it relative to the total number of reads mapping to annotated mature miRNAs in the given sample. We further removed miRNAs with a read count <10 that were detected in >95.0% of samples in each cancer type data set. The microarray and sequencing data expression values were log2 transformed for subsequent analysis. All mRNA expression data sets were obtained from TCGA Data Portal and profiled via RNA sequencing. For mRNA expression, mapped and gene-level-summarized (Level 3, TPM) RNA-seq data sets were used. We further removed mRNAs with a read count <20 that were detected in >95.0% of samples in each cancer type data set. The expression values were log2 transformed for subsequent analysis. DNA methylation data sets were also obtained from TCGA Data Portal (Level 3), and the methylation values for each protein-coding gene were defined as the average beta-values of probes mapping to the corresponding gene promoter (±2 kb of annotated transcription start sites). In addition, DNA copy number data sets were obtained from Firehose (http://gdac.broadinstitute.org/runs/analyses__2015_04_02/). We used Level 4 non-discretized gene-summarized log2-transformed aCGH copy number calls (tumour/normal ratio) computed by the Gistic2 algorithm [83]. Functional annotation of target genes The Biological Process (BP) terms for GO were downloaded from the MSigDB (v5.1) database [84]. As in previous studies, process categories from GO were restricted to BP terms such that the number of genes annotated to a term was at least 5 and no >500. Ultimately, 792 filtered GO BP terms were used for further analysis. Moreover, a list of GO terms determined to be related to the hallmarks of cancer was obtained from a previous study [85]. Protein–protein interaction network of target genes We assembled protein–protein interaction data from HPRD [86] and further removed self-loop interactions. Then, the gene symbols for each interaction were mapped to their corresponding Entrez gene identifiers. Finally, the maximum component of the whole protein–protein interaction network was extracted, which contained 35 865 interactions among 9028 genes. Collection of cancer-associated miRNAs Several database systems have chosen to provide a comprehensive resource indicating miRNA dysregulation in various human diseases. We manually collected the cancer-associated miRNAs from four databases, including HMDD [87], miR2Disease [88], miRCancer [89] and OncomiRDB [44] (Supplementary Table S5). Overview of the identification of miRNA cooperation across cancer types A multiple-step model was proposed to construct miRNA-miRNA functional cooperative networks in each cancer type. First, for each miRNA, a multivariate linear regression model was used to identify cancer context miRNA–target regulation (Supplementary Figure S1A). Then, the functionally cooperative miRNA pairs were identified as follows: for each miRNA pair, we initially identified miRNA pairs that significantly shared target genes; the shared target genes of the miRNA pair were used to identify candidate functional modules; and then, the candidate module sets were further filtered using two topological features in the protein interaction network (Supplementary Figure S1D). Here, we defined a pair of miRNAs as cooperative if they significantly co-regulated at least one functional module. Finally, all cooperative miRNA pairs were assembled into a functional miRNA cooperative network in each cancer type, where nodes represent miRNAs and edges represent their functional cooperative interactions. Identification of miRNA–target regulation in each cancer type For identifying cancer-specific miRNA–target regulation, we first measured the association between miRNA and mRNA expression for each miRNA-mRNA pair across a set of tumours in individual cancer types using a multivariate linear model that also took into account variation in mRNA expression induced by changes in DNA copy number and promoter methylation at the mRNA gene (Supplementary Figure S1A). For each predicted miRNA-mRNA regulation of miRNA u and mRNA i in a specific cancer, mRNA i expression y changes as a linear function of miRNA u expression x, DNA copy number x and DNA methylation x across tumour samples of a given cancer type: yi=βuxu,i+βcnxcn,i+βmexme,i+ β0+εi, i=1,…, n. In this linear model, β, βcn and β are the corresponding regression coefficients for miRNA expression, copy number and methylation variables, respectively, and β is the intercept. All P-values of the miRNA coefficient β in each cancer were subject to FDR correction. Under FDR < 0.05 for the signification of miRNA coefficients, the candidate cancer context miRNA-mRNA interactions were identified. Last, only candidate miRNA–target pairs with both negative correlation and CLIP-supported target sites were considered as cancer context miRNA regulation. Identification of functionally cooperative miRNA pairs across 18 cancer types Next, the functionally cooperative miRNA pairs were identified according to the protocol described in our previous study [9], which briefly included three steps: first, we initially identified miRNA pairs that share target mRNAs. Second, the shared target mRNAs of the miRNA pair were used to identify candidate functional modules. Third, the candidate module sets were further filtered using two topological features in the protein interaction network. A pair of miRNAs was defined as cooperative only if they significantly co-regulated at least one functional module. For a given pair of miRNAs, A and B, we first obtained their co-regulated target subset (A∩B), which was required to have at least three genes. Then, target subset functional enrichment was performed using the hypergeometric test across each selected GO term. At a given significance level, we not only obtained the enriched GO terms but also captured the subset in A∩B annotated to each term, GAB, the set of candidate function modules. We next further filtered these candidate function modules according to two topological features in the protein interaction network: (i) the minimum distance from every gene to others in the subset is no >2; (ii) the characteristic path length (CPL) is significantly shorter than random. The P-value indicating significance was calculated using the edge-switching method and was defined as the fraction of the CPL for the same subset in random networks that was shorter than that in the real network. We generated 1000 random networks using the software Mfinder (available at http://www.weizmann.ac.il/mcb/UriAlon/). After performing function enrichment and two topological restrictions in the network, miRNAs A and B were considered cooperative if they co-regulated at least one functional module. Topological measurements of the functional miRNA cooperative network In this study, we investigated several topological features of the miRNA-miRNA cooperative network using the R package ‘igraph’. First, we evaluated whether the degree distribution of the cooperative network satisfied a power-law model. For each miRNA in a miRNA-miRNA cooperative network, degree was defined as the number of edges linked to it. Clustering coefficient was defined as the ratio between the number of edges linking adjacent nodes and the total number of possible edges among them [90]. Classification of miRNA hubs Furthermore, we selected the top 30.0% of miRNAs with the highest degrees in the pan-cancer cooperative network as the hub miRNAs. To systematically analyse the hubs across cancers, we split the hubs into three groups: (i) cancer-specific hubs, in which the miRNAs were only hubs in the cooperative network of only one cancer type; (ii) common miRNA hubs, in which the miRNAs were hubs in more than one cancer type; and (iii) others. To define the first category, miRNAs belonging to hubs in the pan-cancer network that were ranked in the top 15.0% in one network but were not ranked in the top 15.0% in any other cancers were selected. Conversely, common hubs were required to be hubs in at least one cancer. The remaining miRNA hubs belonged to others. Identification of hub miRNA biomarkers For each cancer, we randomly selected tumour samples with the same number of normal samples ten times. Each time, t-test was used to estimate differential expression. The miRNAs with an FDR < 0.05 all 10 times were considered to be differentially expressed. To evaluate the potential of hub miRNAs with differential expression to be biomarkers, a scoring classifier was constructed. For each cancer, we first performed Z-score transformation on the expression levels across the samples for each hub miRNA with differential expression and then summarized the Z-scores as the integrated expression signature (score). Then, the samples could be divided into two classes (normal and tumour) by choosing a cut-off. In addition, an ROC curve that was drawn by plotting sensitivity against specificity was used for classifier evaluation. This procedure was performed using the R package pROC [91]. In addition, cluster analysis was performed using the R package pheatmap with default parameters. Identification of miRNA cooperative modules Next, miRNA functional cooperative modules, which were defined as cliques, were identified in the pan-cancer miRNA cooperative network using the clique percolation clustering method [92]. Cliques are all of the complete subgraphs, and all cliques in the network were identified using CFinder (2.0.6) [50]. Identification of survival-related miRNA cooperative modules To evaluate the prognostic performance of our miRNA cooperative modules identified above in each cancer, we first used univariate Cox regression analysis to evaluate the association between survival and the expression level of each miRNA in the module. A regression coefficient with a plus sign indicated that increased expression is associated with an increased risk of survival (risk miRNAs), and conversely, a minus sign indicated that increased expression is associated with a decreased risk of survival (protective miRNAs). After selecting the miRNAs that were significantly correlated with survival (P < 0.05) in the module, a mathematical formula for survival prediction was constructed, taking into account both the strength and positive or negative association of each miRNA with survival. More specifically, we assigned a risk score to each patient according to a linear combination of the miRNA expression values ExpmiRNAi weighted by the regression coefficients from the aforementioned univariate Cox regression analysis. The risk score for each patient was calculated as follows: Risk_score=∑i=1nβi*ExpmiRNAi, where βi is the Cox regression coefficient of miRNA i, and n is the number of miRNAs that were significantly associated with survival in the module. All patients were thus dichotomized into high-risk and low-risk groups using the median risk score as the cut-off point. Patients with higher risk scores were expected to have poor survival outcomes. The Kaplan–Meier method was further used to estimate the overall survival time for the two subgroups. The differences in the survival times were analysed using a log-rank test. Randomization test To determine if the miRNA-miRNA cooperative network is a small-world network, we used the duplication model to construct random graphs with the function of sample_gnm() in R. This is a well-known model with power-law degree distributions that provides small-world networks. We generated 1000 instances and computed the average clustering coefficients. The significant P-value was the fraction of the average clustering coefficients in random conditions, which is greater than the value in the real condition. Statistical analysis Because the regression coefficients do not follow a normal distribution, we used a Wilcoxon rank-sum test to evaluate differences in the mean of regression coefficients of cooperative miRNAs from that of non-cooperative miRNAs in each network. The miRNAs were descending sorted by expression variance, and the top 20.0% of the miRNAs were defined as miRNAs with greater variance. Fisher's exact test was used to determine if there were non-random associations between two cooperative properties and variance of expression. Key Points The miRNA-miRNA cooperative landscape was surveyed across cancer types. Pan-cancer miRNA-miRNA cooperative network exhibited a scale-free and modular architecture. Cooperative miRNAs exhibit higher expression, greater expression variation and a stronger regulatory strength towards targets as well as regulated cancer hallmark functions. Hub miRNAs regulate many cancer-associated functions and tend to target known anticancer drug targets. miRNA cooperative modules were found to serve as novel prognostic biomarkers. Supplementary Data Supplementary data are available online at https://academic.oup.com/bib. Funding The National High Technology Research and Development Program of China (863 Program, grant number 2014AA021102), the National Program on Key Basic Research Project (973 Program, grant number 2014CB910504), the funds for Creative Research Groups of the National Natural Science Foundation of China (grant number 81421063), the National Natural Science Foundation of China (grant numbers 91439117, 61473106, 31571331, 61502126, 31601065); the China Postdoctoral Science Foundation (grant numbers 2017M611385); Heilongjiang Province Postdoctoral Fund (grant number LBH-Z16152), Weihan Yu Youth Science Fund Project of Harbin Medical University, Harbin Special Funds of Innovative Talents on Science and Technology Research Project (grant number RC2015QN003080), the Harbin Medical University Innovation Science Research Fund (grant number 2016JCZX45) and the Funds for the Graduate Innovation Fund of Heilongjiang Province (grant number YJSCX2015-7HYD). Tingting Shao is an instructor in the College of Bioinformatics Science and Technology at Harbin Medical University, China. Her research activity focused on ncRNA regulation in complex diseases. Guangjuan Wang is an MS student in the College of Bioinformatics Science and Technology at Harbin Medical University, China. Her research interests focus on bioinformatics methods. Hong Chen is a PhD in the College of Bioinformatics Science and Technology at Harbin Medical University, China. Her research focuses on bioinformatics. Yunjin Xie is an undergraduate in the College of Bioinformatics Science and Technology at Harbin Medical University, China. Xiyun Jin is a master student in the College of Bioinformatics Science and Technology at Harbin Medical University, China. Her research focuses on bioinformatics. Jing Bai is an instructor in the College of Bioinformatics Science and Technology at Harbin Medical University, China. Her research interests focus on method development and ncRNA regulation. Juan Xu is an associate professor in the College of Bioinformatics Science and Technology at Harbin Medical University, China. Her researches focus on computational biology in complex diseases. Xia Li is a professor in the College of Bioinformatics Science and Technology at Harbin Medical University, China. Her research interests focus on ncRNA regulation in complex diseases. Jian Huang is a professor in the Center for Informational Biology at University of Electronic Science and Technology of China, China Yan Jin is a professor in the Department of Medical Genetics at Harbin Medical University, China. Yongsheng Li is an associate professor in the College of Bioinformatics Science and Technology at Harbin Medical University, China. His research interests focus on ncRNA regulation and bioinformatics methods development. References 1 Fujii YR. Quantum Language of MicroRNA: application for new cancer therapeutic targets . Methods Mol Biol 2018 ; 1733 : 145 – 57 . Google Scholar CrossRef Search ADS PubMed 2 He Z , Wang Y , Huang G , et al. The lncRNA UCA1 interacts with miR-182 to modulate glioma proliferation and migration by targeting iASPP . Arch Biochem Biophysics 2017 ; 623–624 : 1 – 8 . Google Scholar CrossRef Search ADS 3 Feng J , Xue S , Pang Q , et al. miR-141-3p inhibits fibroblast proliferation and migration by targeting GAB1 in keloids . Biochem Biophys Res Commun 2017 ; 490 ( 2 ): 302 – 8 . Google Scholar CrossRef Search ADS PubMed 4 Ziebarth JD , Bhattacharya A , Chen A , et al. PolymiRTS database 2.0: linking polymorphisms in microRNA target sites with human diseases and complex traits . Nucleic Acids Res 2012 ; 40 ( D1 ): D216 – 21 . Google Scholar CrossRef Search ADS PubMed 5 Deng Y , Bai H , Hu H. rs11671784 G/A variation in miR-27a decreases chemo-sensitivity of bladder cancer by decreasing miR-27a and increasing the target RUNX-1 expression . Biochem Biophys Res Commun 2015 ; 458 ( 2 ): 321 – 7 . Google Scholar CrossRef Search ADS PubMed 6 Naveed A , Ur-Rahman S , Abdullah S , et al. A concise review of MicroRNA exploring the insights of MicroRNA regulations in bacterial, viral and metabolic diseases . Mol Biotechnol 2017 ; 59 ( 11–12 ): 518 – 29 . Google Scholar CrossRef Search ADS PubMed 7 Kim D , Chang HR , Baek D. Rules for functional microRNA targeting . BMB Rep 2017 ; 50 ( 11 ): 554 – 9 . Google Scholar CrossRef Search ADS PubMed 8 Kozomara A , Griffiths-Jones S. miRBase: annotating high confidence microRNAs using deep sequencing data . Nucleic Acids Res 2014 ; 42 ( D1 ): D68 – 73 . Google Scholar CrossRef Search ADS PubMed 9 Xu J , Li CX , Li YS , et al. MiRNA-miRNA synergistic network: construction via co-regulating functional modules and disease miRNA topological features . Nucleic Acids Res 2011 ; 39 ( 3 ): 825 – 36 . Google Scholar CrossRef Search ADS PubMed 10 Wu S , Huang S , Ding J , et al. Multiple microRNAs modulate p21Cip1/Waf1 expression by directly targeting its 3' untranslated region . Oncogene 2010 ; 29 ( 15 ): 2302 – 8 . Google Scholar CrossRef Search ADS PubMed 11 Bandi N , Vassella E. miR-34a and miR-15a/16 are co-regulated in non-small cell lung cancer and control cell cycle progression in a synergistic and Rb-dependent manner . Mol Cancer 2011 ; 10 ( 1 ): 55 . Google Scholar CrossRef Search ADS PubMed 12 Kandi R , Gutti U , Saladi RGV , et al. MiR-125b and miR-99a encoded on chromosome 21 co-regulate vincristine resistance in childhood acute megakaryoblastic leukemia . Hematol Oncol Stem Cell Ther 2015 ; 8 ( 2 ): 95 – 7 . Google Scholar CrossRef Search ADS PubMed 13 Zhao Y , Cui X , Zhu W , et al. Synergistic regulatory effects of microRNAs on brain glioma cells . Mol Med Rep 2017 ; 16 ( 2 ): 1409 – 16 . Google Scholar CrossRef Search ADS PubMed 14 Borzi C , Calzolari L , Centonze G , et al. mir-660-p53-mir-486 network: a new key regulatory pathway in lung tumorigenesis . Int J Mol Sci 2017 ; 18 ( 1 ): 222 . Google Scholar CrossRef Search ADS 15 Lai X , Gupta SK , Schmitz U , et al. MiR-205-5p and miR-342-3p cooperate in the repression of the E2F1 transcription factor in the context of anticancer chemotherapy resistance . Theranostics 2018 ; 8 ( 4 ): 1106 – 20 . Google Scholar CrossRef Search ADS PubMed 16 Frampton AE , Castellano L , Colombo T , et al. MicroRNAs cooperatively inhibit a network of tumor suppressor genes to promote pancreatic tumor growth and progression . Gastroenterology 2014 ; 146 ( 1 ): 268 – 77.e18 . Google Scholar CrossRef Search ADS PubMed 17 Chen HH , Huang WT , Yang LW , et al. The PTEN-AKT-mTOR/RICTOR pathway in nasal natural killer cell lymphoma is activated by miR-494-3p via PTEN but inhibited by miR-142-3p via RICTOR . Am J Pathol 2015 ; 185 ( 5 ): 1487 – 99 . Google Scholar CrossRef Search ADS PubMed 18 Xu J , Shao T , Ding N , et al. miRNA–miRNA crosstalk: from genomics to phenomics . Brief Bioinform 2017 ; 18 ( 6 ): 1002 – 11 . Google Scholar PubMed 19 Ying H , Lyu J , Ying T , et al. RETRACTED ARTICLE: risk miRNA screening of ovarian cancer based on miRNA functional synergistic network . J Ovarian Res 2014 ; 7 ( 1 ): 226 . Google Scholar CrossRef Search ADS 20 Alshalalfa M. MicroRNA response elements-mediated miRNA-miRNA interactions in prostate cancer . Adv Bioinform 2012 ; 2012 : 839837 . 21 Hua L , Zhou P , Li L , et al. Prioritizing breast cancer subtype related miRNAs using miRNA–mRNA dysregulated relationships extracted from their dual expression profiling . J Theor Biol 2013 ; 331 : 1 – 11 . Google Scholar CrossRef Search ADS PubMed 22 Yang Y , Xing Y , Liang C , et al. Crucial microRNAs and genes of human primary breast cancer explored by microRNA-mRNA integrated analysis . Tumor Biol 2015 ; 36 ( 7 ): 5571 – 9 . Google Scholar CrossRef Search ADS 23 Zhang TF , Cheng KW , Shi WY , et al. MiRNA synergistic network construction and enrichment analysis for common target genes in small-cell lung cancer . Asian Pac J Cancer Prev 2012 ; 13 ( 12 ): 6375 – 8 . Google Scholar CrossRef Search ADS PubMed 24 Li Y , Liang C , Wong K-C , et al. Mirsynergy: detecting synergistic miRNA regulatory modules by overlapping neighbourhood expansion . Bioinformatics 2014 ; 30 ( 18 ): 2627 – 35 . Google Scholar CrossRef Search ADS PubMed 25 Zhao X , Song H , Zuo Z , et al. Identification of miRNA-miRNA synergistic relationships in colorectal cancer . Int J Biol Macromol 2013 ; 55 : 98 – 103 . Google Scholar CrossRef Search ADS PubMed 26 Na YJ , Kim JH. Understanding cooperativity of microRNAs via microRNA association networks . BMC Genomics 2013 ; 14 ( Suppl 5 ): S17 . Google Scholar CrossRef Search ADS PubMed 27 Xiao Y , Ping Y , Fan H , et al. Identifying dysfunctional miRNA-mRNA regulatory modules by inverse activation, cofunction, and high interconnection of target genes: a case study of glioblastoma . Neuro Oncol 2013 ; 15 ( 7 ): 818 – 28 . Google Scholar CrossRef Search ADS PubMed 28 Zhu W , Zhao Y , Xu Y , et al. Dissection of protein interactomics highlights microRNA synergy . PLoS One 2013 ; 8 ( 5 ): e63342 . Google Scholar CrossRef Search ADS PubMed 29 Meng X , Wang J , Yuan C , et al. CancerNet: a database for decoding multilevel molecular interactions across diverse cancer types . Oncogenesis 2015 ; 4 : e177. Google Scholar CrossRef Search ADS PubMed 30 Song R , Catchpoole DR , Kennedy PJ , et al. Identification of lung cancer miRNA-miRNA co-regulation networks through a progressive data refining approach . J Theor Biol 2015 ; 380 : 271 – 9 . Google Scholar CrossRef Search ADS PubMed 31 The Cancer Genome Atlas Research Network , Weinstein JN , Collisson EA , et al. The Cancer Genome Atlas Pan-Cancer analysis project . Nat Genet 2013 ; 45 : 1113 – 20 . Google Scholar CrossRef Search ADS PubMed 32 The International Cancer Genome Consortium , Hudson TJ , Anderson W , et al. International network of cancer genome projects . Nature 2010 ; 464 : 993 – 8 . Google Scholar CrossRef Search ADS PubMed 33 Jacobsen A , Silber J , Harinath G , et al. Analysis of microRNA-target interactions across diverse cancer types . Nat Struct Mol Biol 2013 ; 20 ( 11 ): 1325 – 32 . Google Scholar CrossRef Search ADS PubMed 34 Yang JH , Li JH , Shao P , et al. starBase: a database for exploring microRNA-mRNA interaction maps from Argonaute CLIP-Seq and Degradome-Seq data . Nucleic Acids Res 2011 ; 39 : D202 – 9 . Google Scholar CrossRef Search ADS PubMed 35 Chi SW , Zang JB , Mele A , et al. Argonaute HITS-CLIP decodes microRNA-mRNA interaction maps . Nature 2009 ; 460 ( 7254 ): 479 – 86 . Google Scholar CrossRef Search ADS PubMed 36 Luna JM , Barajas JM , Teng KY , et al. Argonaute CLIP defines a deregulated miR-122-bound transcriptome that correlates with patient survival in human liver cancer . Mol Cell 2017 ; 67 ( 3 ): 400 – 10.e7 . Google Scholar CrossRef Search ADS PubMed 37 Zhang XQ , Yang JH. Discovering circRNA-microRNA interactions from CLIP-Seq Data . Methods Mol Biol 2018 ; 1724 : 193 – 207 . Google Scholar CrossRef Search ADS PubMed 38 Korpal M , Lee ES , Hu G , et al. The miR-200 family inhibits epithelial-mesenchymal transition and cancer cell migration by direct targeting of E-cadherin transcriptional repressors ZEB1 and ZEB2 . J Biol Chem 2008 ; 283 ( 22 ): 14910 – 4 . Google Scholar CrossRef Search ADS PubMed 39 Yoneyama K , Ishibashi O , Kawase R , et al. miR-200a, miR-200b and miR-429 are onco-miRs that target the PTEN gene in endometrioid endometrial carcinoma . Anticancer Res 2015 ; 35 : 1401 – 10 . Google Scholar PubMed 40 Hoadley KA , Yau C , Wolf DM , et al. Multiplatform analysis of 12 cancer types reveals molecular classification within and across tissues of origin . Cell 2014 ; 158 : 929 – 44 . Google Scholar CrossRef Search ADS PubMed 41 Rosenfeld N , Aharonov R , Meiri E , et al. MicroRNAs accurately identify cancer tissue origin . Nat Biotechnol 2008 ; 26 ( 4 ): 462 – 9 . Google Scholar CrossRef Search ADS PubMed 42 Xu Q , Chen J , Ni S , et al. Pan-cancer transcriptome analysis reveals a gene expression signature for the identification of tumor tissue origin . Mod Pathol 2016 ; 29 ( 6 ): 546 – 56 . Google Scholar CrossRef Search ADS PubMed 43 Mar JC , Matigian NA , Mackay-Sim A , et al. Variance of gene expression identifies altered network constraints in neurological disease . PLoS Genetics 2011 ; 7 ( 8 ): e1002207 . Google Scholar CrossRef Search ADS PubMed 44 Wang D , Gu J , Wang T , et al. OncomiRDB: a database for the experimentally verified oncogenic and tumor-suppressive microRNAs . Bioinformatics 2014 ; 30 ( 15 ): 2237 – 8 . Google Scholar CrossRef Search ADS PubMed 45 Wang E , Zaman N , McGee S , et al. Predictive genomics: a cancer hallmark network framework for predicting tumor clinical phenotypes using genome sequencing data . Semin Cancer Biol 2015 ; 30 : 4 – 12 . Google Scholar CrossRef Search ADS PubMed 46 Zotenko E , Mestre J , O'Leary DP , et al. Why do hubs in the yeast protein interaction network tend to be essential: reexamining the connection between the network topology and essentiality . PLoS Comput Biol 2008 ; 4 ( 8 ): e1000140 . Google Scholar CrossRef Search ADS PubMed 47 Lu M , Shi B , Wang J , et al. TAM: a method for enrichment and depletion analysis of a microRNA category in a list of microRNAs . BMC Bioinformatics 2010 ; 11 ( 1 ): 419 . Google Scholar CrossRef Search ADS PubMed 48 Li J , Lenferink AEG , Deng Y , et al. Corrigendum: identification of high-quality cancer prognostic markers and metastasis network modules . Nat Commun 2012 ; 3 ( 1 ): 655 . Google Scholar CrossRef Search ADS 49 Wu G , Stein L. A network module-based method for identifying cancer prognostic signatures . Genome Biol 2012 ; 13 ( 12 ): R112. Google Scholar CrossRef Search ADS PubMed 50 Niklas N , Hafenscher J , Barna A , et al. cFinder: definition and quantification of multiple haplotypes in a mixed sample . BMC Res Notes 2015 ; 8 ( 1 ): 422 . Google Scholar CrossRef Search ADS PubMed 51 Watts DJ , Strogatz SH. Collective dynamics of ‘small-world’ networks . Nature 1998 ; 393 ( 6684 ): 440 – 2 . Google Scholar CrossRef Search ADS PubMed 52 Park S , Lehner B. Cancer type‐dependent genetic interactions between cancer driver alterations indicate plasticity of epistasis across cell types . Mol Syst Biol 2015 ; 11 ( 7 ): 824. Google Scholar CrossRef Search ADS PubMed 53 Xu J , Li Y , Lu J , et al. The mRNA related ceRNA-ceRNA landscape and significance across 20 major cancer types . Nucleic Acids Res 2015 ; 43 ( 17 ): 8169 – 82 . Google Scholar CrossRef Search ADS PubMed 54 Hamilton MP , Rajapakshe K , Hartig SM , et al. Identification of a pan-cancer oncogenic microRNA superfamily anchored by a central core seed motif . Nat Commun 2013 ; 4 : 2730. Google Scholar CrossRef Search ADS PubMed 55 Davoli T , Xu Andrew W , Mengwasser KE , et al. Cumulative haploinsufficiency and triplosensitivity drive aneuploidy patterns and shape the cancer genome . Cell 2013 ; 155 ( 4 ): 948 – 62 . Google Scholar CrossRef Search ADS PubMed 56 Heim D , Budczies J , Stenzinger A , et al. Cancer beyond organ and tissue specificity: next‐generation‐sequencing gene mutation data reveal complex genetic similarities across major cancers . Int J Cancer 2014 ; 135 ( 10 ): 2362 – 9 . Google Scholar CrossRef Search ADS PubMed 57 Zhao X , Modur V , Carayannopoulos LN , et al. Biomarkers in pharmaceutical research . Clin Chem 2015 ; 61 ( 11 ): 1343 – 53 . Google Scholar CrossRef Search ADS PubMed 58 Saito T , Sætrom P. Target gene expression levels and competition between transfected and endogenous microRNAs are strong confounding factors in microRNA high-throughput experiments . Silence 2012 ; 3 ( 1 ): 3. Google Scholar CrossRef Search ADS PubMed 59 Mendell JT. miRiad roles for the miR-17-92 cluster in development and disease . Cell 2008 ; 133 ( 2 ): 217 – 22 . Google Scholar CrossRef Search ADS PubMed 60 Di Lena M , Travaglio E , Altomare DF. New strategies for colorectal cancer screening . World J Gastroenterol 2013 ; 19 ( 12 ): 1855 – 60 . Google Scholar CrossRef Search ADS PubMed 61 Wu HH , Lin WC , Tsai KW. Advances in molecular biomarkers for gastric cancer: miRNAs as emerging novel cancer markers . Expert Rev Mol Med 2014 ; 16 : e1 . Google Scholar CrossRef Search ADS PubMed 62 Chen X , Ba Y , Ma L , et al. Characterization of microRNAs in serum: a novel class of biomarkers for diagnosis of cancer and other diseases . Cell Res 2008 ; 18 ( 10 ): 997 – 1006 . Google Scholar CrossRef Search ADS PubMed 63 Cui L , Zhang X , Ye G , et al. Gastric juice MicroRNAs as potential biomarkers for the screening of gastric cancer . Cancer 2013 ; 119 ( 9 ): 1618 – 26 . Google Scholar CrossRef Search ADS PubMed 64 Zhang X , Cui L , Ye G , et al. Gastric juice microRNA-421 is a new biomarker for screening gastric cancer, Tumour biology: the journal of the International Society for Oncodevelopmental . Biol Med 2012 ; 33 ( 6 ): 2349 – 55 . 65 Yu X , Luo L , Wu Y , et al. Gastric juice miR-129 as a potential biomarker for screening gastric cancer . Med Oncol 2013 ; 30 ( 1 ): 365 . Google Scholar CrossRef Search ADS PubMed 66 Li C , Li JF , Cai Q , et al. MiRNA-199a-3p in plasma as a potential diagnostic biomarker for gastric cancer . Ann Surg Oncol 2013 ; 20 ( S3 ): 397 – 405 . Google Scholar CrossRef Search ADS 67 Li BS , Zhao Yl , Guo G , et al. Plasma microRNAs, miR-223, miR-21 and miR-218, as novel potential biomarkers for gastric cancer detection . PLoS One 2012 ; 7 ( 7 ): e41629 . Google Scholar CrossRef Search ADS PubMed 68 Yun SJ , Jeong P , Kim WT , et al. Cell-free microRNAs in urine as diagnostic and prognostic biomarkers of bladder cancer . Int J Oncol 2012 ; 41 ( 5 ): 1871 – 8 . Google Scholar CrossRef Search ADS PubMed 69 Park NJ , Zhou H , Elashoff D , et al. Salivary microRNA: discovery, characterization, and clinical utility for oral cancer detection . Clin Cancer Res 2009 ; 15 ( 17 ): 5473 . Google Scholar CrossRef Search ADS PubMed 70 Gallardo E , Navarro A , Viñolas N , et al. miR-34a as a prognostic marker of relapse in surgically resected non-small-cell lung cancer . Carcinogenesis 2009 ; 30 ( 11 ): 1903 – 9 . Google Scholar CrossRef Search ADS PubMed 71 Wang Q , Li P , Li A , et al. Plasma specific miRNAs as predictive biomarkers for diagnosis and prognosis of glioma . J Exp Clin Cancer Res 2012 ; 31 ( 1 ): 97 . Google Scholar CrossRef Search ADS PubMed 72 Zhang L , Shan X , Wang J , et al. A three-microRNA signature for lung squamous cell carcinoma diagnosis in Chinese male patients . Oncotarget 2017 ; 8 ( 49 ): 86897 – 907 . Google Scholar PubMed 73 Li JH , Liu S , Zhou H , et al. starBase v2.0: decoding miRNA-ceRNA, miRNA-ncRNA and protein-RNA interaction networks from large-scale CLIP-Seq data . Nucleic Acids Res 2014 ; 42 ( D1 ): D92 – 7 . Google Scholar CrossRef Search ADS PubMed 74 Grimson A , Farh KK , Johnston WK , et al. MicroRNA targeting specificity in mammals: determinants beyond seed pairing . Mol Cell 2007 ; 27 ( 1 ): 91 – 105 . Google Scholar CrossRef Search ADS PubMed 75 Betel D , Koppal A , Agius P , et al. Comprehensive modeling of microRNA targets predicts functional non-conserved and non-canonical sites . Genome Biol 2010 ; 11 ( 8 ): R90 . Google Scholar CrossRef Search ADS PubMed 76 Anders G , Mackowiak SD , Jens M , et al. doRiNA: a database of RNA interactions in post-transcriptional regulation . Nucleic Acids Res 2012 ; 40 ( D1 ): D180 – 6 . Google Scholar CrossRef Search ADS PubMed 77 Kertesz M , Iovino N , Unnerstall U , et al. The role of site accessibility in microRNA target recognition . Nat Genet 2007 ; 39 ( 10 ): 1278 – 84 . Google Scholar CrossRef Search ADS PubMed 78 Miranda KC , Huynh T , Tay Y , et al. A pattern-based method for the identification of MicroRNA binding sites and their corresponding heteroduplexes . Cell 2006 ; 126 ( 6 ): 1203 – 17 . Google Scholar CrossRef Search ADS PubMed 79 Ru Y , Kechris KJ , Tabakoff B , et al. The multiMiR R package and database: integration of microRNA–target interactions along with their disease and drug associations . Nucleic Acids Res 2014 ; 42 ( 17 ): e133 . Google Scholar CrossRef Search ADS PubMed 80 Xiao F , Zuo Z , Cai G , et al. miRecords: an integrated resource for microRNA-target interactions . Nucleic Acids Res 2009 ; 37 : D105 – 10 . Google Scholar CrossRef Search ADS PubMed 81 Hsu SD , Lin FM , Wu WY , et al. miRTarBase: a database curates experimentally validated microRNA–target interactions . Nucleic Acids Res 2011 ; 39(Suppl 1) : D163 – 9 . Google Scholar CrossRef Search ADS 82 Sethupathy P , Corda B , Hatzigeorgiou AG. TarBase: a comprehensive database of experimentally supported animal microRNA targets . RNA 2005 ; 12 ( 2 ): 192 – 7 . Google Scholar CrossRef Search ADS PubMed 83 Mermel CH , Schumacher SE , Hill B , et al. GISTIC2.0 facilitates sensitive and confident localization of the targets of focal somatic copy-number alteration in human cancers . Genome Biol 2011 ; 12 ( 4 ): R41 . Google Scholar CrossRef Search ADS PubMed 84 Liberzon A , Subramanian A , Pinchback R , et al. Molecular signatures database (MSigDB) 3.0 . Bioinformatics 2011 ; 27 ( 12 ): 1739 – 40 . Google Scholar CrossRef Search ADS PubMed 85 Plaisier CL , Pan M , Baliga NS. A miRNA-regulatory network explains how dysregulated miRNAs perturb oncogenic processes across diverse cancers . Genome Res 2012 ; 22 ( 11 ): 2302 – 14 . Google Scholar CrossRef Search ADS PubMed 86 Keshava Prasad TS , Goel R , Kandasamy K , et al. Human protein reference database–2009 update . Nucleic Acids Res 2009 ; 37 : D767 – 72 . Google Scholar CrossRef Search ADS PubMed 87 Li Y , Qiu C , Tu J , et al. HMDD v2.0: a database for experimentally supported human microRNA and disease associations . Nucleic Acids Res 2014 ; 42 : D1070 – 4 . Google Scholar CrossRef Search ADS PubMed 88 Jiang Q , Wang Y , Hao Y , et al. miR2Disease: a manually curated database for microRNA deregulation in human disease . Nucleic Acids Res 2009 ; 37 : D98 – 104 . Google Scholar CrossRef Search ADS PubMed 89 Xie B , Ding Q , Han H , et al. miRCancer: a microRNA-cancer association database constructed by text mining on literature . Bioinformatics 2013 ; 29 ( 5 ): 638 – 44 . Google Scholar CrossRef Search ADS PubMed 90 Hsu CW , Juan HF , Huang HC. Characterization of microRNA-regulated protein-protein interaction network . Proteomics 2008 ; 8 ( 10 ): 1975 – 9 . Google Scholar CrossRef Search ADS PubMed 91 Robin X , Turck N , Hainard A , et al. pROC: an open-source package for R and S+ to analyze and compare ROC curves . BMC Bioinformatics 2011 ; 12 ( 1 ): 77 . Google Scholar CrossRef Search ADS PubMed 92 Palla G , Derenyi I , Farkas I , et al. Uncovering the overlapping community structure of complex networks in nature and society . Nature 2005 ; 435 ( 7043 ): 814 – 18 . Google Scholar CrossRef Search ADS PubMed © The Author(s) 2018. Published by Oxford University Press. All rights reserved. For permissions, please email: journals.permissions@oup.com This article is published and distributed under the terms of the Oxford University Press, Standard Journals Publication Model (https://academic.oup.com/journals/pages/about_us/legal/notices)

Journal

Briefings in BioinformaticsOxford University Press

Published: May 25, 2018

There are no references for this article.

You’re reading a free preview. Subscribe to read the entire article.


DeepDyve is your
personal research library

It’s your single place to instantly
discover and read the research
that matters to you.

Enjoy affordable access to
over 18 million articles from more than
15,000 peer-reviewed journals.

All for just $49/month

Explore the DeepDyve Library

Search

Query the DeepDyve database, plus search all of PubMed and Google Scholar seamlessly

Organize

Save any article or search result from DeepDyve, PubMed, and Google Scholar... all in one place.

Access

Get unlimited, online access to over 18 million full-text articles from more than 15,000 scientific journals.

Your journals are on DeepDyve

Read from thousands of the leading scholarly journals from SpringerNature, Elsevier, Wiley-Blackwell, Oxford University Press and more.

All the latest content is available, no embargo periods.

See the journals in your area

DeepDyve

Freelancer

DeepDyve

Pro

Price

FREE

$49/month
$360/year

Save searches from
Google Scholar,
PubMed

Create lists to
organize your research

Export lists, citations

Read DeepDyve articles

Abstract access only

Unlimited access to over
18 million full-text articles

Print

20 pages / month

PDF Discount

20% off