Functional Stratification of Cancer Drugs through Integrated Network Similarity

doi:10.21203/rs.3.rs-646155/v1

Download PDF

Article

Functional Stratification of Cancer Drugs through Integrated Network Similarity

https://doi.org/10.21203/rs.3.rs-646155/v1

This work is licensed under a CC BY 4.0 License

Journal Publication

published 19 Apr, 2022

Read the published version in npj Systems Biology and Applications →

You are reading this latest preprint version

Heterogeneity across tumors is the main obstacle in developing treatment strategies. Drug molecules not only perturb their immediate protein targets but also modulate multiple signaling pathways. In this study, we explored the networks modulated by several drug molecules across multiple cancer cell lines by integrating the drug targets with transcriptomic and phosphoproteomic data. As a result, we obtained 236 reconstructed networks covering five cell lines and 70 drugs. A rigorous topological and pathway analysis showed that chemically and functionally different drugs may modulate overlapping networks. Additionally, we revealed a set of tumor-specific hidden pathways with the help of drug network models that are not detectable from the initial data. The difference in the target selectivity of the drugs leads to disjoint networks despite sharing the exact mechanism of action, e.g., HDAC inhibitors. We also used the reconstructed network models to study potential drug combinations based on the topological separation, found literature evidence for a set of drug pairs. Overall, the network-level exploration of the drug perturbations may potentially help optimize treatment strategies and suggest new drug combinations.

Network reconstruction

data integration

drug modulation

drug-disease-protein interactions

Transformation of normal cells to tumor cells is a multi-stage process where multiple signaling pathways and biomolecular connections alter^1–7. However, response to the drug treatment is highly dependent on cellular and physiological factors in cancer, and many drugs have multiple targets⁸. The molecular heterogeneity across the tumor types may result in different signaling alterations for the same drug^9–11. Moreover, drugs perturb simultaneously a set of protein networks besides their immediate targets. Therefore, network-based approaches may address several unknowns about the drug efficacy and mechanisms of action across various cancer types^12–14. Integrative multi-omic approaches can provide a realistic view of network-level alterations toward developing better treatment strategies against the complexity and heterogeneity of cancer^15–19.

Different network-based strategies were previously used to investigate protein-drug and protein-protein interactions. The topological separation of disease modules within the human interactome was studied to elucidate disease-disease interactions¹². In a similar line, the network proximity of drug targets in the human interactome was considered to predict drug-disease interactions for drug repurposing¹³. Cheng et al. measured the closest distances between drug targets and disease genes and used patient-specific data to investigate the effects of drugs on different diseases¹⁴. Ritz et al. developed PathLinker to reconstruct human signaling pathways, which finds multiple shortest paths from receptors to transcriptional regulators in a reference protein interactome²⁰. Wu et al. used BioNetGen software²¹ to model a detailed rule-based biochemical network of VEGF-mediated eNOS signaling pathway to interpret the effects of an angiogenic inhibitor Thrombospondin-1 (TSP1)²². Halasz et al. studied signal transduction networks in colorectal cancer by integrated network reconstruction using Bayesian mechanistic modeling algorithm ²³. Naldi et al. combined weighted shortest paths and random walk methods²⁴, and Buffard et al. used this method to identify significant pathways from phosphoproteomic data from cancer cells²⁵. Scoring the topological proximity between disease proteins in the interactome was used to prioritize disease genes, infer new drug targets, identify the efficacy of drugs, and predict phenotypes^26–31. One challenge in network-based approaches is the incompleteness of the human interactome. It has many false negatives^32,33 and is biased to well-studied proteins. Previous works have considered incomplete interactomes, and different link prediction approaches are combined with network-based studies^34–37.

Many cancer drugs eventually lead to resistance and cause adverse effects when applied continuously and with high doses^38,39. The combination of drugs with different mechanisms of action is an alternative approach to improve the treatment strategies and repurpose drugs. Therefore, combinatorial drug treatment approaches are rigorously studied to eliminate or reduce resistance, recurrence, and possible side effects ^40–45. Effective drug combinations can be predicted with the help of network topology-based analysis of drug-disease interactions and inference of affected pathways⁴⁶.

All these reviewed approaches have the overarching aim of better understanding molecular alterations in diseases, how drugs act within the cell and finding the best treatment options. In this study, we further elaborate on perturbed networks to mechanistically understand the similarities and differences of drugs. Conceptually we (i) stratified the drugs at the network level beyond their immediate targets, (ii) evaluated the alterations of drug modulation in different cell lines (iii) suggested potential drug combinations based on topological separation of the networks. For this purpose, we used an integrative approach based on reverse engineering principles that combine a link prediction strategy to modify the underlying interactome and the solution of the prize-collecting Steiner forest problem to reconstruct drug and cell line-specific networks. All-pair comparison of the reconstructed networks shows that chemically and functionally different drugs may modulate shared pathways. We next considered these reconstructed networks for coherence with the available drug response data and possible drug combination prediction. Thus, these network models are rich resources for examining different aspects of drug actions at the pathway level in different cancer types.

Overview of the Method

We used transcriptomic and phosphoproteomic data of five cancer cell lines treated with 89 perturbagens and the associated control treatment (Connectivity Map Project - CMap) to understand the drugs' signaling level differences and commonalities systematically. We obtained the upstream regulators - the set of transcription factors - of the significantly expressed genes for each cell line-drug pair from transcriptomic data. Additionally, we retrieved the targets of each drug from CMAP Drug Repurposing tool⁴⁷, which combines the data from DrugBank, PubChem, and other drug-related databases. Finally, we merged the set of transcription factors, phosphoprotein hits, and drug targets to obtain the list of seed proteins of each cell line-drug pair for the network modeling. In Figure 1, we conceptually illustrated our integrative approach.

Although we started with 89 drugs and five cell lines, the number of proteins in the seed list is very low for some drugs, making the network modeling not feasible. We have also tested less stringent thresholds to enlarge the size of the seed list. However, we still could not overcome this issue for some drugs (i.e., decitabine, ginkgetin in A549 and YAPC, vemurafenib, and tacrolimus in MCF7 and PC3 cells). Therefore, we continued with 70 drug treatments on five cell lines. We calculated Tanimoto similarities and MACCS key distances of these 70 drugs (Supplementary Figure 1, Supplementary Note 1) to understand their chemical similarities and elucidated their mechanism of action (Supplementary Figure 2) to later use as a reference for the analysis.

To discover the hidden mechanisms underlying the effects of drugs, we integrated the seed proteins list (drug targets, transcription factors, and phosphoproteins) with the human interactome using Omics Integrator⁴⁸ software. The edge weights in the human interactome (iRefWeb v13.0) represent the confidence of the interactions and are calculated by the MI-Score function. We additionally refined the human interactome with several approaches to decrease the impact of false positive and false negative interactions. Additionally, the link prediction followed by cellular localization filters enriched the interactome with the missing interactions (see Methods). Omics Integrator solves the prize-collecting Steiner Forest (PCSF) problem to reconstruct optimal subnetworks by integrating the list of seed proteins (terminal set) and a reference interactome. The performance of the network modeling approaches is highly dependent on the interactome, parameter tuning, and inclusion of known biological insights. Omics Integrator has a better performance compared to all pairs shortest paths, page rank, and heat diffusion approaches when rigorous parameter tuning is applied and multiple suboptimal solutions are merged. This modification increases the precision and recall of PCSF⁴⁹. Additionally, PCSF can reconstruct a network from a single list of seed nodes whereas other powerful approaches such as PathLinker²⁰ and ResponseNet⁵⁰ require two sets of inputs. Because of its good performance among other methods requiring only a single seed list, Omics Integrator was selected as the main tool for network reconstruction. In this study, we merged the outputs of multiple parameter sets and refined the underlying interactome to reconstruct better networks. The terminal set size varies between 6 and 214 across the cell line-drug pairs, and the representation of phosphoproteins is relatively limited (max of 20 hits). Proteins have weights reflecting their importance in the terminal set based on their type. Transcription factors have the mean of their target gene expressions, phosphoproteins have the fold change compared to the control as the weight, and drug targets have a uniform weight inferred from the overall weight distribution of all proteins.

Overall, we constructed 236 subnetworks from the combination of 70 drugs and five cell lines. However, not all cell lines have the same number of subnetworks. Out of 236 networks, cell line A375 has 70, A549 has 46, MCF7 has 43, PC3 has 59, and YAPC has 18 drug-specific subnetworks (Supplementary Figure 3A). The topological characteristics of the drug networks can give clues about the perturbation or the connectivity of the modulated proteins. We summarized the topological properties of networks for each cell line (Supplementary Figure 4). We found that the final network size is dependent on the number of altered proteins/genes from multi-omic data rather than the centrality of the drug targets in the interactome. All topological analyses are present in the Supplementary Material (Supplementary Note 2, Supplementary Figure 5).

We initially quantified the node level overlap between network pairs where extensive overlap means similar network neighborhood. 98.4% of network pairs share at least one protein. Therefore, a comparison solely based on the overlapping nodes and node frequencies does not allow us to understand the comprehensive modulation of drugs (Supplementary Figure 3 and Supplementary Note 2).

The reconstructed networks preserve more detailed information about the drug effects beyond the common proteins and their association with cancer pathways. If two drug networks highly overlap, these drugs' action and phenotypic outcomes may be potentially similar. We applied the network-based separation approach, developed by Menche et al. ⁵¹, to reveal disease-disease relations on each drug network pair within and across cell lines to explore the network level separation. The separation score represents how close two networks are to each other, a topology-based comparison calculating the average shortest distances between the nodes in each network in the reference interactome (see Methods). The topological overlap of two drug networks in a cell line reflects their similarities at the pathway level. We found several overlapping subnetworks for drug pairs that do not have common target proteins or chemical structures, such as BIX-01338, entinostat, etoposide. The similarity between two networks represented by the separation scores is found statistically significant based on a hypergeometric test (Supplementary Data 1).

We also compared the application of the separation score method on reconstructed networks against the direct application on the seed proteins. Direct application of the method resulted in a distribution of separation scores mainly within the negative range. Because of the molecular heterogeneity, it is rarely expected for any drugs to modulate the same pathways in all cancer types. This situation implies false positives across cancer types. The difference between the two applications stems mainly from the intermediate (Steiner) nodes found via the network reconstruction method. These intermediate proteins can potentially reveal the off-target effects of drugs (Supplementary Figure 6 and Supplementary Note 3).

Moreover, the effect of link prediction applied on the interactomes used in network reconstruction is investigated to reveal whether the predicted edges cause any misleading results or not. We observed that predicted edges constitute very low percentages of total edge count in the networks, and they do not cause any spurious proteins to be included in the networks. On the other hand, the predicted edges can be considered as the noise introduced to the reference interactome. The very low level contribution of link prediction to the reconstructed networks show the robustness to this noise (Supplementary Note 4, Supplementary Table 1 and Supplementary Data 2).

Chemically and Functionally Different Drugs May Modulate Overlapping Networks

Transforming the networks into a matrix of separation scores (Supplementary Data 3) and clustering represents the overlaps of the networks as drug modules. We observed that the network overlap of some drugs is very high in specific cell lines such as A375, MCF7, and PC3, although their chemical similarity is limited. All-pair separation scores of drug modules in the A375 (skin melanoma) cell line are illustrated in Figure 2A with corresponding MoA similarities (T1=same MoA, T2=different MoA). We noticed that the chemically and functionally different drugs might modulate similar networks. Separation scores of network pairs in other cell lines can be found in Supplementary Figure 7. For example, we observed that tacedinaline and geldanamycin have overlapping networks in the A375 cell line, suggesting similar signaling output despite their different targets or mechanism of action (separation-score = -0.61, Figure 2B).

Tacedinaline is a selective HDAC1 inhibitor, while geldanamycin is a heat-shock protein (HSP) inhibitor targeting HSP90AB1 and HSP90AA1 proteins. Effects of HDAC inhibitor and HSP90 inhibitor drugs on HDACs and HSPs and the potential interplay between HDACs (especially HDAC6) and HSP90 has been highlighted in several studies^52–58. In our analysis, two drugs share several downstream transcription factors (i.e., HIF1A, HSF1, CEBPA, CEBPB, AR, FOXO1/3) in the reconstructed networks. We then linked these drug targets to cancer phenotypes using CancerGeneNet⁵⁹ which calculates the shortest paths between genes and phenotypes. We found that the downstream transcription factors of HDAC1 and HSP90s are shared on the path toward cell death and differentiation phenotypes (Figure 2B sub-panel). The intersection of tacedinaline and geldanamycin networks in A375 is enriched in several pathways such as MAPK, AMPK, and TGF-beta signaling pathways. Moreover, transcriptomic data of the two drugs are highly correlated. Despite the high overlap, they also have some disjoint pathways enriched in each specific network. Tacedinaline affects chemokine, neurotrophin, ErbB, TNF, FoxO, PPAR, prolactin, and thyroid hormone signaling pathways that are not present in the geldanamycin network (Figure 2C).

To systematically evaluate the performance of this method, we classified the drug pairs having the same MoA and similar resulting networks as true positives (TP) and drug pairs having different MoA and overlapping networks as false positives (FP). We defined drug pairs having the same MoA and separated networks as false negatives (FN) and drug pairs having different MoA and separated networks as true negatives (TN). We plotted the ROC curve and Precision-Recall curve using different threshold values (Supplementary Figure 8 and Supplementary Data 4). Because the data is imbalanced (200 positive cases - pairs with same MoA, 6017 negative cases - pairs with different MoA), the best performance is found based on Matthews Correlation Coefficient (MCC). The best performance is achieved when the separation score threshold is selected as -0.45. (precision=0.138, recall=0.400, accuracy=0.900, sensitivity=0.400, specificity=0.917, TPR= 0.400, FPR=0.083, F1-score=0.205, MCC=0.192). Same MoA and overlapping networks are expected, but the interesting cases are the drug pairs with different MoA and overlapping networks. We, therefore, further studied the literature to find evidence about mechanistic similarities among FP drug pairs. This evidence strengthens our claim that even drugs without the same MoA can perturb similar pathways (Supplementary Table 2). Despite being an FP based on the ground truth, the synergistic behavior of decitabine and entinostat was previously shown in the activation of FoxO1⁶⁰.

There are network pairs following the ground truth, i.e., drugs having the same MoA and overlapping networks. One example of drugs with the same MoA is sirolimus and OSI-027, which are selective mTOR inhibitors. Phase I clinical trials of OSI-027 were completed in 2013 for the investigation on patients with advanced solid tumors or lymphoma. It has activity on pancreatic ductal adenocarcinoma as an inhibitor of cell proliferation⁶¹. The networks of sirolimus and OSI-027 in A375 cells highly overlap (ss = -0.562). They have many common proteins and shared signaling pathways such as MAPK signaling, PI3K-Akt signaling, and cGMP-PKG signaling pathways. However, the OSI-027 network has specifically an enrichment of Hippo signaling, Jak-STAT signaling pathways, while sirolimus differs from OSI-027 with the enrichment of Wnt signaling, ErbB signaling pathways, and some immunity-related pathways such as Toll-like receptor signaling pathway (Supplementary Figure 9A-B).

Selumetinib and PD-0325901(Mirdametinib) are MEK inhibitors targeting MAP2K1. We showed a consistent similarity between two drugs for the same cancer types. The separation scores between two drugs reach up to -0.64 in YAPC (pancreatic cancer) (Supplementary Figure 9C-D). This similarity is a supportive result for the clinical studies. However, the network perturbation of these drugs may vary in different cancer types. For example, the separation score between PD-0325901 in the PC3 cell line and selumetinib in the YAPC cell line is 0.48. While the affected pathways are FoxO, VEGF, ErbB cAMP, Rap1, Ras, GnRH signaling pathways in PD-0325901 treated PC3, Jak-STAT, TGF-beta, PI3K-Akt, TNF, and Wnt signaling pathways are enriched in selumetinib treated YAPC (Supplementary Figure 9E-F).

Target selectivity of the drug pairs determines the network separation, despite having the same MoA.

The analysis of drug modulation in a cell line-specific manner helped determine a subset of drugs acting similarly in A375, MCF7, and PC3 (Supplementary Data 5). Afterward, we measured the separation of the drugs across different cell lines to discover if the same subset of drugs similarly acts in various tumor types. We compared the separation scores of drug pairs from 236 networks (Figure 3A) by dividing them into the same cell line and the different cell line groups. Next, we checked if two similar drugs in the same cell line act similarly in different cancer types. The results show that the network level overlap between drug pairs in the same cell line is preserved across different cell lines, but their separation is relatively higher. For example, resveratrol and geldanamycin are highly overlapping in the A375, MCF7, and PC3 cell lines (in A375, separation score=-0.60; in MCF7, separation score=-0.47; in PC3, separation score=-0.54), but their overlap is lower across cell lines such that separation score between A375-resveratrol and MCF7-geldanamycin is -0.40 and separation score between A375-resveratrol and PC3-geldanamycin is -0.43. Another aspect of this difference is the target selectivity of the drugs. Two examples of this phenomenon are the HDAC inhibitors and JAK inhibitors.

RGFP-966 is a slow-on/slow-off, competitive tight-binding selective HDAC3 inhibitor. Tacedinaline selectively targets HDAC1. Other HDAC inhibitors, Belinostat, Entinostat, Trichostatin-a, and Vorinostat, inhibit several HDACs. The network separations of these HDAC inhibitors across all cell lines vary based on the target selectivity of these drugs. We observed both overlapping and non-overlapping networks depending on the cell types. Belinostat, Trichostatin-a, and Vorinostat networks usually have some part overlapping in all cell types, as these three drugs are in the same structural group, and it is expected for them to behave similarly. Tacedinaline treatment activates similar mechanisms on A375, MCF7, and PC3 cells as separation scores between these are negative, but the network of Tacedinaline on YAPC cells is more distant than other cell networks. The separation score of Tacedinaline networks on A375 and YAPC cells is 0.38. YAPC-Tacedinaline network is a small network. It is not significantly enriched for any pathways other than Osteoclast differentiation which is recently shown to be induced in pancreatic cancer-derived exosomes ⁶². In contrast, the A375-Tacedinaline network is enriched for several critical signaling pathways such as MAPK, ErbB, FoxO, TGF-beta, and TNF signaling pathways (Figure 3B-C). However, RGFP-966 treatment on A375, MCF7, and PC3 results in distant networks both from each other and from the networks of other HDAC inhibitors. All pairwise separation scores of RGFP-966 networks are positive (MCF7-PC3:0.42, A375-PC3:0.28, and A375-MCF7:0.45), and the number of shared proteins is minimal. Selectivity of RGFP-966 on HDAC3 among all HDAC proteins results in the induction of various pathways in different cell types. When RGFP-966 is compared to Belinostat, Entinostat, Trichostatin-a, and Vorinostat, the most distant networks have separation scores around 0.4, and these networks generally belong to A375 and PC3 cells. For example, the PC3-RGFP-966 network is separated from A375-Trichostatin-a with a score of 0.46. Two networks are not commonly enriched in any pathway, while PC3-RGFP-966 is only enriched in cell cycle, but A375-Trichostatin-a is enriched in several pathways such as estrogen, neurotrophin, and TNF signaling pathways (Supplementary Figure 10A-B). As Tacedinaline is the other selective HDAC inhibitor, RGFP-966 and Tacedinaline treatments cause more distant PPI networks. RGFP-966 in PC3 and Tacedinaline in A375 have the most different networks with a separation score of 0.62. They only have three proteins in common, and there are no commonly affected pathways (Supplementary Figure 10C-D).

Among three JAK inhibitors, ruxolitinib, tofacitinib, and TG-101348 (fedratinib), only TG-101348 is known for its selectivity to JAK2 protein, and the rest targets JAK1, JAK2, JAK3, and TYK2 proteins with varying affinities^63–68. According to our results, the same drug may modulate different subnetworks in different cell lines, despite a marginal overlap. Within the same cell lines, tofacitinib and ruxolitinib networks have the most distant network pairs with separation scores of 0.68 in A375 and 0.46 in PC3. In A375, they only share the target proteins, JAK1 and JAK2. The tofacitinib network differs from ruxolitinib with the enrichment of Jak-STAT, ErbB, PI3K-Akt signaling pathways. In PC3, they share two more proteins other than JAK1/2 while the ruxolitinib network is enriched in the AMPK signaling pathway, and the tofacitinib network is only enriched in the Jak-STAT signaling pathway. A relatively low separation score (-0.147) is obtained from networks of TG-101348- and tofacitinib-treated A549 cells. Two networks share 23 proteins besides the target proteins, and commonly induced pathways include Jak-STAT and Notch signaling pathways. However, TG-101348 and tofacitinib networks have a separation score of 0.4 in A375, and the Jak-STAT signaling pathway is commonly enriched while the A375-tofacitinib network differs from TG-101348 with ErbB, mTOR, PI3K-Akt signaling pathways (Supplementary Figure 11).

Topologically separated drug networks may guide the use of drug combinations

The most effective strategy for finding drug combinations is to modulate the orthogonal sets of proteins or pathways simultaneously. The modulated pathways should be close enough to influence the disease module commonly. That can be described as targeting a disease through multiple signaling pathways^11,69. Therefore, networks play an essential role in identifying drug combinations. Cheng et al. showed the efficacy of drug combinations using drug targets and disease proteins by calculating their separation mapped on the interactome⁴⁶. Given two drugs, if each drug module overlaps with a disease module and these two drug modules are topologically separated, this is defined as complementary exposure⁴⁶. They demonstrated six classes of drug combination-disease interactions from FDA-approved or experimentally validated pairs and found one class correlates with therapeutic effects (complementary exposure). In our study, we found the intersection between our network models and the drug combinations provided in Cheng et al. Our reconstructed drug networks include many hidden proteins and modulated omic hits besides the drug targets. Therefore, the coverage of the networks is higher than the drug modules described in previous studies. In our analysis, we considered the disease module as the set of cancer driver genes (obtained from Cancer Genome Interpreter⁷⁰). We defined a rule where if two drug networks have at least one cancer driver protein in common and each has a topologically disjoint region, then the combination of these drugs can be effective. (Figure 4A). We only have four drug pairs intersecting with Cheng et al., but when we consider them per cell type, we could analyze seven drug pairs in total.

The separation scores of seven cases vary between -0.09 and 0.675. The network pairs in this analysis share a limited number of nodes (min three, max 16 proteins). KEGG pathway enrichment of these drug networks shows that each drug modulates a different region of the interactome so that different signaling pathways are perturbed, and their combination may be effective (Supplementary Figure 12-13).

Considering these seven cases and the concept of complementary exposure as a guideline, we rated other drug pairs in our dataset with similar criteria. If two networks are separated (ss>0.0) and share less than two common enriched pathways and the networks are not in small size (number of nodes > 40), and at least one is a large network (number of nodes > 100), then these drug pairs have potential to be used in combination. Given these criteria, we found 1,441 potential drug combinations out of 6,217 possible cell line–drug network pairs. Some of these new drug pairs were previously supported in the literature for their use in combination therapy. The top-ranking drug pair was I-BET-762 and lenalidomide in the A375 cell line. Two networks are separated with a score of 0.693, and they do not share any modulated signaling pathways. Lenalidomide network is enriched in several essential pathways such as PI3K-Akt, ErbB, MAPK, Ras, Rap1, VEGF, and FoxO signaling pathways, while the I-BET-762 network is enriched in Jak-STAT signaling and transcriptional misregulation in cancer pathways (Figure 4B). Three proteins are shared between these networks. Among them, "BRAF" is a cancer driver protein. These networks also contain different disease-related proteins (Figure 4C). The combination of lenalidomide with CPI-203, a primary amide analog of (+)-JQ1 and having the same mechanism of action as I-BET-762, is shown to synergistically induce cell death in bortezomib-resistant mantle cell lymphoma⁷¹. Similar activities of BET inhibitors and lenalidomide are reported in different studies of multiple myeloma^72,73.

Other examples within the top-ranking 100 predicted drug pairs are geldanamycin and tofacitinib (A375, ss= 0.656), trichostatin-a and tofacitinib (A375, ss= 0.562), and resveratrol and ginkgetin (PC3, ss= 0.53). HSP90 and JAK2 inhibition was shown to synergistically overcome resistance to JAK2-TKI in human myeloproliferative neoplasm cells⁷⁴. Tofacitinib targets JAK proteins nonspecifically; however, selective JAK2 inhibitor TG-101348 was also predicted in combination with tofacitinib. Moreover, HDAC and JAK dual inhibition were previously studied to improve treatment strategies^75–78. Resveratrol and ginkgetin combination was also shown to suppress VEGF-induced angiogenesis in colorectal cancer synergistically⁷⁹.

The separation between networks of a drug across cell lines gives clues about their sensitivity levels.

Next, we associated the network models with the drug sensitivity of the cell lines. For this purpose, we collected drug sensitivity scores of cell lines deposited in the Genomics of Drug Sensitivity in Cancer (GDSC) Database⁸⁰. The intersection of GDSC and CMap datasets results in five drugs that at least one cell line is significantly resistant or sensitive to. A375 is sensitive to CHIR-99021 and PD-0325901. CHIR-99021 has reconstructed networks in three cell lines (A375, MCF7, and PC3). PD-0325901 has reconstructed networks in four cell lines (A375, A549, PC3, and YAPC). Trichostatin-a has networks in two cell lines (MCF7 and A375), and MCF7 is significantly sensitive to Trichostatin-a. Thus, we can compare the networks to better understand the changes in sensitivity of cell lines to these drugs. PC3 is resistant to Staurosporine and YAPC is resistant to Dinaciclib. Since there may be several processes underlying the resistance of cells and we can not directly infer these mechanisms from our reconstructed networks, we focused only on the drugs that cell lines are sensitive to (CHIR-99021, PD-0325901, and Trichostatin-a).

If a drug has dissimilar network models across different cell lines, sensitivity to that drug may also vary (Figure 5A). Since Trichostatin-a(TSA) has highly overlapping networks across three cell lines, we could not observe this pattern for sensitive MCF7 across A375 and PC3. MCF7-TSA network differs from others with active TGF-beta signaling and cell cycle pathways, while A375 is differently enriched in MAPK signaling, neurotrophin signaling, estrogen signaling, TNF signaling pathways, and PC3 differs with Jak-STAT signaling and AMPK signaling pathways. TSA targets class I and class II HDACs. Treatment with HDAC inhibitors is reported to restore TGF-beta signaling in breast cancer⁸¹, so it is expected to observe that the TGF-beta signaling pathway is enriched in MCF7 cells (Figure 5B). The higher difference between z-scores implies the more separated networks in CHIR-99021 and PD-0325901 modulation. Therefore, analyzing the networks of these drugs and exploring enriched pathways may give clues about the resistance to these drugs. PD-0325901, to which A375 cells are sensitive, is a selective MAP2K1 (MEK1) inhibitor directly related to cell proliferation. Given the MEK is in the middle of the RAS/RAF/MEK/ERK pathway, it is on the upstream of several cellular mechanisms for cell proliferation and cell survival. In networks of two cell lines whose z-scores are close to the resistance threshold, YAPC and PC3, we observed several modulated pathways as expected. While the neurotrophin signaling pathway is enriched in all cells, PC3 cells differ with several pathways such as cGMP-PKG signaling, sphingolipid signaling, Rap1, and Ras signaling, VEGF signaling, Oxytocin signaling, FoxO signaling, and Fc epsilon RI signaling. Sphingolipid signaling is shown to be involved in the resistance of prostate cancer cell lines to antineoplastic drug treatment (z-score of PC3 is 1.95, very close to the resistance threshold)⁸² (Figure 5C).

We also analyzed all possible relationships between z-score differences and separation scores. We could find z-scores of 30 drugs and calculated pairwise z-score differences of those with cell line-drug networks and their associated separation scores. After filtering for pairs with one negative z-score and one positive z-score, and separation scores higher than -0.45, we performed a regression analysis and observed a moderately positive correlation such that separation score increases with increasing delta z-scores (R=.29, p-value=0.011) (Figure 5D).

Modulation of drug targets is not local when the complex interactions between molecules are considered within the cell. The impact can diffuse to distant parts or modulate several pathways simultaneously. Each drug has its specific modulated network, which can be used to detect the similarity between drugs. We, therefore, reconstructed the networks of cancer drugs by integrating transcriptomic, phosphoproteomic, and drug targets with protein-protein interactions using Omics Integrator. We modified the underlying reference interactome with the gene expression data and the predicted interactions to overcome its incompleteness and make it cell line-specific. The resulting drug networks cover the hidden proteins besides the initial data to explore connected mechanisms. We hypothesize that the topological overlap between each drug network pair on the interactome can be a solid basis for identifying their pharmacodynamic similarity. Our topological and pathway-level analysis of 236 reconstructed networks from five cell lines and 70 drugs demonstrated that chemically and functionally different drugs may modulate overlapping networks, which conventional comparisons cannot reveal. For example, tacedinaline and geldanamycin target different proteins, but their networks in the A375 cell line share a module consisting of common transcription factors leading to cell death and differentiation phenotypes. Another example is the drugs targeting the same proteins and modulating separated networks i.e., sirolimus and OSI-027. Sirolimus network is enriched in Wnt, cAMP, ErbB, and FoxO signaling pathways, while OSI-027 network is differently enriched in Jak-STAT and Hippo signaling pathways in A375 cells. Besides exploring the different drugs in the same tumor type, we can also compare the network-level impact of a drug across multiple cell lines. Tacedinaline, one of the HDAC inhibitors, modulates MAPK, ErbB, FoxO, TGF-beta pathways in the A375 cell line in the reconstructed network while these pathways are absent in the network of YAPC cell line. Additionally, the target selectivity of the drugs is a strong signature of the separation of the networks. We, therefore, specifically investigated the networks of HDAC inhibitors and JAK inhibitors. As an example, RGFP is the selective member of HDAC inhibitors. Networks of drugs, such as vorinostat, belinostat, trichostatin-a, targeting the same HDAC proteins have overlapping regions in almost all cancer cells. On the other hand, the network of RGFP-966, which targets HDAC3 selectively, is topologically distant from other HDAC inhibitor networks in different cancer types.

As evidenced by Menche et al. ⁵¹ that the overlap of the disease modules can quantify the molecular similarity of the diseases, we applied this quantification to interpret the drug similarities and eventually the effects of potential drug combinations. Treatment strategies combining drugs are used to decrease adverse effects and to prevent resistance to drugs^83–86. Therefore, we compared our drug network models with experimentally validated drug combinations. Recently, the DREAM AstraZeneca-Sanger drug combination prediction challenge provided a large combinatorial cell line screening dataset. The challenge was open to methods to predict synergistic drug combinations⁸⁷, i.e., machine learning or network-based approaches^88–91. In this study, we found that the experimentally synergistic drug pairs modulate topologically separated networks. Using their separation scores, common pathway enrichments, and network sizes, we rated cell line–drug network pairs and found 1,441 drug pairs out of 6,217 possible pairs. These predictions can be further used to experimentally study their synergy on specific conditions and provide guidance to finding alternatives to currently used drug pairs. In the future, we will integrate our network-level results with learning-based approaches to better understand the pharmacodynamics of drugs.

Another use of drug network models is to understand the drug resistance mechanisms. To reveal the differences between sensitive and resistant cells to a given drug, we compared the networks of which one is significantly sensitive, and the other is not to the given drug. We found that the networks of cell lines having different responses to the same drug have limited or no overlap. Therefore, pathway-level comparison and the topological analysis of drug networks provided the potential to guide us better understand drug response.

In summary, we leveraged networks to find pharmacodynamic similarities of drugs within and across multiple cancer cell lines. Further, these reconstructed networks are used for understanding drug response and combinations. We believe that network-based approaches are fundamental in interpreting the drug actions in tumors and integrating multi-omic data. This study is directly devoted to this task.

Data

The 'omic' data is obtained from ConnectivityMap (CMap), which contains 1.3 M L1000 data for nine different core cell types, treated with 27,927 perturbagens. Within CMap, the Touchstone V1.1 dataset contains 8,388 perturbagens (well-annotated genetic and small-molecular perturbagens) for nine cell lines for three time points and three replicates; and P100 dataset contains proteomic data for six different core cell types, treated with 90 perturbagens^92,93. In addition to transcriptional data (L1000), proteomic and phosphoproteomic data for a subset of these perturbagens have been recently released (P100, released on 4/13/18). P100 dataset was used in our study, which contains already filtered transcriptional and phosphoproteomic data for the same set.

The P100 dataset represents three different read-outs (phosphosignaling (P100), chromatin modifications (GCP), and transcriptional changes(L1000)) of the same 90 small-molecule perturbations in six cell lines. The duration of treatment was 3 hours for P100, 24 hours for GCP, and 6 hours for L1000. These data are available at multiple levels of processing: level 1 is fluorescence intensity (for L1000) or mass spectrometry extracted ion chromatogram traces (for P100, GCP); level 2 is gene expression or proteomic values without normalization; level 3 is normalized; and level 4 is differential (i.e., each sample is compared to all other samples on a plate). We used normalized (level3) data of L1000 and P100 for five cancer cell lines in our studies.

L1000 assay measures transcription levels of 978 genes (referred to as Landmark genes) and 80 control transcripts directly. It then infers the expression levels of the remaining 11,350 genes via Ordinary Least Square (OLS) regression. From 11,350 inferred genes, 9,196 genes are considered well inferred and called Best Inferred Genes (BING). We used landmark genes in the terminal set preparation procedure. Phosphoproteomic data generated by the P100 assay consists of 96 phosphosites and is a reduced representation of phosphoproteomics, targeting common signaling pathways.

The drug-target interactions are retrieved mainly from CLUE Drug Repurposing tool⁴⁷, which curates information from several databases such as DrugBank & PubChem (https://clue.io/repurposing-app). The reconstructed networks have been oriented through these targets for each drug.

We also have the edge-weighted protein-protein interaction network for the network modeling part retrieved from iRefWeb, v13.0 reference human interactome, which has 15,404 nodes (proteins) and 175,820 weighted edges (protein interactions) without self-loops. The weights of edges represent how real an interaction is based on the MIScore function.

Calculation of Tanimoto and MACCS key distance similarities

For 70 drugs listed with network reconstruction studies, SMILES signatures are collected from chemical databases such as PubChem and DrugBank. All pairwise Tanimoto and MACCS key similarities are calculated with python using the open-source RDkit⁹⁴ module. Hierarchical clustering on similarity matrices is performed with the Morpheus tool of BROAD institute with the parameters as the metric is 'euclidean distance' and linkage method is 'average' (https://software.broadinstitute.org/morpheus/ ).

Network Modeling Approach

For network reconstruction purposes, a solution to the prize-collecting Steiner tree problem was searched through Omics Integrator software¹⁹ to find an optimum tree. Nodes obtained from the experiments (terminal nodes) and nodes not detected in experiments and obtained by the algorithm (Steiner nodes) were determined by this process. Let G = (V, E, c, p) be an undirected graph, with the vertices/nodes V associated with non-negative prizes p(v), and with the edges E associated with non-negative costs c(e). The Prize-Collecting Steiner Tree problem (PCST) consists of finding a connected subgraph T = (V', E' ) of G, that minimize the weight of T, which is the sum of its edge costs plus the sum of the penalties of the vertices of G outside of the solution.

Omics Integrator is a software package that applies the prize-collecting Steiner forest (PCSF) approach to construct the most biologically relevant protein-protein interaction network. This tool is efficient in the integration of different omics data using interactome data.

There are two distinct tools within the Omics Integrator package: Garnet and Forest. To integrate different experimental results derived from mRNA, proteins, metabolites, etc., measurements, Garnet and Forest complement each other and use omics data either gathered from experiments or derived from several databases. In this study, only the Forest tool was used.

Forest tool constructs the interaction network by using user-defined omic data that is the list of proteins/genes together with their importance(prizes), and by using interactome data together with their significance levels. Each protein/gene given as input to the Forest tool is defined as 'terminal'. If necessary, Forest can add extra nodes from the interactome data called 'Steiner' nodes. When constructing the network, the algorithm optimizes the score by calculating the sum of prizes of nodes not included plus the costs of edges included. The algorithm seeks to minimize the score to find the most optimum and biologically relevant protein-protein interaction network. There is always a possibility to include 'hub' nodes, highly-connected nodes. Forest uses a generalized prize function to decide whether these 'hub' nodes are, in fact, essential and should be included in the network or not. This generalized prize function assigns negative weights to nodes according to the number of connections they have in the interactome. The function is:

(1)\(p\text{'} \left(v\right)= \beta \bullet p\left(v\right)- \mu \bullet degree\left(v\right)\)

where β and µ are scaling parameters and degree(v) is the number of connections of node v in the interactome. β is used to calibrate the effect of terminal nodes, and µ is used to calibrate the effect of hub nodes. µ is 0 in the default parameters where hub correction is disabled; if it is increased, the algorithm tries to exclude these hub nodes. β enables the algorithm to include terminal nodes, and increasing β facilitates more terminal nodes to be included in the final network.

Given a directed, partially directed, or undirected network G(V, E, c(e), p' (v)), it is aimed to find a forest F(V_F, E_F) that minimizes the objective function:

(2) \({f}^{\text{'}}\left(F\right)= {\sum }_{v\notin {V}_{F}}^{}{p}^{\text{'}}\left(v\right)+ {\sum }_{eϵ{E}_{F}}^{}c\left(e\right)+\) ω∙κ

where p' (v) is defined in Eq. 1, c(e) is the cost of each edge, κ is the number of trees in the forest, and ω is another scaling parameter which is a uniform edge cost of each node connected to dummy node (explained below).

The forest tool has six PCSF parameters; however, ω, β and D are the ones that at least need to be defined by the user. D is the depth parameter which is the maximum path length from v₀ to terminal nodes. The other three optional parameters are µ, g (reinforcement parameter, default is 1e-3), and garnetBeta (scales the Garnet output prizes relative to the provided protein prizes, default is 0.01).

Transcriptional and phosphoproteomic datasets (P100) were used to model protein-protein interaction networks. Gene expression data allows detecting transcription factors while phosphoproteomic data allows differentiating active and inactive proteins. By integrating these two data types and drug target data, multilayer networks were constructed. For parameter tuning purposes, we run the algorithm for each combination of D=10, µ=[0.00, 0.005, 0.01, 0.015, 0.02, 0.025, …, 1.0], β=[2, 3, 4, 5, 10] and ω=[1, 2, 3]. After collecting all the networks, we chose the ones with the highest number of terminal nodes and the minimum number of 'hub' proteins (degree>100) and then merged all of these networks to get the final cell line-drug network.

Preparation of Terminal Sets as an Input for Network Reconstruction Method

All statistical analyses were performed with the python scipy⁹⁵ module. L1000 data includes fluorescence values of several replicas for transcriptional measurements of both drug-treated samples and DMSO treated samples as control. These two conditions for each landmark (genes directly measured by L1000 assay) gene were compared with one way-Anova method after verifying the assumptions and p-values for each landmark gene are collected. The p-values for each gene as drug-treated condition compared against control condition is generated, genes that have lower p-value than the selected threshold were collected and labeled as 'significantly transcribed genes'. These 'significantly transcribed genes' are listed for each cell line-drug pair. Also, log2 fold changes of average transcription values of each gene compared to the control condition are calculated and stored to be later used in the prize designation of terminals.

From L1000 data, by using the 'significantly transcribed genes' list and transcription factor regulatory network⁹⁶, we found out transcription factors that are regulating any of the genes in our list. If a transcription factor regulates at least 3 of the genes in the significant gene list, we used it as a terminal, and mean log2 fold changes of interactors of these transcription factors are used as their prizes. From P100 data, 'significantly phosphorylated proteins' are collected using the phosphosites passing the p-value threshold(p<0.05).

Lastly, these two lists were joined together, prioritizing those coming from proteomic data. Also, the targets of the drug of interest were appended to the terminal set with a uniform weight inferred from the overall weight distribution of all proteins.

Application of link prediction strategy on the interactome

Link prediction is another step of the input file preparation. Firstly, interactome was processed to exclude hub proteins that have degree > 900 and more than ten standard deviations away from the mean ('UBC', 'APP', 'ELAVL1', 'SUMO2', 'CUL3') as defined in ⁹⁷. These hub proteins constitute 13,738 interactions in the reference interactome. All edges that include at least one of these hub proteins were excluded. After, the interactome was processed for each drug treatment in which the edges having at least one protein with transcription level, fluorescence measurement, below a threshold were excluded. L1000 data have expression values for each gene that is inferred from direct measurements of landmark genes. For each drug treatment, expression values lie between 0.0 and 15.0. The threshold was set as 2.0, and edges with at least one protein whose expression value is below 2.0 were excluded. iRefWeb interactome and L1000 genes share 11,002 genes in common. After excluding low expressed proteins, condition-specific interactomes included ~160,800 interactions.

After processing with expression values, the link prediction approach is applied. Adamic/Adar scoring method was selected because this link prediction method is known for its ability to weigh rarer features more heavily and exclude hubs. The scoring formula of Adamic/Adar is as below:

(3)\(Score(x,y)={\sum }_{\left(w ϵ Neighbors\right(x)\cap Neighbors(y\left)\right)}^{}\frac{1}{log\left|Neighbors\right(w\left)\right|}\)

Edge predictions were scored based on the processed interactome (interactome without hub proteins and low expressed genes). Then the same number of best scoring predictions as the edges found in the processed interactome were taken for further filtering based on subcellular location information. Localization information was gathered from the Human Protein Atlas. Using this information, predictions in which two proteins that do not have any common location were filtered out. The predictions in which at least one of the proteins do not have any localization information were collected since it is currently unknown and it has a possibility to be found in the same location with other protein and the predictions in which two proteins having at least one common location were also included in the artificial edge list. After filtering according to the subcellular localization, the remaining predicted edges were appended to the processed interactome by scaling their scores between 0.0 and 0.5. Their scaled scores were assigned as edge weight.

Comparison of Reconstructed Networks

Topology-based comparison

Menche et al., 2015 reported a pairwise network similarity scoring method based on mapping drug networks on the interactome and scoring their overlap⁵¹. They calculate the relative distances of two networks in the interactome. The separation score is measured based on the below formula:

(4)\(sAB\equiv \langle dAB \rangle - \frac{\langle dAA \rangle +\langle dBB \rangle }{2}\)

Where s_AB is the separation score of networks A and B, d_AB is the shortest distances between A-B proteins, and d_AA and d_BB are the shortest distances between proteins within A and B, respectively. If a protein is common between both networks (A and B), the distance d_AB of that protein is 0.

For our comparison purposes, this separation scoring method is used. Several separation score matrices were prepared to be able to compare drug networks on different levels. As there are different cell lines, there are matrices for each cell line that effects of drug treatments were compared within a cell type. Separation score matrices for each drug were also prepared to compare the drug's effect based on the cell type. Finally, a separation score matrix containing scores of pairwise network comparisons of all cell line-drug pair combinations was prepared to have a broader idea. All matrices are subjected to hierarchical clustering using python scipy.cluster.hierarchy module.

Pathway-based Comparison

Reconstructed networks are subjected to functional analyses based on their KEGG pathway enrichments. Pathway enrichments are calculated using DAVID source code⁹⁸. All significantly enriched pathways (p<.05) are collected for each network and signaling pathways are filtered to be later used for comparison purposes.

Calculation of topological features of the reconstructed networks

The number of nodes, number of edges, average degree, average shortest path lengths, density, and diameter of the networks are calculated with python using Networkx⁹⁹ module.

Hypergeometric Tests

Hypergeometric tests are used to analyze the significance of the network overlaps. For this purpose, hypergeom method of python scipy.stats module is used. The parameters of this method are x, M, n, N that are defined as the number of nodes common in both network, number of nodes in the larger network, number of nodes in the smaller network and number of total nodes found in the human interactome, respectively.

Data Availability

All datasets used in this work are publicly available from the following sources: The CMap data was downloaded using accession code GSE101406. Installation details and documentation for Omics Integrator software can be found on this link: https://github.com/fraenkel-lab/OmicsIntegrator. We used the PPI network from https://github.com/fraenkel-lab/OmicsIntegrator (IRefIndex Version 13.0), and drug-target data is curated from Cmap Drug Repurposing Tool (https://clue.io/repurposing-app ). The script for parameter tuning performed for network reconstruction can be downloaded from https://github.com/gungorbudak/forest-tuner. The separation score method is adapted from the source code provided in the supplementary material of the study of Menche et al.,2015⁵¹ (http://science.sciencemag.org/content/suppl/2015/02/18/347.6224.1257601.DC1 ).

Code Availability

All code used to conduct this study can be found from the GitHub repository on this link: https://github.com/metunetlab/IntegrativeNetworkModeling

Acknowledgements

N.T. has received support from UNESCO-L'Oreal National for Women in Science Fellowship, UNESCO-L'Oréal International Rising Talent Fellowship, and TUBA-GEBIP. S.U.B would like to thank Mahmut Beyge for his help on drug classification based on known mechanisms of action and valuable discussions on the challenges of drug classification methods.

Competing Interests

The authors declare that they have no competing interests.

Author Contributions

Conceptualization: SUB, NT

Data curation: SUB

Formal analysis: SUB, NT

Methodology: SUB, NT

Project administration: NT

Supervision: NT

Visualization: SUB.

Lunt SY, Vander Heiden MG. Aerobic glycolysis: meeting the metabolic requirements of cell proliferation. Annu Rev Cell Dev Biol. 2011;27:441–464. doi:10.1146/annurev-cellbio-092910-154237
Ramundo V, Giribaldi G, Aldieri E. Transforming Growth Factor-β and Oxidative Stress in Cancer: A Crosstalk in Driving Tumor Transformation. Cancers. 2021;13(12). doi:10.3390/cancers13123093
Zhan T, Rindtorff N, Boutros M. Wnt signaling in cancer. Oncogene. 2017;36(11):1461–1473. doi:10.1038/onc.2016.304
Liu X, Zhang Y, Wang Y, Yang M, Hong F, Yang S. Protein Phosphorylation in Cancer: Role of Nitric Oxide Signaling Pathway. Biomol. 2021;11(7). doi:10.3390/biom11071009
Bhat AA, Uppada S, Achkar IW, et al. Tight Junction Proteins and Signaling Pathways in Cancer and Inflammation: A Functional Crosstalk. Front Physiol. 2019;9:1942. https://www.frontiersin.org/article/10.3389/fphys.2018.01942.
Ward PS, Thompson CB. Metabolic Reprogramming: A Cancer Hallmark Even Warburg Did Not Anticipate. Cancer Cell. 2012;21(3):297–308. doi:https://doi.org/10.1016/j.ccr.2012.02.014
Hsu PP, Sabatini DM. Cancer Cell Metabolism: Warburg and Beyond. Cell. 2008;134(5):703–707. doi:https://doi.org/10.1016/j.cell.2008.08.021
Lee YT, Tan YJ, Oon CE. Molecular targeted therapy: Treating cancer with specificity. Eur J Pharmacol. 2018;834:188–196. doi:https://doi.org/10.1016/j.ejphar.2018.07.034
Panda AK, Chakraborty D, Sarkar I, Khan T, Sa G. New insights into therapeutic activity and anticancer properties of curcumin. J Exp Pharmacol. 2017;9:31–45. doi:10.2147/JEP.S70568
Clara JA, Monge C, Yang Y, Takebe N. Targeting signalling pathways and the immune microenvironment of cancer stem cells — a clinical update. Nat Rev Clin Oncol. 2020;17(4):204–232. doi:10.1038/s41571-019-0293-2
Dent P, Curiel DT, Fisher PB, Grant S. Synergistic combinations of signaling pathway inhibitors: mechanisms for improved cancer therapy. Drug Resist Updat Rev Comment Antimicrob Anticancer Chemother. 2009;12(3):65–73. doi:10.1016/j.drup.2009.03.001
Menche J, Sharma A, Kitsak M, et al. Uncovering disease-disease relationships through the incomplete interactome. Science (80-). 2015;347(6224):841. doi:10.1126/science.1257601
Guney E, Menche J, Vidal M, Barábasi AL. Network-based in silico drug efficacy screening. Nat Commun. 2016;7(May 2015):1–13. doi:10.1038/ncomms10331
Cheng F, Desai RJ, Handy DE, et al. Network-based approach to prediction and population-based validation of in silico drug repurposing. Nat Commun. 2018;9(1):1–12. doi:10.1038/s41467-018-05116-5
Irmisch A, Bonilla X, Chevrier S, et al. The Tumor Profiler Study: integrated, multi-omic, functional tumor profiling for clinical decision support. Cancer Cell. 2021;39(3):288–293. doi:https://doi.org/10.1016/j.ccell.2021.01.004
Turanli B, Karagoz K, Bidkhori G, et al. Multi-Omic Data Interpretation to Repurpose Subtype Specific Drug Candidates for Breast Cancer. Front Genet. 2019;10:420. https://www.frontiersin.org/article/10.3389/fgene.2019.00420.
Dhanasekaran R. Deciphering Tumor Heterogeneity in Hepatocellular Carcinoma (HCC)-Multi-Omic and Singulomic Approaches. Semin Liver Dis. 2021;41(1):9–18. doi:10.1055/s-0040-1722261
González-Reymúndez A, Vázquez AI. Multi-omic signatures identify pan-cancer classes of tumors beyond tissue of origin. Sci Rep. 2020;10(1):8341. doi:10.1038/s41598-020-65119-5
Ortega MA, Poirion O, Zhu X, et al. Using single-cell multiple omics approaches to resolve tumor heterogeneity. Clin Transl Med. 2017;6(1):46. doi:10.1186/s40169-017-0177-y
Ritz A, Poirel CL, Tegge AN, et al. Pathways on demand: automated reconstruction of human signaling networks. NPJ Syst Biol Appl. 2016;2:16002. doi:10.1038/npjsba.2016.2
Harris LA, Hogg JS, Tapia J-J, et al. BioNetGen 2.2: advances in rule-based modeling. Bioinformatics. 2016;32(21):3366–3368. doi:10.1093/bioinformatics/btw469
Wu Q, Finley SD. Mathematical Model Predicts Effective Strategies to Inhibit VEGF-eNOS Signaling. J Clin Med. 2020;9(5). doi:10.3390/jcm9051255
Halasz M, Kholodenko BN, Kolch W, Santra T. Integrating network reconstruction with mechanistic modeling to predict cancer therapies. Sci Signal. 2016;9(455). doi:10.1126/scisignal.aae0535
Naldi A, Larive RM, Czerwinska U, et al. Reconstruction and signal propagation analysis of the Syk signaling network in breast cancer cells. PLoS Comput Biol. 2017;13(3):1–27. doi:10.1371/journal.pcbi.1005432
Buffard M, Naldi A, Radulescu O, Coopman PJ, Larive RM, Freiss G. Network Reconstruction and Significant Pathway Extraction Using Phosphoproteomic Data from Cancer Cells. Proteomics. 2019;19(21-22). doi:10.1002/pmic.201800450
Mazza A, Klockmeier K, Wanker E, Sharan R. An integer programming framework for inferring disease complexes from network data. Bioinformatics. 2016;32(12):i271-i277. doi:10.1093/bioinformatics/btw263
Gottlieb A, Magger O, Berman I, Ruppin E, Sharan R. Principle: A tool for associating genes with diseases via network propagation. Bioinformatics. 2011;27(23):3325–3326. doi:10.1093/bioinformatics/btr584
Macropol K, Can T, Singh AK. RRW: Repeated random walks on genome-scale protein networks for local cluster discovery. BMC Bioinformatics. 2009;10(Mc):283. doi:10.1186/1471-2105-10-283
Erten S, Bebek G, Ewing RM, Koyutürk M. DADA: Degree-aware algorithms for network-based disease gene prioritization. BioData Min. 2011;4(1):19. doi:10.1186/1756-0381-4-19
Kim YA, Wuchty S, Przytycka TM. Identifying causal genes and dysregulated pathways in complex diseases. PLoS Comput Biol. 2011;7(3). doi:10.1371/journal.pcbi.1001095
Madeddu L, Stilo G, Velardi P. A feature-learning-based method for the disease-gene prediction problem. Int J Data Min Bioinform. 2020;24(1):16–37. doi:10.1504/IJDMB.2020.109502
Madeddu L, Stilo G, Velardi P. Network-based methods for disease-gene prediction. 2019. http://arxiv.org/abs/1902.10117.
Wang X, Gulbahce N, Yu H. Network-based methods for human disease gene prediction. Brief Funct Genomics. 2011;10(5):280–293. doi:10.1093/bfgp/elr024
Turki T, Wei Z. A link prediction approach to cancer drug sensitivity prediction. BMC Syst Biol. 2017;11(Suppl 5). doi:10.1186/s12918-017-0463-8
Kovács IA, Luck K, Spirohn K, et al. Network-based prediction of protein interactions. Nat Commun. 2019;10(1):1–8. doi:10.1038/s41467-019-09177-y
Cannistraci CV, Alanis-Lobato G, Ravasi T. From link-prediction in brain connectomes and protein interactomes to the local-community-paradigm in complex networks. Sci Rep. 2013;3. doi:10.1038/srep01613
Lu Y, Guo Y, Korhonen A. Link prediction in drug-target interactions network using similarity indices. BMC Bioinformatics. 2017;18(1):1–9. doi:10.1186/s12859-017-1460-z
Lippert TH, Ruoff H-J, Volm M. Intrinsic and acquired drug resistance in malignant tumors. The main reason for therapeutic failure. Arzneimittelforschung. 2008;58(6):261–264. doi:10.1055/s-0031-1296504
Kelderman S, Schumacher TNM, Haanen JBAG. Acquired and intrinsic resistance in cancer immunotherapy. Mol Oncol. 2014;8(6):1132–1139. doi:10.1016/j.molonc.2014.07.011
Bayat Mokhtari R, Homayouni TS, Baluch N, et al. Combination therapy in combating cancer. Oncotarget. 2017;8(23):38022–38043. doi:10.18632/oncotarget.16723
Hanahan D, Bergers G, Bergsland E. Less is more, regularly: metronomic dosing of cytotoxic drugs can target tumor angiogenesis in mice. J Clin Invest. 2000;105(8):1045–1047. doi:10.1172/JCI9872
Wang X, Zhang H, Chen X. Drug resistance and combating drug resistance in cancer. Cancer Drug Resist. 2019;2(2):141–160. doi:10.20517/cdr.2019.10
Sun X, Vilar S, Tatonetti NP. High-Throughput Methods for Combinatorial Drug Discovery. Sci Transl Med. 2013;5(205):205rv1 LP-205rv1. doi:10.1126/scitranslmed.3006667
Jia J, Zhu F, Ma X, Cao ZW, Li YX, Chen YZ. Mechanisms of drug combinations: interaction and network perspectives. Nat Rev Drug Discov. 2009;8(2):111–128. doi:10.1038/nrd2683
Ali MA, Rizvi S, Syed BA. Trends in the market for antihypertensive drugs. Nat Rev Drug Discov. 2017;16(5):309–310. doi:10.1038/nrd.2016.262
Cheng F, Kovács IA, Barabási AL. Network-based prediction of drug combinations. Nat Commun. 2019;10(1). doi:10.1038/s41467-019-09186-x
Corsello SM, Bittker JA, Liu Z, et al. The Drug Repurposing Hub: a next-generation drug library and information resource. Nat Med. 2017;23(4):405–408. doi:10.1038/nm.4306
Tuncbag N, Gosline SJC, Kedaigle A, Soltis AR, Gitter A, Fraenkel E. Network-Based Interpretation of Diverse High-Throughput Datasets through the Omics Integrator Software Package. PLoS Comput Biol. 2016;12(4). doi:10.1371/journal.pcbi.1004879
Kaan Arici NT. A Critical Assessment of the Performance of Network Reconstruction Approaches Based on Interactome Data Sets. Front Mol Biosci. 2021. doi:doi: 10.3389/fmolb.2021.666705
Lan A, Smoly IY, Rapaport G, Lindquist S, Fraenkel E, Yeger-Lotem E. ResponseNet: revealing signaling and regulatory networks linking genetic and transcriptomic screening data. Nucleic Acids Res. 2011;39(Web Server issue):W424-9. doi:10.1093/nar/gkr359
Menche J, Sharma A, Kitsak M, et al. Disease networks. Uncovering disease-disease relationships through the incomplete interactome. Science. 2015;347(6224):1257601. doi:10.1126/science.1257601
Krämer OH, Mahboobi S, Sellmer A. Drugging the HDAC6-HSP90 interplay in malignant cells. Trends Pharmacol Sci. 2014;35(10):501–509. doi:10.1016/j.tips.2014.08.001
Glozak MA, Sengupta N, Zhang X, Seto E. Acetylation and deacetylation of non-histone proteins. Gene. 2005;363:15–23. doi:https://doi.org/10.1016/j.gene.2005.09.010
Pinzi L, Benedetti R, Altucci L, Rastelli G. Design of Dual Inhibitors of Histone Deacetylase 6 and Heat Shock Protein 90. ACS Omega. 2020;5(20):11473–11480. doi:10.1021/acsomega.0c00559
Kovacs JJ, Murphy PJM, Gaillard S, et al. HDAC6 regulates Hsp90 acetylation and chaperone-dependent activation of glucocorticoid receptor. Mol Cell. 2005;18(5):601–607. doi:10.1016/j.molcel.2005.04.021
Chai RC, Vieusseux JL, Lang BJ, et al. Histone deacetylase activity mediates acquired resistance towards structurally diverse HSP90 inhibitors. Mol Oncol. 2017;11(5):567–583. doi:10.1002/1878-0261.12054
Li X, Yin S, Sakr W, Sarkar F, Heath E, Sheng S. Maspin inhibits HDAC1-mediated Hsp90 deacetylation. Cancer Res. 2008;68(9 Supplement):327 LP - 327. http://cancerres.aacrjournals.org/content/68/9_Supplement/327.abstract.
Zhang J, Xu E, Chen X. TAp73 protein stability is controlled by histone deacetylase 1 via regulation of Hsp90 chaperone function. J Biol Chem. 2013;288(11):7727–7737. doi:10.1074/jbc.M112.429522
Iannuccelli M, Micarelli E, Surdo P Lo, et al. CancerGeneNet: linking driver genes to cancer hallmarks. Nucleic Acids Res. 2020;48(D1):D416-D421. doi:10.1093/nar/gkz871
Wang C, Hamacher A, Petzsch P, et al. Combination of Decitabine and Entinostat Synergistically Inhibits Urothelial Bladder Cancer Cells via Activation of FoxO1. Cancers (Basel). 2020;12(2):337. doi:10.3390/cancers12020337
Zhi X, Chen W, Xue F, et al. OSI-027 inhibits pancreatic ductal adenocarcinoma cell proliferation and enhances the therapeutic effect of gemcitabine both in vitro and in vivo. Oncotarget. 2015;6(28):26230–26241. doi:10.18632/oncotarget.4579
Zhou Y, Zhu Y, Dong X, et al. Exosomes Derived from Pancreatic Cancer Cells Induce Osteoclast Differentiation Through the miR125a-5p/TNFRSF1B Pathway. Onco Targets Ther. 2021;14:2727–2739. doi:10.2147/OTT.S282319
Quintás-Cardama A, Vaddi K, Liu P, et al. Preclinical characterization of the selective JAK1/2 inhibitor INCB018424: therapeutic implications for the treatment of myeloproliferative neoplasms. Blood. 2010;115(15):3109–3117. doi:10.1182/blood-2009-04-214957
Verstovsek S, Kantarjian H, Mesa RA, et al. Safety and efficacy of INCB018424, a JAK1 and JAK2 inhibitor, in myelofibrosis. N Engl J Med. 2010;363(12):1117–1127. doi:10.1056/NEJMoa1002028
Clark JD, Flanagan ME, Telliez J-B. Discovery and development of Janus kinase (JAK) inhibitors for inflammatory diseases. J Med Chem. 2014;57(12):5023–5038. doi:10.1021/jm401490p
Changelian PS, Flanagan ME, Ball DJ, et al. Prevention of organ allograft rejection by a specific Janus kinase 3 inhibitor. Science. 2003;302(5646):875–878. doi:10.1126/science.1087061
Geyer HL, Tibes R, Mesa RA. JAK2 inhibitors and their impact in myeloproliferative neoplasms. Hematology. 2012;17 Suppl 1:S129-32. doi:10.1179/102453312X13336169156375
Wernig G, Kharas MG, Okabe R, et al. Efficacy of TG101348, a selective JAK2 inhibitor, in treatment of a murine model of JAK2V617F-induced polycythemia vera. Cancer Cell. 2008;13(4):311–320. doi:10.1016/j.ccr.2008.02.009
Zou J, Ji P, Zhao Y-L, et al. Neighbor communities in drug combination networks characterize synergistic effect. Mol Biosyst. 2012;8(12):3185–3196. doi:10.1039/c2mb25267h
Tamborero D, Rubio-Perez C, Deu-Pons J, et al. Cancer Genome Interpreter annotates the biological and clinical relevance of tumor alterations. Genome Med. 2018;10(1):25. doi:10.1186/s13073-018-0531-8
Moros A, Rodríguez V, Saborit-Villarroya I, et al. Synergistic antitumor activity of lenalidomide with the BET bromodomain inhibitor CPI203 in bortezomib-resistant mantle cell lymphoma. Leukemia. 2014;28(10):2049–2059. doi:10.1038/leu.2014.106
Siu KT, Eda H, Santo L, et al. Effect of the BET Inhibitor, Cpi-0610, Alone and in Combination with Lenalidomide in Multiple Myeloma. Blood. 2015;126(23):4255. doi:10.1182/blood.V126.23.4255.4255
Tania Díaz, Vanina Rodríguez, Ester Lozano, et al. The BET bromodomain inhibitor CPI203 improves lenalidomide and dexamethasone activity in in vitro and in vivo models of multiple myeloma by blockade of Ikaros and MYC signaling. Haematologica. 2017;102(10 SE-Articles):1776–1784. doi:10.3324/haematol.2017.164632
Fiskus W, Verstovsek S, Manshouri T, et al. Heat shock protein 90 inhibitor is synergistic with JAK2 inhibitor and overcomes resistance to JAK2-TKI in human myeloproliferative neoplasm cells. Clin Cancer Res. 2011;17(23):7347–7358. doi:10.1158/1078-0432.CCR-11-1541
Bae D, Choi Y, Lee J, et al. M-134, a novel HDAC6-selective inhibitor, markedly improved arthritic severity in a rodent model of rheumatoid arthritis when combined with tofacitinib. Pharmacol Reports. 2021;73(1):185–201. doi:10.1007/s43440-020-00188-x
Liang X, Liu H, Zhang Y. Novel-targeted therapy for hematological malignancies with JAK and HDAC dual inhibitors. Future Med Chem. 2019;11(15):1849–1852. doi:10.4155/fmc-2019-0168
Liang X, Zang J, Li X, et al. Discovery of Novel Janus Kinase (JAK) and Histone Deacetylase (HDAC) Dual Inhibitors for the Treatment of Hematological Malignancies. J Med Chem. 2019;62(8):3898–3923. doi:10.1021/acs.jmedchem.8b01597
Yao L, Ramanujulu PM, Poulsen A, Ohlson S, Dymock BW. Merging of ruxolitinib and vorinostat leads to highly potent inhibitors of JAK2 and histone deacetylase 6 (HDAC6). Bioorg Med Chem Lett. 2018;28(15):2636–2640. doi:10.1016/j.bmcl.2018.06.037
Hu W-H, Chan GK-L, Duan R, et al. Synergy of Ginkgetin and Resveratrol in Suppressing VEGF-Induced Angiogenesis: A Therapy in Treating Colorectal Cancer. Cancers (Basel). 2019;11(12). doi:10.3390/cancers11121828
Yang W, Soares J, Greninger P, et al. Genomics of Drug Sensitivity in Cancer (GDSC): a resource for therapeutic biomarker discovery in cancer cells. Nucleic Acids Res. 2013;41(D1):D955-D961. doi:10.1093/nar/gks1111
J Biol Chem. 2004;279(31):32620-32625. doi:10.1074/jbc.M402691200
Wang XZ, Beebe JR, Pwiti L, Bielawska A, Smyth MJ. Aberrant sphingolipid signaling is involved in the resistance of prostate cancer cell lines to chemotherapy. Cancer Res. 1999;59(22):5842–5848.
Sun X, Bao J, Shao Y. Mathematical Modeling of Therapy-induced Cancer Drug Resistance: Connecting Cancer Mechanisms to Population Survival Rates. Sci Rep. 2016;6(1):22498. doi:10.1038/srep22498
Johnson SF, Cruz C, Greifenberg AK, et al. CDK12 Inhibition Reverses De Novo and Acquired PARP Inhibitor Resistance in <em>BRCA</em> Wild-Type and Mutated Models of Triple-Negative Breast Cancer. Cell Rep. 2016;17(9):2367–2381. doi:10.1016/j.celrep.2016.10.077
Fabbro D, Manley PW, Jahnke W, et al. Inhibitors of the Abl kinase directed at either the ATP- or myristate-binding site. Biochim Biophys Acta - Proteins Proteomics. 2010;1804(3):454–462. doi:https://doi.org/10.1016/j.bbapap.2009.12.009
Griffin M, Scotto D, Josephs DH, et al. BRAF inhibitors: resistance and the promise of combination treatments for melanoma. Oncotarget; Vol 8, No 44. 2017. https://www.oncotarget.com/article/19836/text/.
Menden MP, Wang D, Guan Y, et al. A cancer pharmacogenomic screen powering crowd-sourced advancement of drug combination prediction. bioRxiv. January 2018:200451. doi:10.1101/200451
Menden MP, Wang D, Mason MJ, et al. Community assessment to advance computational prediction of cancer drug combinations in a pharmacogenomic screen. Nat Commun. 2019;10(1):2674. doi:10.1038/s41467-019-09799-2
Celebi R, Bear Don’t Walk O, Movva R, Alpsoy S, Dumontier M. In-silico Prediction of Synergistic Anti-Cancer Drug Combinations Using Multi-omics Data. Sci Rep. 2019;9(1):8949. doi:10.1038/s41598-019-45236-6
Li H, Li T, Quang D, Guan Y. Network Propagation Predicts Drug Synergy in Cancers. Cancer Res. 2018;78(18):5446 LP - 5457. doi:10.1158/0008-5472.CAN-18-0740
Sharma A, Rani R. An integrated framework for identification of effective and synergistic anti-cancer drug combinations. J Bioinform Comput Biol. 2018;16(05):1850017. doi:10.1142/S0219720018500178
Abelin JG, Patel J, Lu X, et al. Reduced-representation Phosphosignatures Measured by Quantitative Targeted MS Capture Cellular States and Enable Large-scale Comparison of Drug-induced Phenotypes * □. 2016:1622–1641. doi:10.1074/mcp.M116.058354
Lamb J, Crawford ED, Peck D, et al. The Connectivity Map: Using Gene-Expression Signatures to Connect Small Molecules, Genes, and Disease. Science (80-). 2006;313(5795):1929 LP - 1935. doi:10.1126/science.1132939
Landrum G. RDKit: Open-Source Cheminformatics Software. 2016. https://github.com/rdkit/rdkit/releases/tag/Release_2016_09_4.
Virtanen P, Gommers R, Oliphant TE, et al. {SciPy} 1.0: Fundamental Algorithms for Scientific Computing in Python. Nat Methods. 2020;17:261–272. doi:10.1038/s41592-019-0686-2
Liu Z, Wu C, Miao H, Wu H. Original article RegNetwork: an integrated database of transcriptional and post-transcriptional regulatory networks in human and mouse. 2015;(September). doi:10.1093/database/bav095
Hristov BH, Singh M. Network-Based Coverage of Mutational Profiles Report Network-Based Coverage of Mutational Profiles Reveals Cancer Genes. Cell Syst. 2017;5(3):221-229.e4. doi:10.1016/j.cels.2017.09.003
Huang DW, Sherman BT, Lempicki RA. Systematic and integrative analysis of large gene lists using DAVID bioinformatics resources. Nat Protoc. 2009;4(1):44–57. doi:10.1038/nprot.2008.211
Hagberg AA, Schult DA, Swart PJ. Exploring Network Structure, Dynamics, and Function using NetworkX. In: Varoquaux G, Vaught T, Millman J, eds. Proceedings of the 7th Python in Science Conference. Pasadena, CA USA; 2008:11-15.

(Not answered)

Download PDF

Journal Publication

published 19 Apr, 2022

Read the published version in npj Systems Biology and Applications →

Editorial decision: accept
21 Jan, 2022
Submission checks completed at journal
13 Jan, 2022
Unknown event
10 Jan, 2022
Editor assigned by journal
08 Jan, 2022
First submitted to journal
08 Jan, 2022

You are reading this latest preprint version

Functional Stratification of Cancer Drugs through Integrated Network Similarity

Status:

Journal Publication

Version 1

Abstract

Figures

Introduction

Results

Overview of the Method

Chemically and Functionally Different Drugs May Modulate Overlapping Networks

Topologically separated drug networks may guide the use of drug combinations

Discussion

Methods

Data

Calculation of Tanimoto and MACCS key distance similarities

Network Modeling Approach

Preparation of Terminal Sets as an Input for Network Reconstruction Method

Application of link prediction strategy on the interactome

Comparison of Reconstructed Networks

Topology-based comparison

Pathway-based Comparison

Calculation of topological features of the reconstructed networks

Hypergeometric Tests

Declarations

References

Additional Declarations

Supplementary Files

Status:

Journal Publication

Version 1