Identification of novel prognostic indicators for triple-negative breast cancer patients through integrative analysis of cancer genomics data and protein interactome data

Triple negative breast cancers (TNBCs) are highly heterogeneous and aggressive without targeted treatment. Here, we aim to systematically dissect TNBCs from a prognosis point of view by building a subnetwork atlas for TNBC prognosis through integrating multi-dimensional cancer genomics data from The Cancer Genome Atlas (TCGA) project and the interactome data from three different interaction networks. The subnetworks are represented as the protein-protein interaction modules perturbed by multiple genetic and epigenetic interacting mechanisms contributing to patient survival. Predictive power of these subnetwork-derived prognostic models is evaluated using Monte Carlo cross-validation and the concordance index (C-index). We uncover subnetwork biomarkers of low oncogenic GTPase activity, low ubiquitin/proteasome degradation, effective protection from oxidative damage, and tightly immune response are linked to better prognosis. Such a systematic approach to integrate massive amount of cancer genomics data into clinical practice for TNBC prognosis can effectively dissect the molecular mechanisms underlying TNBC patient outcomes and provide potential opportunities to optimize treatment and develop therapeutics.


S1 Text. Biological insights from network-based prognostic biomarkers
Among the subnetworks detected from the three different interaction networks, most strikingly, one subnetwork(inthint1/iref3/multinet2) was shared by the three PPI networks and two were shared by iRefIndex and MultiNet (iref2/multinet1, iref4/multinet6) ( Figure 2). To further understand how these subnetwork biomarkers may be related to patient survival, we investigated individual gene using NCBI Entrez Gene [1] and literature search and performed functional enrichment analysis based on the known pathways or ontologies using Enrichr [2].
Such a thorough analysis identified several pathways known to be involved in TNBC such asiref4/multinet6, which contain BRCA1 interacting partners (RBBP7 by high mRNA, low methyl; BAP1 by high mRNA and CNV) together with transcriptional repressors (DHX30 by high CNV; HSPA2 by high mRNA, low methyl) and PIK3A signaling pathway (discuss in detail later). BRCA1 inactivation and PIK3CA high activity have been identified as key features in basal-like breast cancers (often referred to as TNBCs because of 75% overlap of mRNA profiles) [3]. Several metabolism-related subnetworks have also been identified, such as multinet3 (discuss in detail later) and multinet4/8/9/12. Members in multinet4/8/9/12 are mainly involved amino acid metabolism and enriched in mitochondria from Enrichr analysis. We found that elevated level of several synthetases involved in one carbon unit (MAT2B, MTHFD1, MTHFS, MTR), high presence of GPX1 (one of the most important antioxidant enzymes), and low presence of various enzymes involved in biosynthesis (CBS, PC, SPTLC2, NADSYN1, PAICS) were linked to better prognosis. Though we cannot reach a conclusion of which particular product/metabolic or activating or repressing which particular pathway is linked to better prognosis, such a phenomenon leads us to propose that these subnetworks together with various acetyl-CoA related enzymes in multient3 may be related to the Warburg effect in tumors, such as reactive oxygen (ROS) and biosynthesis [4]. Most importantly, we found quite large amount of a subnetworks related to GTPase, endoplasmic reticulum (ER)-Golgi-cell surface trafficking, ubiquitin/proteasome system, and complement system. Here, we discuss these top-ranked and shared subnetworks in detail ( Figure 3).

iref12/inthint2/multinet11 and inthint3/ multinet14
Complement related subnetworks stood out from our analysis. Complement system is well known to play key roles in innate and adaptive immunity [5]. Its role in cancer is a bit controversial: its has anti-tumor effect by direct complement attacked or complement-dependent cytotoxicity; it has tumor promoting effect by C5a induced ROS; tumors can also escape complement attack by expressing or recruiting complement regulators or inhibitors to inhibit the formation of terminal complement complex [6]. Previous studies have shown that that low activity of the classical complement pathway predicts short survival of patients with chronic lymphocytic leukaemia [7] and activation of complement and immune response pathways is associated with good prognosis in ER-negative breast cancer [8]. In fact, high presence of the majority of the genes in these subnetworks is linked to better survival (high mRNA_C1R, C1S, SERPING1, SERPINI1, CYTH1, VNN2; high CNV_ARF1, low methyl_CFB, CR1, CYTIP, ICAM3, ITGAM, ITGB2, PLAUR). In contrast, low presence of C1QA, C3, CFD, SURF2 by high methylation is also linked to better survival. Particularly, though high mRNA_C1R, C1S, SERPING1, SERPINI1 are linked to better survival, serpin family proteins (SERPING1, SERPINI1) apparently work as inhibitors of C1R and C1S [9]. Compared with inthint4/iref3/multinet12, which are more related to the peptidase activity, inthint3/multinet14 are more related to glycoprotein binding and guanyl-nucleotide exchange factor activity by Enrichr analysis. Therefore, our results may not fully support the idea that high complement activity is linked to better survival, but a sensitive and tightly regulated complement system is associated with patient survival.

Supplementary Materials
C-terminal cleaved UFM1 by forming a high-energy thioester bond, and UFM1 was further transferred to E2like enzyme UFC1 via thioester linkage [10]. Knockdown or depletion of UFL1 or UBA5 could lead to unfolded protein response (UPR), ER stress, and inhibited vesicle trafficking [11][12][13]. Knockdown of UBA5 has been shown to inhibit breast cancer cell growth [14], implying the low activity of this ubiquitination system is beneficial for patient survival. APEH (acylaminoacyl-peptide hydrolase) encodes the enzyme called acylpeptide hydrolase, which plays an important role in destroying oxidatively damaged proteins in living cells [15]. Deletion of APEH locus have been found in small cell lung carcinoma [16] and renal cell carcinoma [17]. Therefore, high activity of APEH is beneficial to patient survival by reducing the oxidative proteins. CBR1 (carbonyl reductase 1) has protective roles against xenobiotic carbonyls [18], thus it may reduce the effect of many anti-cancer drugs by metabolization, such as reducing the effect of doxorubicin in breast cancer patients [19,20]. It could also attenuate the effect of arsenic trioxide in leukemia patients [21]. Thus, low level of CBR1 (CNV) is linked to better prognosis. As a mammalian lectin, LGALS8 (galectin 8) can inhibit cell adhesion and induce apoptosis by binding to integrin a4αβ3 to modulate cell-matrix interactions [22], which agrees with our prediction that high LGALS8 is linked to better survival. In summary, this shared subnetwork implied that low activity of ubiquitin/proteasome degradation by low presence of UBA5/UFM1, and effective removal of oxidatively damaged proteins by APEH are linked to better prognosis.

iref2/multinet1
RAB6A is a member of RAS oncogene family localized at Golgi apparatus [23], and RABGAP1 is a RAB6A activating protein playing a role in the coordination of microtubule and Golgi dynamics during the cell cycle [24]. CHML (CHM like, Rab escort protein 2) can facilitate the attachment of geranylgeranyl groups to several Rab proteins [25], while free CHML inhibits the prenylation reaction on RABs by RabGGTase by forming a CHML-RabGGTase complex [26]. TMF1 (TATA element modulatory factor 1) is a conserved Golgi protein binds to RAB6 and influences Golgi morphology [27], depletion of which blocks both retrograde membrane transport from endosome to trans-Golgi and RAB6dependent cargo-specific transport from Golgi to ER (i.e. retention of glycosylation enzymes at Golgi) [28]. Downregulation of TMF has been shown in solid tumors while overexpression of TMF significantly attenuated the development and growth of xenograft tumors by decreasing proangiogenic genes in metabolic stress state (i.e. direct the ubiquitination and proteasomal degradation of p65/RelA) [29]. Therefore, our prediction that low RAB oncogenic activity of by low CNV_RAB6A and mRNA_ RABGAP1 and high level of their inhibitory factors CHML (mRNA, CNV) and CNV_TMF1 are beneficial for patient survival agrees well with previous studies.

iref4/multinet6
This network contains a tightly regulated RAS signaling pathway. PIK3A is a well-known oncogene contributing to Warburg effect and affecting various pathways related to cancer (i.e. glucose transport and catabolism, cell adhesion, RAS signaling) and high CNV_PIK3A is commonly found in ovarian cancer biosynthesis [4,30,31]. Meanwhile, we reveal that low methyl_DGKZ and high mRNA_PPM1B are linked better survival. DKGZ (diacylglycerol kinase zeta) can inhibit both Ras-ERK and PKC pathways by regulating diacylglycerol (DAG) level [32,33]. DKGZ can also bind to retinoblastoma protein (pRB) or act as a downstream effecter of pRB to reconstitute cell cycle arrest in pRBnull fibroblasts [34]. As a member of PP2C family of Ser/Thr protein phosphatases, PPM1B (phosphatase, Mg2+/Mn2+ dependent 1B) can dephosphorylate cyclindependent kinases (CDKs) to inhibit cell growth [35], IkappaB kinases (IKKs) to down-regulate cytokineinduced NF-kappaB activation [36], and stress-activated protein kinase (a member of MAPKs) to attenuate stress response [37]. Therefore, we revealed a subnetwork that high PIK3A by high mRNA and CNV is linked to worse survival, while high presence of DKGZ and PPM1B, two inhibitory proteins, are linked to better survival.

iref14/multinet3
Besides sharing three high-heat score genes, SRPRB, UQCR1, and USP19, iref14 is more transcription related while mMultinet3 is related to CoA metabolism. SPRRB (SRP receptor beta subunit) is a GTPase associated with ER membrane [38] and highly expressed in various tumors [39]. UQCRC1 (ubiquinol-cytochrome c reductase core protein I) is a critical component in mitochondrial respiratory electron transport chain [40]. It also plays a role in mitochondria-to-nucleus retrograde response invoked by depletion of mitochondrial DNA (mtDNA) and its expression is positively correlated with cytochrome c-oxidase encoded by mtDNA [41]. USP19 (ubiquitin specific peptidase 19) is an ER membrane-anchored deubiquitinating enzyme to rescue ER-associated degradation (ERAD) stubstrates from proteasomal degradation [42] and mount an appropriate response to hypoxia by rescuing hypoxia-inducible factor 1α (HIF-1α) from degradation [43]. Our prediction is low activity of SRPRB (mRNA, CNV) and high activity of UQCRC1 (CNV) and USP19 (mRNA, CNV) are linked to better survival, implying general low activities of GTPase (SRPRB) and ubiquitin/proteasome degradation (by high activity of deubiquitinating enzyme USP19) as well as normal mitochondrial respiratory electron transport chain (UQCRC1) are linked to better prognosis. Furthermore, we would like to point out two one pathways highly associated with TNBC patient survival: iref56, semaphorins-plexins signaling pathway with RAS family small GTPases as downstream signaling [44] and iref18, nidogen proteins, members of extracellular matrix (ECM) which can interact with integrins and their expression is reduced in tumors [45]. Additionally, we would like to highlight some of the genes that we do not explain the relative subnetworks in detail but are recently reported as novel tumor suppressors linked to cancer prognosis, including PHLDA3 (a p53 target gene) by repressing AKT [46,47], GGT7 by regulating anti-oxidative damage [48], and CTNNA1 (α-catenin) by inhibiting Wnt/β-catenin, NF-κB, Hippo-YAP, and Hedgehog signaling pathways [49].