Network pharmacology strategies toward multi-target anticancer therapies: from computational models to experimental design principles.

Polypharmacology has emerged as novel means in drug discovery for improving treatment response in clinical use. However, to really capitalize on the polypharmacological effects of drugs, there is a critical need to better model and understand how the complex interactions between drugs and their cellular targets contribute to drug efficacy and possible side effects. Network graphs provide a convenient modeling framework for dealing with the fact that most drugs act on cellular systems through targeting multiple proteins both through on-target and off-target binding. Network pharmacology models aim at addressing questions such as how and where in the disease network should one target to inhibit disease phenotypes, such as cancer growth, ideally leading to therapies that are less vulnerable to drug resistance and side effects by means of attacking the disease network at the systems level through synergistic and synthetic lethal interactions. Since the exponentially increasing number of potential drug target combinations makes pure experimental approach quickly unfeasible, this review depicts a number of computational models and algorithms that can effectively reduce the search space for determining the most promising combinations for experimental evaluation. Such computational-experimental strategies are geared toward realizing the full potential of multi-target treatments in different disease phenotypes. Our specific focus is on system-level network approaches to polypharmacology designs in anticancer drug discovery, where we give representative examples of how network-centric modeling may offer systematic strategies toward better understanding and even predicting the phenotypic responses to multi-target therapies.

INTRODUCTION

Over the past decades, there has been a massive progress in many scientific and technological developments in pharmaceutical research; yet, over the same time period, the number of new drugs approved or successfully translated into clinical use has significantly declined, despite of massive investment on drug research and development by the global biotechnology and pharmaceutical industries [1]. This decline in pharmaceutical efficiency can be attributed to many factors, including lack of efficacy due to drug resistance and individual variation in treatment responses, as well as clinical safety or toxicology observed for the candidate drug compounds in pre-clinical or clinical studies in vivo [1-3]. It has also been increasingly understood that most drug molecules elicit their bioactivities by modulating multiple cellular targets and that such polypharmacological effects are behind many of the adverse side effects observed in clinical practice. However, polypharmacology can also be seen as part of the solution to the rather modest progress made so far in pursuing the expensive and suboptimal route of the current drug discovery. In particular, rather than trying to design selective ligands that target individual proteins only, polypharmacology aims to modify multiple cellular targets either by multi-target drugs or targeted drug combinations. Such multi-target treatments are being considered as a promising strategy to tackle the compensatory mechanisms and robustness of cellular systems, as well as to reduce unwanted off-target effects that often limit the clinical utility of many conventional drug treatments [3-5].

The potential of drug target combinations is perhaps best appreciated in anticancer research, where both genetic and non-genetic bypass mechanisms have led to inherent redundancy and robustness of compensatory signaling pathways in many cancer phenotypes [3; 6]. While most cancer cells show initial sensitivity to single-targeted drugs, their molecular heterogeneity often results in secondary outgrowth of new clones of rare cells that are resistant to the same therapy [7]. In contrast, combinations of drugs that target each individual clone and cancer escape pathway have the potential to kill cancers even at their advanced stages. However, one of the bottlenecks in the development of safe and effective anticancer drugs lies in the current inability to identify targeted compounds that will kill cancer cells at doses low enough to avoid severe side-effects. Therefore, much of the contemporary anticancer research is aimed at identifying specific genetic dependencies associated with cancer cells, with the hope of using such addictions or vulnerabilities to target directly cancer cells, while simultaneously reducing or even eliminating any unwanted side effects. However, large-scale cancer genome sequencing efforts have revealed tremendous mutational heterogeneity and clonal evolution, which renders it difficult to translate the genetic information into clinically actionable treatment strategies [8]. In particular, although cancer cells may harbor hundreds of genomic alterations in various biological pathways, only a subset of these alterations are driving the cancer initiation or progression in the different clones. Accordingly, the phenotypic response of even single drugs is often hard to predict because many of the compensatory cross-talk and feedback loops are still poorly understood in most cancer-related signaling pathways [9]. A global view of the interconnectivity of the signaling proteins and their functional contribution to cancer growth is therefore critical for the success of targeted single or multi-drug anticancer therapies.

Network graphs provide a convenient conceptual framework for system-level modeling, integrating and mining of high-throughput experimental datasets for understanding and gaining insights into different types of molecular relationship, such as how genes are linked to various diseases or interactions between drugs and their cellular targets. Systems or has proven useful for deciphering fundamental research questions, such as how perturbations in the cellular networks lead to certain phenotypes, including human diseases [9-16], whereas more recent field of aims at applying network modeling to tackle treatment-oriented drug discovery questions, such as where in the disease networks one should target in order to inhibit the disease phenotypes, for instance, cancer progression [9; 17-x003B; 17-23]. In such network context, the efficacy of multi-target therapy can be understood from a robustness of disease networks to deal with single node perturbations, due to inherent diversity and redundancy of compensatory signaling pathways that result in highly resilient network architecture with modular and interconnected topology [3]. Therefore, network-based drug discovery seeks for drug target combinations to perturb a specific subset of nodes in the disease-associated networks to inhibit the bypass mechanisms at systems level [4]. These recent developments have resulted in moving away from the traditional 'one target' strategy toward sub-network targets and systems or , a novel paradigm with potential to provide more global understanding of the mechanism behind drug action, resistance and side effects by considering drugs and their targets in the context of biological networks and interconnected pathways [2; 24-28]. However, even if sounding theoretically ideal, there is a need for efficient approaches to prioritizing most effective target combinations, while maintaining their safety and drug-like properties [2].

Rational and systematic design of network polypharmacology at systems level faces considerable challenges due to a number of experimental, modeling and computational challenges. Even with modern genome-level technologies for molecular profiling and drug screening, such as those based on next-generation sequencing or high-throughput screening of chemical compound libraries [29; 30], the number of possible drug and target combinations leads to a combinatorial explosion in the pharmacological and molecular spaces. Furthermore, the efficacy and side effects are also dependent on several other factors, such as the dose levels and time order of drug administration [31], which should be considered when predicting various drug response phenotypes. Therefore, the exponentially increasing number of possible combinations makes the pure experimental approach quickly unfeasible, and translates into a critical need for computationally efficient algorithms that can take into account the effect of partly overlapping target sets of promiscuous compounds to effectively reduce the search space for prioritizing most promising combinations for experimental evaluation. Ideally, given an input set of molecular or drug profiling measurements, the computational approaches should be able to model the dynamic information flow through cellular networks, starting from the individual ligand-receptor interactions and resulting in combined effects on biological sub-system or sub-networks [4; 32]. Moreover, there is an increasing evidence that non-genetic mechanisms may in some instances be even more important than the genetic factors when explaining, for instance, cancer evolutionary dynamics or intra-tumor heterogeneity [33; 34]. While more comprehensive computational modeling frameworks are still relatively far from the today’s reality and clinical practice, there are already a number of potential developments that are addressing specific steps in the modeling pipeline and information flow toward in silico prediction of drug treatment phenotypes in vitro or even in vivo [35-37].

Towards these goals, this review describes a number of computational-experimental systems pharmacology approaches, in which the concept of a network graph is a key component for searching potential drug target combinations, and how such approaches can be used especially for designing effective and safe anticancer treatments. As substantial efforts have recently been devoted to developing network-based methods for disease modeling and drug discovery, only representative examples of different approaches can be surveyed here, with an emphasis on methods related to concrete treatment applications, rather than more general network medicine approaches that have been reviewed elsewhere [11; 12; 14; 15; 17; 19; 20; 22]. Although our specific focus here is on global system-level modeling approaches, we will also briefly cover some related formulations, such as those based on differential equation or Boolean logic modeling, where typically only the model structure is represented in the form of a network graph, with connections describing either kinetic parameters or logic gates, respectively. There are also a number of recent excellent reviews and books that have touched upon some of these concepts, either in terms of the mechanistic or biological rationale for existing drug combinations, or from the experimental design or computational algorithm points of view [4; 38-43]. Rather than providing mere perspectives, however, our aim here is to focus on concrete network-centric approaches toward developing computational-experimental polypharmacology designs.

The remaining sections of the review follow a conceptual workflow of multi-target drug discovery process (see Fig. ()). We start by reviewing the current state of network analysis methods for in silico prediction of drug-target interactions. To understand the drug-target networks in the cellular context, systematic assessment of polypharmacological effects in the biological networks is next described and its implications for anticancer therapeutic design are highlighted. We further provide an overview of current resources and tools for integration of pharmacological and biological data, with representative examples of their successful applications to discovery of novel anticancer drug combinations. Finally, various models for scoring and optimization of potential drug combinations in the experimental validation phase are introduced as an integral part of the computational-experimental design strategy.

GLOBAL PREDICTION OF DRUG-TARGET INTERACTION NETWORKS

In the conventional paradigm for rational drug discovery, special emphasis is placed on the molecular mechanism of a particular disease to first pinpoint ‘druggable’ proteins, followed by finding or designing suitable drug compounds that interact with the desired target proteins. An alternative to this ‘target-based’ approach is so called ‘phenotype-based’ drug discovery, which starts from large-scale screening of chemical libraries to select those compounds that elicit desired phenotypes, and then tries to probe the targets of the candidate lead compounds in a specific disease setting to identify their (MoA) [44]. Regardless of the drug discovery paradigm, however, identification of cellular targets for the bioactive compounds is a common and essential step in the drug discovery pipeline, especially when going to multi-target designs. Global prediction of drug-drug and drug-target interaction networks provides a systematic means to explore the pharmacological space. In particular, network analyses at systems level allow us to link the pharmacological signatures of drugs to the molecular context of underlying cellular and tissue environment, from where the physiological consequences of drug perturbations can be elucidated from a more global perspective. This will hold great promise to the discovery of specific chemical compounds which modulate such cellular components that are therapeutically important in the disease processes, while better understanding and controlling their adverse effects prior to the actual clinical trials.

A rational design of anticancer treatments greatly benefits from quantitative proteome-wide characterization of drug-target interactions at the molecular level using both experimental and computational approaches. High-throughput experimental techniques for rapid evaluation of drug-target interactions have been developed, including cell-based phenotypic screens [45], binding affinity assays [46], in vitro modeling of ADMET (absorption, distribution, metabolism, excretion and toxicity) [47], activity-based probes [48] and transcriptional profiling [49]. However, the current screening capacities are still limited in testing the huge amount of possible compounds against hundreds of therapeutically relevant protein targets. Recent technological improvements especially in chemo-proteomic strategies are gradually improving the coverage and accuracy of the experimental drug-target mappings [50]. On the other hand, in silico target profiling methods have been considered as powerful complementary tool for providing efficient and systematic solutions to deduce relationships between chemicals and their cellular targets at the proteome scale. The predicted drug-target interactions can provide mechanics insights into drug MoA to facilitate decision making for designing selective lead compounds during experimental studies. In this section, we focus on recent computational techniques for global prediction of target-target, drug-target or drug-drug networks. These networks among drugs and targets can be constructed based on similarities defined either on the basis of the chemical structure of drugs, the phenotypic responses of drugs or the molecular information of the targets [51].

Ligand-based Prediction

Ligand is a substance (usually a small molecule) that binds to a site on a target protein (usually a receptor) to alter its chemical conformation leading to functional changes. Ligand-based approaches are based on the assumption that the bioactivity of a compound against a target can be predicted from other targets that are similar in their ligand chemistry. If two drugs share similar ligand-binding properties, then their targets may be related. Although this assumption might not hold true in specific cases, where a small change in ligand structure leads to a dramatic change of bioactivity [52], ligand-based prediction of drug targets has been successfully used in quantitative structure-activity relationship (QSAR) models, mainly for virtual screening of a large amount of compounds to identify those that have high probability of binding to a specific target [53]. Ligands are usually represented using chemical descriptors about their substructural fragments, such as 2D topological fingerprints or 3D pharmacophore structures. The metrics for the chemical similarity between a pair of ligands is often based on the conventional Tanimoto coefficient, which is defined by the fraction of descriptors shared between the ligand pair [52]. By utilizing efficient statistical correlation-based methods, network models that relate targets to each other can be built based on the selected similarity metrics. For instance, Keiser et al. recently developed a ligand set comparison method, named similarity ensemble approach (SEA) [54]. The SEA method is adapted from statistical techniques of the basic local alignment search tool (BLAST), which were originally designed for sequence similarity scoring [55]. More specifically, SEA calculates a raw similarity score by summing the Tanimoto coefficients of the ligand pairs across two sets that are known to bind to their targets. To access the significance of the ligand set similarity, the raw score is normalized to a Z-score and an expectation value (E-value) is derived by fitting the Z-scores of random sets to an extreme value distribution. The resulting matrix of SEA E-values among the ligand sets thus enables construction of a similarity network among targets either by sequential linkage of edges with significant E-values, or by constructing a minimum spanning tree to connect only the most similar neighbors [54].

It has been shown that the ligand-based target networks defined by SEA tend to connect the targets that are most likely to be inhibited by a common class of chemical compounds [56]. On the other hand, many targets that are related in biological functions also are clustered together in terms of their ligand chemistry. Applications of SEA to the ligand sets in the MDDR [57] and WOMBAT [58] databases have further revealed that the ligand-based chemoinformatics networks appear to be highly self-organized in a way that resemble networks of many natural phenomena and human activities [59]. This observation has implied that polypharmacology might be an intrinsic property of many chemical agents, in the sense that their ligands can be related to multiple targets through only a few connections (so-called ‘small-world’ network property). This polypharmacological nature of chemical compounds would offer great promises for discovering new therapeutic indications of existing drugs, a problem commonly referred to as [60]. Among some of the recent efforts to predict potential off-target affinities, Diller exploited pharmacological data stored in WOMBAT and BindingDB [61] to derive a ligand-based Shannon Entropy Descriptor (SHED) model that predicted the target profiles for a set of 767 drugs against a panel of 684 therapeutically-relevant targets [62]. The polypharmacological network was established by linking all targets sharing at least ten drugs. A fundamental discovery from such a target-target network analysis was that aminergic G protein-coupled receptors (GPCRs) appear as the most connected hub in the network, suggesting that the drugs targeting GPCRs tend to be the most promiscuous ones. The GPCRs hub is also closely related to opioid, sigma, NMDA and 5-HT3 receptors [62], providing polypharmacology information for antidepressants and antipsychotics drugs on their targets that were previously unknown to bind.

Target-based Prediction

The ligand sets collected for a given target provide an empirical description of potential drug-target interactions, where the information comes directly from compound libraries typically synthesized by combinatorial chemistry. Even though one might argue that ligand-based drug target discovery is more pharmacologically relevant, the mechanisms of drug-target binding are often difficult to understand without knowing the molecular information of the drug targets. With the rapid technical advances in molecular biology, the genome sequences coding for target proteins have been deciphered. Functional annotations of drug targets can therefore be made based on their amino-acid sequences or 3D structures. The sequence similarity reflects evolutionary relationship between ligands and associated proteins and thus provides a homology view of functional and structural links between proteins that are targeted by the same ligand [63]. In particular, recent advances in the determination of 3D structures of receptor proteins have revealed many factors that may contribute to drug target binding, such as the concavity of ligand binding sites and polar interactions [64]. Given that the structure and sequence-based information of target proteins is known, one can also predict the drug-target binding interactions by analyzing the molecular similarities of protein targets [65].

Statistical and machine learning approaches, such as Bayesian methods and kernel-based methods, have proven useful for binary classification of drug-target interactions on the basis of sequence and/or structure-related information [66-69]. In these models, proteins are usually encoded into a number of discrete vectors in terms of their biochemical and physiochemical features including, e.g. hydrophobicity, polarity and secondary structure. For example, Li et al. applied support vector machine (SVM) algorithms to classify target proteins using a 146-dimension vector of physicochemical features based solely on protein sequence information [69]. Along the same lines, Yamanishi et al. developed a bipartite graph model to map drugs and proteins into a unified pharmacological space, where the distance is determined by a kernel regression model [70]. The method was applied to integrating genomic sequence of proteins and chemical structures of ligands, with the aim to derive target interaction networks separately for four target classes including human involving enzymes, ion channels, GPCRs and nuclear receptors. A similar bipartite graph approach was also constructed for the FDA-approved drugs based on their known binary associations [71]. As an alternative approach, He et al. introduced a feature selection method called Maximum Relevance Minimum Redundancy (mRMR) algorithm, which ranks the sequence-based features according to their biochemical and physicochemical properties [72]. They showed an average success rate of 82.6% in the prediction of the drug-target benchmark data used by [70], albeit that many of predictions were not experimentally validated.

Mathematical modeling that utilizes 3D structures of targets, such as those based on reverse docking, have also become available [73; 74]. Reverse docking predicts those drugs that fit into the ligand binding site of a given target. Identification of ligand binding sites on proteins has recently become an area of tremendous interest, including a number of methods, such as pocket detection based on protein geometry, the energy contours on the protein surface, or the sequence and structure similarity with proteins with known functional sites [75]. Toward side effect prediction, Xie et al. showed that off-target binding can be also predicted on a proteome-wide scale using structure-based ligand binding site models [76].

Phenotype-based Prediction

Drug-target interactions often lead to phenotypic responses that can be measured using advanced pharmacological and molecular profiling methods. Information on drug target binding potency can thus be inferred indirectly by detecting the sets of proteins perturbed in response to the treatment. These phenotypic response profiles provide a new data layer that captures functional consequences of drug action beyond the conventional ligand-based or target-based information, and thus can help to expand our understanding of drugs’ MoA at a systems pharmacology level. The sour-ces of phenotypic measurements can be divided into several major classes: transcriptomics, proteomics, or side-effect phenotypes.

Transcriptomics-based methods assume that drug-target interactions are the primary causes of transcriptional changes in cells as a response to the drug perturbations. The underlying idea is to identify gene-expression signatures for drug treatments, and then the similarity between signatures can be related to predict new targets. A notable transcriptomics-based ongoing effort is the Connectivity Map (CMap), where the transcriptional profiles in five human cancer cell lines were originally profiled to assess the phenotypic changes brought about by more than 1000 bioactive small molecules [77]. Utilizing the CMap data, Iorio et al. developed a systematic approach to score the similarities of drugs based on their gene expression profiles and then they went on and constructed a drug-drug interaction network for the prediction of therapeutic and off-target effects [78]. Using this approach, they correctly predicted the MoA of 9 anticancer drugs, including HSP90 inhibitors and cyclin-dependent kinases (CDKs) inhibitors, which were not included in the original CMap dataset. A related approach was introduced by Hassane et al., where they combined the CMap data and the Gene expression omnibus (GEO) database [79] to construct a drug-disease network; they discovered two new compounds, celastrol and 4-hydroxy-2-nonenal, which produce a similar transcriptional response in acute myelogenous leukemia (AML) cancer cells to that of a known agent parthenolide (PTL) [80].

Proteomics-based approaches utilize drug response measurements at protein level, such as those based on bioactivity experiments using affinity chromatograph or protein expression measurements using protein microarrays [81]. Once the enzyme activities affected by drug treatment have been profiled, similar computational tools that are used for characterizing transcriptional responses of drug treatment can be applied. For example, Chen et al. built target-target interaction networks from the protein target similarities based on their binding affinity profiles computed from drug dose-responsive assays [82]. These bioactivity-based networks were overlapped with the ligand-based network determined by SEA to capture common connecting edges that are indicative of strong pairwise relations between targets captured by the two network perspectives. In addition to many promiscuous targets, the networks also revealed unexpected links between targets that exhibit similar binding properties while being unrelated in their biological functions.

Drug target identification can be also approached by side-effect phenotypic information. It has been observed that unrelated drugs can cause similar side effects due to their common off-target binding. Based on this concept, novel drug target associations can be found by comparing side effects of drugs [83]. In a seminal work, Campillos et al. exploited similarity in side effect to infer the drug target interactions and observed a relatively small overlap in the drug-drug relations predicted based on the side-effect similarity and chemical similarity, respectively [84]. This suggested that side effect phenotypes contain added complementary information which is not encoded by the drugs’ chemical similarity or the sequence similarity of their targets.

DISCOVERING DRUGGABLE TARGETS USING NETWORK APPROACHES

The drug-target networks identified using the computational tools, such as those reviewed in the previous section, addresses the drug-target interaction problem aiming at discovering new targets for given drugs or new drugs targeting given proteins. As the anticancer drugs’ MoA depends on many cellular factors and context, in-depth understanding of the effects of drug-target interactions need to be approached by the integration of disease network and drug-target interaction network to link the drug response phenotypes with the genetic makeup of the cancer cells [85]. High-throughput ‘-omics’ techniques have led to a better characterization of the global landscape of cancers, enabling e.g. a network-centric approaches to identification of druggable targets at the systems level. The integration of such ‘omics’ profiling with chemoinformatics to associate drugs, targets and disease outcomes has been a crucial step for studying polypharmacological effects. The challenge here is how to relate disease phenotypes and drugs in the cellular context by organizing both genotypes and phenotypes into a complex network through connecting proteins that are either involved in the disease development or serve as targets of drug treatments.

System-level identification of druggable targets requires two critical steps. The first step involves identification of causal relationship between genotypes and phenotypes that are implicated in the disease under analysis. The disease-associated genes might be either causal or non-causal, and also their contributions to the disease process may be dependent on multiple other factors, such as genetic interactions and protein interactions [8; 11]. For instance, prediction of genetic interactions that are implicated in cancer tumorigenesis has become an important area of anticancer research for discovering more effective drug targets based on the concept of synthetic lethality [8]. Pathway analysis can further convert the network interactions into functional sub-networks, reflecting cross-talk between biological processes and molecular-level mechanisms. In particular, metabolomics play an essential role in characterization of cancer progression and it has provided rich links into the metabolite changes after drug treatments. Network modeling for understanding and controlling of side effects is also being developed. The second step is to identify the drug-associated genotypes and phenotypes [44]. The challenge here is to integrate the chemoinformatics and bioinformatics analyses for rational strategy to pinpoint potential targets. An ultimate goal of such integrative analysis is to construct a drug-disease network, where the MoA of drugs and their biological consequences for human physiology and pathophysiology are explicitly represented.

Mathematical Modeling of Cell Signaling Pathways

Signaling pathways are particular types of canonical regulatory sub-networks that are involved in the transmission of cellular information about growth factors, nutrients and chemical perturbations. Signaling pathways play key roles in governing many cellular functions and coordinating cell actions. Knowledge of the defects in signaling transduction underlying cancer tumorigenesis and tumor growth can thus provide valuable information for discovering effective anticancer targets. For example, dys-regulation of protein tyrosine kinase signaling by mutations and other genetic alterations is known to lead to many malignancies [86]. The proteins involved in these dys-regulated pathways are often the intended targets of small-molecule drugs. A number of computational models have been introduced to characterize the functions of drug targets on signaling pathways, or predict drug effects on cellular phenotypes, such as migration and apoptosis from the signaling information [87]. The construction of a network model of signaling pathways can roughly be classified into three major model types: mass action-based, statistical association-based and Boolean logic-based. These modeling frameworks often capture causal protein interactions that can be predictive of cell behavior under a wide range of disease conditions, and thus can help in silico testing of potential drug combinations.

Mass-action based modeling typically begins with a collection of molecular interactions represented as ordinary or partial differential equations (ODEs and PDEs), with network topology describing the underlying kinetic rate parameters and their dynamic properties. A reliable construction of such detailed models often requires sensible selection of model parameters that are experimentally measured, as well as a control of model complexity to avoid overfitting to the data. Even after rigorous computational model validation, such as cross-validation, the model predictions should be subjected with experimental validation [88]. As an example, Birtwistle et al. developed a comprehensive ODE model that describes the ligand-dependent activation kinetics in both extracellular and cytoplasmic compartments, including a total of 117 species, 235 parameters and 96 reactions [89]. The mass-action model predicted the responses of two kinases, ERK and Akt, in the ErbB signaling pathways to the stimulation with EGF and HRG ligands in MCF-7 breast cancer cells. Mechanistic understanding of the ligand-dependent signaling is a crucial step for elucidating the ErbB network’s dys-regulation in many cancers. A similar study on the ErbB pathway was carried out, where two other ligands, TGFalpha and HRG, were modeled after their stimulations using H292 lung cancer cells [90]. These computational modeling frameworks of the ErbB signaling pathways may provide a critical assessment of the factors that influence drug efficacy, and thus have the potential to give insight into the effective combinations with ErbB-targeted drugs, such as ErbB1 kinase inhibitors erlotinib and gefitinib and monoclonal ErbB2 antibody, trastuzumab.

The logic-based models for signaling networks are usually represented as a regulatory graph where the response of the network is given by Boolean logic functions [91]. Nodes in the network represent components of the signaling pathways, such as receptors, ligands and transcription factors. Each node in the network takes value of 0 (inactive) or 1 (active), depending on the values of its regulators. The logic networks are usually structured into three layers: input, intermediate and output layers. Simulation of the network model starts with implementing a set of logical rules for the activation of each node, and the biological output of such a perturbation event can be derived by propagating input signals according to the logical connections in the network. The logical rules can be derived either from stoichiometric models or from the literature. Logic-based models are ideal for evaluating combinatorial effects of activation or inactivation, and thus provide insights into potential drug combinations that are targeting therapeutic components in cancer signaling pathways. Logic-based model can be applied to evaluate the phenotypic outcomes of knockdown of multiple protein components. For example, Sahin et al. employed a Boolean logic model to evaluate gene regulatory interactions, and predicted that combinatorial targeting of ERBB2 and EGFR may not inhibit the cancer growth in trastuzumab resistant breast cancer cells [92]. Instead, c-MYC was identified as a new potential drug target for breast cancer cells. A recent work by Morris et al. has extended Boolean logic network models to quantitative data by using constrained fuzzy transfer functions to provide quantitative relationships between the input and output components [93].

In comparison to the logic-based and mass-action-based models, association-based models for cell signaling are at the other end of the modeling spectrum, where data-driven statistical or machine learning methods are being applied. Network models in this category do not necessarily rely on any prior knowledge of molecular interactions among the entities; rather they utilize the observed correlations and other patterns in the experimental data when making inferences regarding network structure. Statistical approaches have widely been used in cancer research for identifying specific nodes in signaling networks that are amiable for therapeutic intervention. For instance, probabilistic graphical models (PGMs) provide a natural representation of signaling networks that can be used to also model cross-talk between the signaling components [94]. In this setting, the signaling components are considered as random variables and their values are dependent on the related components subject to uncertainty. Edges in the network can be either directional or non-directional, representing e.g. a measure of correlation or causality depending on the model assumptions. Due to the high dimensionality of omics datasets, the challenge has been in the learning of the underlying network topology as well as in reliable estimation of the model parameters [95]. For instance, Yörük et al. utilized dynamic Bayesian networks on protein array data to capture the time progression of protein signaling dynamics [96]. The network topology determined by the Bayesian networks can also be complemented with more detailed ODE-based approaches. In the application to the breast cancer cell line MDA-MB-468, the model discovered possible interaction between the MAPK and JAK/STAT pathways. Notably, the cross-talk was further validated by observing a reduced phosphorylation of STAT3p (S727) after treatment of an MEK inhibitor (MEKi), which inhibited MAPK as well.

System-level Metabolic Modeling of Drug Mode of Action

According to the systems biology view, most of the genetic components of complex disease phenotypes, such as cancer susceptibility, are not based on individual genes, but rather their interactions with other genes as well as with the environmental factors. In this context, the measurement of traits that are modulated but not encoded by the DNA sequence, commonly referred to as intermediate phenotypes, is of particular interest. One important class of intermediate phenotypes is metabolites. Metabolites are small-molecules present in a biological system, such as metabolic intermediates, hormones and other signaling molecules involved in breaking down and synthesis of human biological processes. Techniques for global metabolite profiling include liquid chromatography-mass spectrometry (LC-MS), gas chromatography-mass spectrometry (GC-MS), shot-gun mass spectrometry and nuclear magnetic resonance (NMR) [97]. Metabolomics as the systematic study of metabolic compounds has recently emerged as a powerful tool for the characterization of complex phenotypes, as well as for the identification of biomarkers for the onset and progression of cancers. During cancer development, cellular metabolism is often altered to adapt the requirements of excessive proliferation. Drugs also influence the metabolism of living cells by targeting enzymes which catalyze metabolic reactions where the metabolite concentrations and fluxes are affected. Therefore, metabolic networks lay down an important foundation for computational tools for disease modeling, and evaluation of metabolites that are perturbed by intended targets has suggested a potential way to drug discovery [98]. In particular, the reconstruction and analysis of system-level metabolic networks allow us to systematically explore the metabolic behavior of cells, and may thus increase the coverage of drug-target identification for cancers.

Reconstruction of metabolic networks usually starts by collecting from the integrated databases or literature the knowledge of numerous metabolic reactions that have been individually studied, and manually collecting the metabolic components in a bottom-up manner. Large-scale experimental mapping efforts for human metabolic networks are also under progress. One notable work is the first generic genome-wide human metabolic network, namely Human Recon 1, where stoichiometric interactions between metabolites were mapped and visualized [99]. An alternative is the Edinburgh Human Metabolic Network (EHMN), which has a higher coverage in the number of reactions than in the Human Recon 1 [100]. A third human metabolic network reconstruction is HumanCyc [101], which is derived computationally by matching the annotated human genome with the metabolic pathways provided in MetaCyc [102]. It was found that the consensus between the generic human metabolic networks was surprisingly low, especially at the reaction level, where about 18% of the reactions were commonly identified in a pairwise comparison [103]. The low level of overlap may be resolved by an integration effort, where useful complementary information can be discerned from inconsistent database biases, so that the accuracy and completeness of the human metabolic network can be improved. Tissue-specific human metabolic networks, such as the ones built for liver, kidney and brain using the INIT (Integrative Network Inference for Tissues) algorithm are also available, which enable the study of particular tissues under various genetic and physiological states [104]. The INIT algorithm starts by constructing the Human Metabolic Reaction (HMR) database, which incorporates existing genome-scale metabolic models, Recon1, EHMN and HumanCyc, as well as the KEGG database. Tissue-specific metabolic enzyme profiles were then obtained from the Human Protein Atlas (HPA) [105], where cell type specific high quality proteomic data are being stored. Using the INIT algorithm, the metabolic networks for 16 different cancer types were reconstructed, from which a number of such reporter metabolites were identified that are significantly enriched in the metabolism of cancer cells. Recently, another method was introduced that utilizes tissue-specific transcriptional profiles and metabolic network structure to infer genome-scale metabolic models for 126 human tissues and cell types [106]. The network reconstruction is based on ranking of gene-associated reactions according to a combination of expression-based and connectivity-based evidence using an algorithm called metabolic Context-specificity Assessed by Deterministic Reaction Evaluation (mCADRE). By comparing the networks built for cancer tissues and normal tissues, such metabolic reactions can be identified that are occurring significantly more often in tumor including those in the ecosanoid metabolism pathway.

Constrained-based models have recently gained considerably popularity for downstream analysis of the system-level network constructions for cancer specific metabolic features. Unlike mechanistic simulations of biological processes, which require measurements of kinetic rate parameters, constraint-based models aim at a steady-state characterization of the genome-scale metabolic networks by utilizing stoichiometric network connectivity structure. The models enforce physico-chemical constraints on biological networks including stoichiometric, mass balance and other constraints to define the space of allowed metabolic fluxes [107]. Metabolic fluxes are the intermediates of many biological signaling pathways in cancer physiology, where the cell genotypes and phenotypes meet. Constraint-based models, such as those based on flux balance analysis (FBA), can be used for identifying an optimal reaction flux distribution in order to achieve a biologically relevant objective function, such as ATP production or cell growth rate [108]. Depending on the application area, FBA can be classified as cell-specific, tissue-specific or context-specific. As the constraints are usually much easier to determine than the kinetic parameters, FBA has enabled near genome-scale models for several organisms. In human diseases, FBA-based modeling has typically been based on integration of genomic and bibliomic data. Using the Recon 1 as a template, a recent study further integrated a human metabolism modeling with gene expression data to predict flux distributions for metabolic reactions of NCI-60 cancer cell lines [109]. By reasoning that a drug’s target can catalyze its associated metabolic reactions, a Drug-Reaction Network (DRN) was constructed for evaluation of the drug influences on the metabolic network. FBA-based modeling has also been used in drug-target identification in cancer. For instance, Folger et al. reported the development of a genome-scale cancer metabolic network to predict characteristic alterations in cancer metabolism [110]. They first constructed a cancer metabolic model to capture major metabolic functions common across different cancer types. The model contained a total of 772 reactions and 683 genes. The FBA predictions on the network model were found to agree well with essential genes and mutations in cancers [110].

Network-based Modeling and Controlling for Side Effects

Designing effective and safe drug combinations is a fine tradeoff between optimizing drug efficacy and minimizing adverse side effects. It has been observed in many cases that complex interactions between on-target and off-target binding of chemical compounds is responsible not only for the side effects but is also critical for therapeutic efficacy of many drugs, including unintended off-target binding of statins and multi-target modulation of protein kinase inhibitors [76]. In anticancer applications, drugs with different toxic effects that attack cancer cells through distinct molecular mechanisms and cancer driving pathways should ideally be combined. However, our rudimentary knowledge of the MoA of many anticancer drugs is hindering such rational drug combination strategies [7]. Moreover, while intelligently selected synergistic drug combinations may allow a therapeutic effect to be achieved with lower doses of administered medicine, enthusiasm for this approach has been tempered by concerns that the therapeutic synergy of a combination will be accompanied by synergistic side effects if these are not properly taken into account in the design phase [5]. Therefore, global understanding of drug-target and drug-drug networks is essential not only for designing maximally effective drug combinations, but also for controlling their potential side effects. However, the efficacy and toxicity of multi-target drugs is not only attributable to complex interactions between drugs and their cellular targets, but also to many pharmacodynamic, pharmacokinetic, genetic, epigenetic, and environmental factors that affect the drug response phenotypes in vivo [4]. Given the great complexity of these interactions, it is understandable that the computational approaches to predicting compound toxicity — let alone adverse side effects in patients — is currently lagging behind the models for predicting drug efficacy [23].

Systems pharmacology is playing also an important role in understanding drug side effects by means of studying the adverse events as complex network responses [111]. As an example, more abstract network-based studies have focused on the global characterization of the correlation between interaction network topology and drug side effects. For instance, Brouwers et al. quantified the contribution of protein interaction network neighborhood on the observed side-effect similarity of drugs, and found out that drugs targeting proteins that are close in the network explain much less fraction of side-effect similarities, compared to side-effect similarities caused by overlapping drug targets; moreover, those targets that cause similar side-effects were more frequently in a linear part of the network, e.g. same pathway, than drug targets in general [112]. Similarly, Mizutani et al. observed in a comprehensive analysis of correlated sets of targeted proteins and side effects that most of these correlated sets were significantly enriched with proteins that are involved in the same biological pathways, even if their molecular functions were different [113]. Along the same lines, Wang et al. observed that while the close distance between the drug targets and disease genes seem to improve the efficacy of the targeted drugs, this may also increase incidence of side effects for drugs with too small distances; in particular, those drugs that have failed in clinical trials due to severe side effects showed smaller network distances than approved drugs [114]. These pathway and network topology characterizations not only provide systematic and even mechanistic interpretation regarding the relationship between drug-targeted proteins and known side effects, but may also be useful for predicting potential side effects of new multi-target drugs or drug combinations based on their protein-binding profiles and other information.

While screening candidate drugs for binding against every protein encoded in the human genome is not experimentally possible, computational approaches to predicting drug-target binding can facilitate also prioritization of follow-up mechanistic or pre-clinical studies for off-target binding or possible adverse effects. Most of these approaches are based on analyzing features in the chemical and molecular structures of drugs and their potential targets. Such structure-based approaches for global analysis of off-target binding have been successfully applied for specific drug chemicals [76]. For instance, Xie et al. used computational docking methods to predicted protein-ligand networks for a set of Cholesteryl Ester Transfer Protein (CETP) inhibitors, with or without known side effects, and suggested that adverse drug effects might be minimized by fine-tuning multiple off-target interactions using single or combinatorial therapies along multiple interconnected pathways [115]. Beyond studying specific therapeutic targets alone, certain structure-based approaches may be systematically applied to a wide range of drug classes. Recently, Lounkine et al. carried out a large-scale evaluation of their SEA approach, in terms of predicting the activity of 656 marketed drugs approved for human use on 73 unintended ‘side-effect’ targets [116]. Since SEA utilizes only chemical similarity in the target prediction, it can be applied systematically for all those targets that have known ligands. Based on a guilt-by-association metric that linked the targets to adverse drug reactions (ADRs), they further constructed a drug-target-ADR network, which may prove useful in prioritizing and streamlining the drug discovery process.

Genome-wide metabolic models are also providing a rich source of system-scale information for controlling potential side effects on a network-level. For instance, Li et al. formulated the problem of detecting optimal drug targets as an integer linear programming model, which finds such sets of targeted enzymes that provide maximal inhibition efficacy and minimal side effects originating from non-target compounds in the context of metabolic networks [117]. Similarly, Facchetti et al. developed an algorithmic solution, which uses genome-scale metabolic networks for systematic investigation of synergistic drug effects, and applied it to finding anticancer drug combinations with minimal side effects on the normal human metabolism [118]. Using a metabolism of the human kidney as a model system, Chang et al. evaluated metabolic drug response phenotypes in silico for a specific CETP inhibitor, torcetrapib, in the context of human renal function [119]. A number of causal drug off-targets were predicted to impact renal function, as well as genetic risk factors for drug treatment, which may play a role in the adverse side effects observed in clinical trials. These studies demonstrate possibilities of integration of structural and systems biology toward computational systems medicine strategies for personalized medicine. In many cases, drug candidates are found to be unsafe only late in the drug discovery process and clinical trials [120]. In fact, some adverse effects are not observed until a drug is on the market and widely used in genetically diverse populations. The field of pharmacogenetics focuses on the role of genetic factors in differential drug efficacy and toxicity in individuals with different genetic backgrounds; however, this field is out of the scope of the present review, and the readers are referred to recent excellent reviews [121; 122].

In general, prediction of drug off-targets and side effects can utilize many different sources of information about drugs and their targets, such as those from signaling, metabolic and protein interactions among targets. Moreover, it has also been shown that computational algorithms can relatively accurately predict side effects of a new drug, given the information on other drugs, such as their structural similarity and known side effects [123]. There are also dedicated databases, such as SIDER [124], which connect drugs to known side effects and ADRs. Towards integrative analyses, a number of recent computational approaches have combined a wide variety of molecular and pharmacological information to predicting drug combinations and interactions. For instance, Huang et al. showed that protein interaction network combined with drug structure can facilitate the prediction of ADR profiles over a number of side effect categories [125]. Further, Zhao et al. combined multiple data sources of drug features, such as their medical indications, molecular targets, toxicity profiles or anatomical therapeutic and chemical classifications to predict effective and safe drug combination [126]. Recently, Gottlieb et al. developed an algorithmic framework, called INDI, which calculates the likelihood that that a query drug pair interacts based on different drug-drug similarity measures, including those based on chemical similarity and registered or predicted side effects, as well as similarity measures constructed between drug targets, including sequence similarity, distance on a protein interaction network and Gene Ontology (GO) semantic similarity [127]. Notably, INDI is capable of handling both pharmacokinetic and pharmacodynamic interactions, and it was shown to provide relatively accurate prediction of adverse drug-drug interactions, as well as the severity level of co-administration of drugs used in the clinical practice.

DATABASE RESOURCES AND KNOWLEDGE DISCOVERY FOR NETWORK PHARMACOLOGY

Regardless of the approach taken, the common objective in multi-target drug design for anticancer treatments is to develop optimized computational-experimental approaches to systematically explore how drugs and their cellular targets interact to modulate cancer phenotypes on a global-network-level. The aim of such network pharmacology analyses is to identify molecular pathways behind drug action, as well as to identify key set of vulnerabilities in cancer networks and suggest effective and safe combinatorial treatment strategies that can block the cancer survival pathways. With the massive amounts of data from pharmacology and molecular biology, together with the development of public and proprietary resources, bioinformatics data integration tools have become crucial for data modeling and mining in a cost-effective manner. The first objective is to integrate multiple pharmacological databases to relate drug target profiles and the therapeutic effects in specific cancer types. The challenge here is that data acquisition from different research groups using different measurements need to be made comparable. A number of recent studies have revealed that drugs in general tend to bind to multiple targets involved in drug efficacy as well as in promiscuous off-target binding [4]. Furthermore, recent results of binding affinity assays of marketed drugs indicate that their therapeutic efficacy is not necessarily associated with high binding affinity only [46]. To predict the effect of a drug combination, it is essential to capture both high- and low-affinity binding of drug-target interactions on a proteome-wide scale.

Integrative drug-target interaction analysis and network construction can be done by combining data from comprehensive public data resources that focus on drug-target relationships. Experimental binding-affinity data can be retrieved from multiple databases, including BindingDB [61], ChEMBL [128] and canSAR [129] databases. DrugBank [130], Therapeutic Target Database [131], SuperTarget [132] and Matador [133] provide additional resources related to drug-target interactions by text mining the literature. The side effects of drugs are also accessible in some of these databases, while SIDER [124] focuses specifically on drug adverse effects. The PubChem Bioassay database is public information resource, consisting of bioactivity data generated by high-throughput screenings and chemical functional assays [134]. It also contains high-throughput siRNA screens targeting genes in human genome. A number of Web portals have also been developed for integrating multiple drug-target databases. For instance, STITCH is a public tool that integrates information about protein-ligand interactions from metabolic pathways, crystal structures, binding experimental data and text mining; it further provides a convenient network representation of protein-chemical relations for over 68,000 chemicals and the numbers are constantly increasing [135]. ChemProt integrates ChEMBL, DrugBank, PubChem, BindingDB and STITCH [136]; its current version includes more than 700,000 unique chemicals with bioactivities for 30,578 proteins, and more than 2 million chemical-protein interactions are represented in the context of protein-protein interaction networks [136]. Data mining can be implemented on the available databases for extracting and standardizing quantitative binding affinity data for anticancer drugs. Particularly, data normalization approaches are needed to allow direct comparison between the binding affinity from different measurement readouts [137].

The continuously growing ‘omics’ databases are commonly being used in drug discovery by deducing information on targeted pathways of drug treatments in terms of related changes in gene expression, protein abundance or metabolic concentration levels. The gene expression signatures of drug action can be retrieved, for example, from the CMap database [77]. One of the limitations is that the data in CMap is based on cancer cell lines and thus may not reflect the pathophysiology at the tissue or organism levels. Gene Expression Omnibus (GEO) is a generic public repository for gene expression profiles, where specific tissue types can be studied [79]. Protein interaction data can be extracted, for example, form MIPS [138], BIND [139] and PRIDE [140] databases. PINA is one of the most comprehensive efforts to integrate protein-protein interaction data from six databases including IntAct [141], MINT [142], BioGRID [143], DIP [144], HPRD [145] and MIPS [146], where a meta-database is provided with a set of web-based tools for network analysis and visualization [147; 148]. Knowledge about metabolic networks is commonly retrieved from KEGG [149] and BiGG [150]. Further integration of the data resources for heterogeneous phenotypes of drug response holds great promise for understanding drug MoA on a global scale. A remarkable examples of database integration are Chem2Bio2RDF [151] and PROMISCUOUS [152], where data of different types such as phenotypic data and drug side effects were linked to chemical data. These databases provide multiple data sources for understanding the relationships between drug-target interactions and the disease physiology from the perspective of the whole organism. Table summarizes a selection of currently available databases for integration of pharmacological and biological information, as well as some of the application areas in cancer polypharmacology.

EXPERIMENTAL TESTING OF THE POLYPHARMACOLOGY PREDICTIONS

It has been observed in large-scale experimental studies that synergistic drug combinations are relatively rare [160]. Therefore, computational models for prioritization of the most potential combinations for experimental testing can speed-up the drug discovery process and save experimental efforts, since many of the combinations in exhaustive testing would end up being negative hits. However, also in the experimental exploration of the potential drug combinations, which correspond to the selected multi-target modulations, one need to define appropriate models for synergy scoring and also efficient experimental strategies for exploring even the prioritized drug-dose combinatorial space.

Computer Aided Design for Drug Combination Screens

Computational search algorithms have successfully been applied in the past to optimize drug combinations prior to their experimental screening in the lab. For example, in smaller drug screening setups, including six drugs to lymphoma cancer, the search space of all possible drug combinations was represented graphically using a hierarchical tree, a special type of network graph [161]. The root level of the tree consisted of individual drugs and at the next levels the size of combination increased by adding one drug. The search for the best drug combinations was achieved by stack sequential algorithm aiming for a step-by-step optimal path through the tree that maximizes the effectiveness of a drug combination. Although showing good results in smaller drug setup, the computational complexity of this algorithm may become a limiting factor, as the maximal number of drugs that can be tested is only nine.

A similar effort to search for minimal drug combinations has been proposed using data from the NCI60 cell line collection [162]. Their goal was to determine a minimal set of drugs which can effectively treat all of the 60 cell lines of various cancer types. The method is mainly applicable to conventional cytotoxic chemotherapeutic drugs, as it considers the drug combination as a universal therapy that will kill both cancer and normal cells, which is seldom the purpose for targeted kinase inhibitors in personalized medicine applications. However, the method should be useful as a prioritization method, which narrows down the list of potential drug combinations, without further implication in their interaction patterns that may vary depending on the cancer type.

Other types of advanced algorithms have been proposed for prioritizing drug combinations in a cost-effective and timely manner using more complex cellular phenotypes. For instance, a multi-objective evolutionary algorithm was proposed for analyzing drug treatment-induced gene expression changes in a selected IL-1b cancer pathway [163], and a stochastic search algorithm has been developed for determination of optimal concentrations of drug combinations [164]. However, in the absence of link between the detailed pharmacological profiles of drugs and the underlying cellular contexts, the synergistic mechanisms behind these proposed combinations may be relatively difficult to understand.

Synergy Scores for Characterization of Drug Combination Effects

The drug perturbation effects are commonly measured using either cell proliferation or apoptosis phenotypes. These phenotypic outcomes can be quantified using simple measures, such as IC50 (defined as the dose of an antagonist that causes half-maximal inhibition) or EC50 (defined as the dose of an agonist that gives half maximal activation), or more recently introduced Activity Area, which corresponds to the area under the non-linear drug-dose response curve. In drug combinatorial screening, the drug synergistic effects are generally manifested either as potency shifts or efficacy boosts [39]. The commonly used synergy metrics are based on models such as Loewe additivity, Bliss independence or Highest single agents model [165]. When pairwise drug combinations are largely available for testing, Tan et al. used so-called S-score for each drug pair by averaging the interactions between replicates normalized by the variance of drug pairs as controls [166]. Similar scoring model has been used previously for detecting positive and negative genetic interactions in model organisms [167; 168].

While the assumptions behind these models are different, the antagonistic and agonistic synergy is generally interpreted as the discrepancy between the observed outcome and the effect that is expected in the cases when the individual drugs act independently (null hypothesis). The common feature about these synergy scores is that they measure the relative synergy compared to single agent effects. The synergy scores can be divided into those that do not take into account the nonlinearity in the dose response curves of the single agents (such as Bliss independence), and those that model also the joint dose-response curves of the combinations (such as Loewe additivity). These synergy metrics have been mainly designed for double drug combinations, but it is relatively straightforward to extend these also to higher combination orders. When the combination order is increased, one may also consider the differential synergy score to measure the gain by adding one more drug [169].

Visually perhaps more informative way to assess pairwise drug synergy is to look for response patterns in the three-dimensional dose-response surfaces where the mechanistic interactions between drugs can be visually compared [170]; however, applying such a method requires a detailed single agent dose-response matrix which may be experimentally difficult to obtain. Based on the Loewe additivity model, Cokol et al. developed a drug interaction score to quantify the concavity of the isophenotypic contours in a combinatorial dose response grid map, and they identified 38 novel drug synergies in the yeast S.cerevisiae which can be explained by genetic interactions on which they act [171]. A network representation of drug-drug interactions was also constructed from which a cluster of 6 drugs were found to be highly synergistic. Drug classes that have intrinsic tendencies towards synergy were also revealed. Moreover, the Loewe synergy score has also been extended to incorporate the structure information of the underlying signaling pathways related to the combinatorial drug targets [172; 173].

A recent study applied a synergy score based on a combination index by deriving a dose matrix for pairwise combinations for 13 experimental screens including anticancer screens [5]. The inhibitory effects of drug combinations were further illustrated in terms of their mechanisms using a FBA on E.coli metabolism models. A strong synergy between LY 294002 and camptothecin was found in H460 lung cancer cells, while not in Colo-205 colon cancer cells, suggesting that the interaction of drug combinations might be cancer-specific. Another study derived an analytical synergy score by first formulating a probabilistic model for the single-agent dose response curve [174]. The method was applied to chemical perturbation data of two drugs on 65 non-small lung cancer cell lines, where the synergy of EGFR/ERBB2 inhibitor BIBW-2992 and a PI3K/mTOR inhibitor PI-103 were found in EGFR and ERBB2 mutant cell lines.

In addition to applying synergy scores on static drug responses, Lopez et al. provided a model-based approach for evaluation of drug combinations using the time progression of tumor volume as drug response data [175]. The novelty here is that the in vivo time trajectories of tumor volume were explicitly represented in a differential equation growth model. The method also considered the differences in timing and order of drug administration. The method was applied to MCF-7 breast cancer cell line data, where the combinatorial effects were evaluated for two drugs doxorubicin and rhuMAb HER-2 that were given at different times and in different orders.

CONCLUSION

Although many current anticancer drug combination designs focus on targeted polypharmacology of a priori selected pathways or protein families, we argue that global network models may provide more comprehensive and unbiased insight into the underlying molecular mechanisms and pathways cross-talks behind cancer development, which are eventually required for the systematic identification of most critical pathways that go awry in disease and identification of optimal therapeutic strategies for controlling the dys-regulated sub-networks. Accordingly, although detailed models, such as those based on ODEs or Boolean logic formulations, have provided valuable quantitative or discrete insights into system behavior for predicting drug response phenotypes [163; 176-178], these modeling frameworks require rather detailed experimental characterization of the interaction kinetics and/or structure of disease-related key pathways, which are typically unavailable for many cancer phenotypes and may also bias the modeling results. In contrast, global models, based on e.g. metabolic network constructions or other global experimental measurements, should prove useful for predicting most effective multi-target drugs or drug combinations for anticancer therapies.

There are many potential directions how the experimental-computational approaches can be improved in the future. From the disease biology point of view, improved understanding of molecular vulnerabilities of cancer cells, in relation to normal cells, using concepts such as signal addictions or synthetic lethality, will likely to provide invaluable additional insights into the models of cancer progression and treatment. In theory, these concepts can address the fundamental challenges of anticancer therapy by optimally targeting differential features in each cancer type while sparing normal cells [8]. From the technology point of view, it seems necessary to develop more comprehensive drug screening panels and high-throughput screening platforms, in addition to more in-depth genomic and other ‘omics’ profiling technologies. Ideally, the experimental setup should be both economical and practical, utilizing large-scale functional measurements and such phenotypic readouts that are readily available in typical drug screening experiments. The experimental improvements, together with development of novel computational approaches for making most of these exciting measurements, should synergize the experimental-computational studies aiming at designing more effective and safe combinatorial drug treatments in the future.

While most current approaches to predicting cancer drug targets focus on large-scale profiling of genetic dependencies, the large number of genetic alterations present in tumor cells makes the discrimination of the cancer type specific driver mutations and pathways highly challenging. Even when genetic aberrations with pathogenetic importance can be identified, directly targeting these is often challenging. Furthermore, genes that are not altered at the genomic level may play essential roles in cancer development and treatment [21]. There is also substantial heterogeneity in response; even patients who share the targeted mutation and cancer type may show drastically different responses to the same treatment in vivo [9]. Therefore, functional screening of genes for their contribution to cancer progression has to go along with the structural characterization of the cancer genome to provide complementary insight into the molecular mechanisms and pathways behind various cancer types and to pinpoint the druggable drivers and other clinically actionable vulnerabilities as targets for personalized therapies. In particular, network-based approaches have the capacity to go toward predicting tumor-specific treatment responses [18; 21].

In more general terms, the emerging paradigms of network medicine and systems pharmacology have the potential to offer holistic information on disease networks and drug responses, enabling identification of more effective drug targets and their combinations tailored for safe and personalized cancer medicine. However, while such approaches to rational selection and evaluation of effective drug target combinations may offer possibilities to move beyond empirical, and often painstaking clinical trial and error, it should be appreciated that network medicine is still an emerging area of research and therefore all these methods should be considered as experimental. Model predictions needs to be carefully validated in experimental settings, something that has already been started in the context of metabolic network models in vitro [110]. However, it will likely take many years of extensive experimentation and iterations between computational and experimental cycles before fulfilling the promises of the systems pharmacology, especially in clinical settings in vivo. However, system-level computational-experimental approaches, such as those discussed here, could eventually offer improved means to facilitate the transition toward network pharmacology and systems medicine in the coming years.

Table 1.

Representative Examples of Database Resources and their Application to the Multi-Scale Modeling in Polypharmacology

Databases	URL	Description	Applications
ChEMBL [128]	http://www.ebi.ac.uk/chembl	A bioactivity database for over 1 million drug-like bioactive compounds and 5400 protein targets	Drug-target interaction predictions [112] Structure-activity relationships [47]
canSAR [129]	http://cansar.icr.ac.uk/	A repository of cancer specific biological data including gene expression, protein-protein interaction and RNAi screens together with chemical screening and pharmacological data	Identification of potential druggable targets from protein interactions [153] Polypharmacology map showing the shared compounds between queried targets [129]
STITCH [135]	http://stitch.embl.de/	A chemical-protein interaction database to query chemicals or proteins for their known and predicted relations using combined evidence from literature, experimental data and other databases	Benchmark for validation of in silico prediction of drug-target interactions [154]
PINA [148]	http://cbg.garvan.unsw.edu.au/pina/	An integrative platform collecting protein-protein interaction data from six manually curated public databases	Network construction, filtering and visualization for protein functional modules for six model organisms [147; 155; 156]
CMap [77]	http://www.broadinstitute.org/cmap/	A database of publicly available genome-wide gene expression profiles of five cancer cell lines in response to over 1300 bioactive small molecule treatments	Drug repurposing by linking drugs to each other or to diseases according to their gene expression signatures [157; 158]
BiGG [150]	http://bigg.ucsd.edu/	A knowledge-based reconstruction of genome-scale metabolic networks including human	Prediction of downstream effect of a drug perturbation in a disease network [109; 159]
SIDER [124]	http://sideeffects.embl.de	A database to connect marketed drugs to their recorded side effects and adverse drug reactions obtained from public resources using text mining	Linking side effects to drug-target interactions and pathways [123; 126]