Systems biology approaches for advancing the discovery of effective drug combinations
Ryall and Tan Journal of Cheminformatics
Systems biology approaches for advancing the discovery of effective drug combinations
Karen A Ryall 1
Aik Choon Tan 0 1 2
0 Department of Biostatistics and Informatics, Colorado School of Public Health, University of Colorado Anschutz Medical Campus , Aurora, CO , USA
1 Translational Bioinformatics and Cancer Systems Biology Laboratory, Division of Medical Oncology, Department of Medicine, School of Medicine, University of Colorado Anschutz Medical Campus , 12801 E.17th Ave., L18-8116, Aurora, CO 80045 , USA
2 Department of Computer Science and Engineering, Korea University , Seoul , South Korea
Complex diseases like cancer are regulated by large, interconnected networks with many pathways affecting cell proliferation, invasion, and drug resistance. However, current cancer therapy predominantly relies on the reductionist approach of one gene-one disease. Combinations of drugs may overcome drug resistance by limiting mutations and induction of escape pathways, but given the enormous number of possible drug combinations, strategies to reduce the search space and prioritize experiments are needed. In this review, we focus on the use of computational modeling, bioinformatics and high-throughput experimental methods for discovery of drug combinations. We highlight cutting-edge systems approaches, including large-scale modeling of cell signaling networks, network motif analysis, statistical association-based models, identifying correlations in gene signatures, functional genomics, and high-throughput combination screens. We also present a list of publicly available data and resources to aid in discovery of drug combinations. Integration of these systems approaches will enable faster discovery and translation of clinically relevant drug combinations.
Drug combinations; Systems biology; Computational modeling; Cancer; Drug discovery
Despite increasing investments in pharmaceutical
research and development, the rate of introduction of
successfully translated drugs has decreased . Reasons for
increased attrition rates in drug development include
toxicity and inadequate efficacy due to individual
variation in therapeutic response and development of drug
resistance . This increased attrition rate has coincided
with increased interest in seeking highly specific ligands
affecting single targets for treatment of disease .
Pharmaceutical research has increasingly relied on
reductionist approaches, even though systemic diseases
such as cancer and heart disease are managed by large,
interconnected networks with many pathways affecting
pathological signaling [4,5]. The redundancy and
feedback in these networks allows for robustness of phenotype
and maintenance of homeostasis [6,7]. This network
complexity has hindered development of new therapies and
indicates a need for more integrative systems approaches to
make better predictions of drug responses [8,9].
The failure of single targets to successfully translate
into clinical practice and the problem of development of
drug resistance with single target cancer therapies has
increased interest in discovery of effective drug
combinations. Administering drug combinations has been effective
in overcoming resistance to anti-microbial therapies for
treatment of infectious diseases such as HIV and
tuberculosis . In cancer, drug resistance can occur through
mutation of the drug target , amplification of an
alternate pathway , or intrinsic resistance of a subset of the
cancer cells . Combinations of drugs could potentially
overcome these resistance strategies by limiting the
potential of escape mutations and pathways .
While combination therapies may dramatically improve
efficacy of cancer therapies, the discovery of effective
combinations is a challenging endeavor. With over 1,500 FDA
approved compounds, experimentally testing every
possible combination of these drugs would be unfeasible, even
with high-throughput experimental methods .
Therefore, new systems approaches are needed to reduce the
search space and prioritize combinations for experimental
testing (Figure 1) .
Figure 1 Diagram depicting estimated ratio of computational and experimental requirements for various methods in this review.
For example, mass-action/kinetic modeling has higher experimental requirements than logic-based and normalized-Hill-based modeling due to
its need for many abundance and rate parameters. Unbiased high-throughput screening of drug combinations has the highest experimental
requirement. Many of the systems biology methods in this review aim to use publicly available data and computational approaches to reduce
the need for exhaustive screens and prioritize combinations for experimental validation.
Here, we review computational and experimental methods
for accelerating the discovery of effective drug
combinations for complex diseases, with special focus on cancer.
In addition, we include a list of publicly available
resources as a reference for future drug combination
studies (Table 1).
Quantifying synergistic drug combinations
When presenting the results of drug combination studies,
it is important to have a standard to statistically define
synergistic drug pairs. Two commonly used methods for
quantifying synergy between drug combinations are Loewe
additivity and Bliss independence. Loewe additivity is
based on the assumption that the two inhibitors act
through a similar mechanism while Bliss independence
assumes independent mechanisms .
Using Loewe additivity, the concentration of two
inhibitors (A and B) which alone results in X% inhibition of
the target ([IA]X%, [IB]X%) can be used to calculate the
theoretical concentrations of each inhibitor needed to
achieve the same X% inhibition when combined ([CA]X%,
The Loewe additivity applies the isobologram analysis
to evaluate the combination effects of two drugs at a
given effect. For example, in a Cartesian coordinate plot
where x and y-axes represent concentrations of drugs A
and B to achieve a defined effect X% (e.g., X = 50% for
half maximal inhibitory concentration (IC50) of [IA]50%
and [IB]50%), respectively. The coordinates ([IA]50%,0) and
(0, [IB]50%) represent the concentration for drugs A and
B, respectively. The line of additivity is constructed by
connecting these two points for a 50% effect
isobologram plot. The concentrations of the two drugs used in
combination to provide the same effect X% (e.g. X =
50%) will be denoted by point ([CA]50%,[CB]50%) and are
placed in the same plot. Synergy, additivity, or antagonism
will be determined when this point ([CA]50%,[CB]50%) is
located below, on, or above the line, respectively. More
generally, linear, concave, and convex isoboles represent
non-interacting, synergy, and antagonistic drug
combination, respectively (Figure 2A).
This approach led to the development of the
combination index (CI) popularized by Chou and Talalay .
Here, the CI provides a quantitative measure of the
extent of drug interaction at a given effect. It measures the
combination concentrations of drugs A and B to
produce a effect X%, [CA] and [CB], normalized by their
corresponding concentrations that produces the same
effect as a single agent, [IA] and [IB], respectively. CI value
is calculated by:
where CI value <1, =1, and >1 represent synergy,
additivity, and antagonism, respectively.
Bliss independence is based on probability theory and
assumes the two inhibitors are working through
independent mechanisms . The inhibitors do not
interfere with each other, but contribute to a common result.
Unlike Loewe additivity, calculating Bliss independence
does not require determination of doseresponse curves
for the individual compounds to determine the theoretical
Table 1 Resources for systems analysis of drug combinations
Database of biological activities of millions
of small molecules.
Database of target, chemical, pharmacological,
and interaction data for 7739 drugs.
Chemical-protein interaction database containing
300,000 small molecules and 2.6 million
proteins from 1133 organisms.
Database of adverse drug reactions from
Comparative Toxicogenomics Manually curated database of over a million interactions
Database (CTD) between chemicals and genes and over 1.6 million
associations between chemicals and diseases and over
15 million associations between genes and diseases.
Drug Gene Interaction
Connectivity Map (CMap)
Database of drug information including dosing guidelines, www.pharmgkb.org
drug labels, signaling pathway diagrams, drug-gene
associations, and drug-phenotype relationships.
Database and web tool for mining over 14,000
Data from 1363 drug combinations.
Database of over 720,000 protein and genetic interactions http://thebiogrid.org/
from model organisms and humans from over
Database of known and predicted protein interactions,
including both direct and functional associations.
It currently covers 5,214,234 proteins from
Gene expression profiles from 1309 FDA
approved small molecules tested in 5
human cell lines.
Public repository of gene expression data.
Web tool that identifies kinase inhibitors
for a set of query kinases.
Pathway database with visual representation
for 21 organisms, which includes over
1500 human pathways.
Large collection of manually drawn pathway
maps of molecular interaction networks for
various biological processes.
Open source software platform for network
analysis and visualization.
Modeling and simulation tool for construction of
normalized-Hill models of signaling networks
from user defined species interactions.
Free software for creating logic-based models
of signaling networks.
Repository of computational models of biological
processes. Includes both peer-reviewed models
and models produced automatically using pathway
resources like KEGG.
Detailed genetic characterization
of ~1000 cancer cell lines.
Table 1 Resources for systems analysis of drug combinations (Continued)
Genomics of Drug Sensitivity
in Cancer (GDSC)
Drug sensitivity data from hundreds of
genetically characterized cancer cell
lines perturbed with a wide
variety of anti-cancer agents Part of an ongoing
project to discover therapeutic biomarkers.
Drug screen data from a diverse panel of 60
human cancer cell lines with extensive
Drug sensitivity data of 242 genetically
characterized cancer cell lines treated with
354 different small molecule probes and drugs.
Each compound selectively targets a distinct
part of cell wiring and collectively affect a vast
array of cell processes.
results, making it easier to compute . Bliss
independence models the combined effect (ET) as the product of
the individual effects with drugs A (EA) and B (EB). The
predicted combined effect (ET) is computed by:
where each effect (E) is expressed as fractional activity
compared to control between 0 (100% inhibition) and 1
(0% inhibition). For example, if drug A and drug B each
result in 40% tumor growth compared to control, then
the predicted tumor growth when combined would be
(0.4)*(0.4)*(100%) = 16% of control according to Bliss
Independence. The predicted combined inhibition level
would therefore be 100%-16% = 84% inhibition of tumor
growth. If the actual tumor growth when drug A and B
are combined is less than 16% of control (greater than
84% growth inhibition), then the compounds would be
synergistic by Bliss Independence. If the tumor growth
level is greater than 16% of control (less than 84%
growth inhibition), then the compounds would be
defined as antagonistic (Figure 2B).
Figure 2 Examples of Loewe Additivity and Bliss Independence in defining drug interactions. A) Additivity, synergy and antagonism of
drug combination as defined by Loewe Additivity. Let x and y-axes represent concentrations of drugs A and B to achieve a defined effect X% (e.g., X =
50% for half maximal inhibitory concentration (IC50) of [IA]50% and [IB]50%), respectively. The coordinates ([IA]50%,0) and (0, [IB]50%) represent the concentration
for drugs A and B, respectively. The line of additivity is constructed by connecting these two points for a 50% effect isobologram plot. The concentrations
of the two drugs used in combination to provide the same effect X% (e.g. X = 50%) will be denoted by point ([CA]50%,[CB]50%) and are placed in the same
plot. Synergy, additivity, or antagonism will be determined when this point ([CA]50%,[CB]50%) is located below, on, or above the line, respectively. More
generally, linear, concave, and convex isoboles represent non-interacting, synergy, and antagonistic drug combination, respectively. B) Additivity, synergy
and antagonism of drug combination as defined by Bliss Independence. For example, if two non-interacting drugs (A and B) each result in 40% tumor
growth compared to control (EA = 0.4, EB = 0.4), then the predicted tumor growth when combined would be EC = (0.4 x 0.4) = 0.16, (16% of control). If the
observed combined (A + B, red bar) tumor growth is similar to, less than, or greater than 16% of control, then the combination would be deemed as
additive, synergistic, or antagonistic, respectively. N.D. denotes no drug (control).
These two methods produce different results, and it is
uncertain which method performs better with uncertainty
of mechanism and noisy data [17,20]. Drugs inhibiting
parts of the same linear pathway may act according to
Loewe additivity . Drugs nonexclusively affecting
parallel pathways may act according to Bliss independence.
Experimental characterization of drug combinations
typically involves generating dose response curves with the
inhibitors separately and combined. The experimental dose
response curve data can then be compared to the
predictions of Loewe additivity or Bliss independence to
determine if the drugs are acting synergistically.
Computational models of signaling networks
Given the complexity of the signaling networks
controlling systemic diseases such as cancer, computational
models of cell signaling pathways are important tools for
increasing understanding of pathological signaling and
prioritizing targets to test experimentally . Models
can be used to quantify systems properties that are often
not apparent in individual experiments. Through model
simulations, one can predict the relative importance of
various proteins in the network, the presence of signal
amplification, and the role of feedback and cross-talk
. These features will be important in the prediction
of viable drug combinations. While model predictions
require experimental validation, they are useful tools for
prioritizing targets for experimental planning.
Mass action and enzyme kinetics-based models
Three predominant signaling network modeling
approaches are mass action and enzyme kinetics-based,
logic-based, and statistical association-based models.
Mass action models are biochemically detailed kinetic
models that typically represent interactions between
molecular species in the signaling network as ordinary
differential equations (ODEs) and require selection of
parameter values for concentrations of species in the
network and rate constants controlling protein-protein
associations . Many of these parameters may be unavailable
in the literature and can either be measured experimentally
or fit to the data by minimizing an objective function such
as the sum of squared error.
When training parameters to data, it is important to
determine the importance of parameter selection on
model predictions. For example, Chen et al. measured
parameter sensitivity for several independent fits and
saw that the rank order of the most sensitive parameters
was nearly the same across the fits for a given output,
therefore parameter uncertainty did not affect major
model predictions . In another approach, Iadevaia
et al. developed a mass-action model of the IGF-1
signaling network in a breast cancer cell line with 161
unknown parameters and fit the model to the time courses
of six proteins measured with reverse-phase protein
array . Given the uncertainty in parameter
estimation with so many unknown parameters, they identified
ten sets of parameters using particle swarm optimization
that equally fit the experimental data. Model predictions
were averaged from three randomly sampled sets of the
ten parameter sets. The trained model was then used to
identify beneficial drug combinations in a breast cancer
Mass-action network models have been used to
predict new beneficial drug combinations for cancer. As an
example, Faratian et al. used a mass-action model of
heregulin-induced HER2/3 signaling through MAPK and
PI3K to study the role of PIK3CA activation in Receptor
Tyrosine Kinase (RTK) inhibitor resistance . Model
results demonstrated that the ratio of PTEN to activated
PIK3CA predicted resistance to RTK inhibitors. This
finding could therefore be used to predict patient
response to anti-HER2 therapies based on clinical
measurements of PTEN. It predicts that PIK3CA inhibition should
be paired with RTK inhibitors in patients with tumors
with low PTEN, a negative regulator of PI3K signaling.
Another group developed a mass action kinetics model of
PI3K signaling by ERBB receptors including ligand
binding, dimerization, internalization, recycling, and
degradation . Sensitivity analysis of this model predicted an
important role of ERBB3 in AKT activation, which was
then validated in mice xenografts. Sensitivity analysis
could be used in future work to find drug combinations
that may work synergistically with ERBB3 inhibition.
A limitation of mass-action modeling approaches is the
amount of data required to generate specific values for
the abundance and rate constant parameters, which can
be prohibitive for large scale network reconstructions
(Figure 1) . Logic-based models use network
topology without the need for specific parameter values.
Network interactions are modeled with OR, AND, and
NOT Boolean logic gates. Each species in the network
takes a value of 0 (inactive) or 1 (active) based on the
state of its effectors . As an example, Sahin et al.
developed a Boolean model of ERBB signaling of G1/S cell
cycle transition . The group used computational
knockouts of network proteins, validation experiments
with RNAi, and model revision based on proteomic data,
to predict the effects of combined inhibition of ERBB2
and c-MYC or EGFR. A combination therapy targeting
c-MYC and ERBB2 was predicted to improve treatment
for breast cancer that is de novo resistant to ERBB2
inhibition. Another group developed a Boolean logic
model of apoptosis signaling in Leukemic T-Cell large
granular lymphocytes . The authors used the model
to determine species that controlled apoptosis and
experimentally validated two of these species, sphingosine
kinase 1 and NFB. Given the limitations in representing
species as either on or off, this modeling approach has
been extended to accommodate intermediate activity
states using fuzzy logic .
Normalized hill differential equation modeling approach
While logic-based modeling approaches benefit from
simple construction using network topology, results can
be difficult to interpret due to assignment of discrete
values to continuous variable such as concentration of
active species, sensitivity to temporal node-updating
schemes, and incompatibility with many systems analysis
tools such as quantitative sensitivity analysis . To
address the limitations of mass-action and logic-based
models, Kraeutler et al. developed the normalized Hill
differential equation modeling approach, which uses
logic-based differential equations to represent activation
or inhibition by molecular species in the network .
Cross-talk is represented with AND and OR gates and
species activation is continuous over time and in units
of fractional activation instead of concentration.
Therefore protein abundance parameters are not required like
with mass-action models. Interactions between species
in the network are modeled with normalized Hill
equations with 3 parameters: reaction weight, half maximal
effective concentration (EC50), and Hill coefficient.
While these parameters can be fit to data, using default
values generated highly similar quantitative predictions
as a previously constructed detailed biochemical model
of the same pathway which used 88 parameters from
literature [33,34]. Therefore, this approach allows for
straightforward model construction of a known network
topology even if kinetic and abundance parameters are
unknown, like with logic-based modeling, while also
allowing for prediction of dynamics and systems analysis
tools such as quantitative sensitivity analysis.
The normalized-Hill modeling approach is a valuable
tool for model construction of larger networks with more
unknown parameters. As an example, Ryall et al. used this
approach to model the cardiac hypertrophy signaling
network, which contained 106 species and 193 reactions .
Since cardiac myocytes have minimal capacity for
proliferation, many of these pathways also regulate proliferation
in cancer cells . Quantitative systems analysis revealed
the most prevalent species involved in growth of cardiac
myocytes, prioritizing future experimental targets .
While Ras, the largest signaling hub, was the highest
influencer on cell size, the correlation between the
number of connections a species has and its influence was
low. Moreover, highly influential species were at many
levels in the network, not just close to the output level.
These findings demonstrate the need for model
reconstructions to predict important drug targets in cell
signaling networks. Highly influential species are not
obvious from intuition alone or data from gain or loss
of function studies of single genes .
Ryall et al.s analysis of the hypertrophy signaling
network also looked at the presence of different signaling
motifs such as bi-fan and feed-forward loops. Motifs can
affect network properties such as signal filtering,
acceleration, pulse generation, ultra-sensitivity, stability, and
robustness [38-40]. Yin et al. modeled three-node enzymatic
networks with many different topologies to study the
effect of topology on drug combinations . Model
simulations were conducted to identify motifs that could result
in synergy. Most of the combinations were not dependent
on parameter selection, demonstrating that network
topology can be used to predict synergistic combinations.
Moreover, synergistic drug combinations were found in
both parallel and series drug combinations. In a similar
study, Zhang et al. made reduced models of the
convergence of two signaling pathways on a target and observed
synergy in only a subset of the motifs . Synergy had a
greater likelihood in motifs with negative feedback
between the target and an upstream effector or mutual
inhibition between parallel signaling pathways. These
findings suggest that searching for synergistic motifs
within a cancer signaling network topology can be a
useful strategy in prioritizing drug combinations to test
experimentally. Networks exported into Cytoscape , a
open source software platform for network visualization,
can use the Netmatch plug-in  to quickly search for
motifs of interest.
Statistical association-based modeling approach
Network modeling approaches are useful when network
topology is known, but these approaches can be biased
towards established pathways and may miss novel
interactions. Statistical association-based models do not
depend on prior knowledge of pathways and instead use
correlations and patterns in experimental data to predict
network structure. As an example, Ryall et al. used data
of correlations among cell shape features and expression
of 12 genes relevant to cardiac hypertrophy to identify a
network map linking input modules to output modules
. Drug combinations could then be prioritized by
selecting targets that enabled adaptive module signaling
and prevented maladaptive module signaling. Molinelli
et al. developed a network inference algorithm based on
Belief Propagation  to construct networks from
phenotypic screen data . They applied their method
to screen data from a melanoma cell line and identified
both new and established pathway interactions and then
used the network to predict efficacious drug targets.
Another useful approach for network inference is Bayesian
network computational methods .
Signature-based approaches for predicting drug
Many effective drug combinations have been discovered
using correlations in gene expression signatures. One
useful tool for this is the Connectivity Map (CMap)
database . The first-generation CMap contained gene
expression profiles from three cancer cell lines perturbed
by 164 distinct small-molecule compounds. The second
generation of CMap (CMap 2.0) includes gene
expression profiles from 1,309 small molecules including FDA
approved drugs tested in five human cancer cell lines
[49,50]. This method assumes gene expression changes
can be used as a universal language to connect distinct
biological states (e.g. diseases), allowing for the
successful repurposing of compounds [51,52]. In short, drugs
known to be effective in one disease can serve as
candidates for use in other diseases marked by similar gene
expression changes. Users query the database to
compute similarity metrics between a test gene expression
signature and each reference set. Similarity metrics are
scaled between 1 and +1, where a positive score indicates
positive correlation and a negative score indicates negative
correlation. An advantage of this approach is the ability to
query CMap with publicly available gene expression data
from sources such as Gene Expression Omnibus (GEO)
, therefore facilitating rapid drug combinations
prediction for experimental validation .
As an example, Riedel et al. applied CMap to predict
drugs that would prevent resistance to chemotherapy
agents in lung cancer cell lines . Genes with the
highest changes after treatment with docetaxel were
analyzed using CMap to identify drugs with negative
connectivity scores, indicating these drugs had antagonistic
effects on the genes associated with docetaxel resistance.
PI3K inhibitor LY294002, which was highly ranked among
these antagonistic compounds, was tested in vitro with
docetaxel and found to synergistically increase cytotoxicity.
Wei et al. used a similar approach to predict drugs to
overcome resistance to glucocorticoid treatment in acute
lymphoblastic leukemia . Microarrays of pre-treated
cell lines either sensitive or resistant to glucocorticoid
in vitro were used to define a sensitive/resistant gene
signature. Using CMap, mTOR inhibitor rapamycin was
found to induce a highly similar signature, leading to
the hypothesis that mTOR inhibition could induce
glucocorticoid sensitivity. Follow-up experiments
supported this hypothesis and showed that rapamycin
conferred sensitivity through down-regulation of MCL1.
Therefore, CMap is a useful tool for using gene
signatures to predict drug combinations that may overcome
Inspired by the CMap concept, Kim et al. recently
developed K-Map (Kinase Inhibitor Connectivity Map) that
systematically connects a set of query kinases to kinase
inhibitors based on quantitative profiles of the kinase
inhibitor activities . Instead of gene expression
signatures, Kim et al. used the kinase activity profiles as the
language for connecting kinases and small molecules
in K-Map to reveal the complex interactions of kinases
and inhibitors. By querying K-Map with the essential
kinases mediating resistance to EGFR-inhibitor gefitinib in
an EGFR mutant non-small cell lung cancer (NSCLC)
cell line, bosutinib was predicted to be a more effective
drug for killing EGFR mutant cancer cells. Follow up
in vitro experiments confirmed that bosutinib alone is a
more effective agent than gefitinib, and that the
combination of bosutinib and gefitinib had synergistic effects in
EGFR mutant NSCLC cells . This demonstrates the
utility of K-Map in connecting kinases with kinase
inhibitors and suggesting candidates for drug combinations.
Network-based approaches for predicting drug
Other computational approaches have been developed
to predict drug combinations using data from
highthroughput screens and drug databases. Pal and Berlow
developed an algorithm based on set theory that uses
tumor drug sensitivities and kinase inhibition profiles for
a set of individual drugs to predict the tumor sensitivity
to new drugs or drug combinations . The algorithm
applies the following rules to generate circuit
representations of tumor pathways: 1) drugs that inhibit a superset
of an effective set of inhibited kinases will also be
successful in inhibiting tumor growth and 2) drugs
inhibiting subsets of ineffective sets of inhibited kinases will
also be unsuccessful. These circuits reveal a set of
kinases that are most predictive of drug sensitivity and
depict combinations of kinases that need to be inhibited to
prevent tumor growth. This analysis is helpful for
identifying drug combinations that inhibit a minimal set of
kinases with as few off target effects as possible to minimize
negative side effects. This approach was validated using
data from four canine cancer cell lines given 60 different
drugs at four different concentrations to generate IC50
values . Tang et al. expanded on this algorithm to
improve the computational cost and accuracy with drug
screen data with little overlap between drug target profiles
. While this approach requires a lot of experimental
data from drug screens to be useful (Figure 1),
technological advancements are enabling larger drug screens at
Given the high cost of exhaustive drug screens, Gujral
et al. exploited the polypharmacology of kinase inhibitors
by developing an approach to select the most predictive
kinase targets from a smaller scale drug screen of
multitarget drugs . They performed a phenotypic screen to
identify kinases regulating cell migration using an optimal
set of 32 kinase inhibitors. Elastic net regularization was
then used to deconvolute the polypharmacology of the
kinase inhibitors, identifying kinases with the greatest
explanatory power for the phenotype. Elastic net
regularization regresses an output variable against a set
of predictor variables (kinase activity) and invokes a
penalty on the number of variables in order to eliminate
kinases with insignificant contributions. This is a useful
approach for reducing experimental time and cost by
extracting more information from a smaller set of
compounds in drug screens.
Huang et al. used drug genomic profiles from the
Connectivity Map database to construct a drug functional
network and then grouped drugs into modules with
similar transcriptional responses . They then built
disease-signaling networks highlighting defective
signaling modules based on patient genomic profiles and
protein interaction data. The DrugComboRanker algorithm
ranked potential drug combinations by selecting drugs
with high overlap in the disease network, affecting
multiple key signaling modules. In a similar approach, Pang
et al. developed an algorithm to identify combination
therapies by building a network of drug-target interactions
from the DrugBank database . Given an input disease
gene set, the algorithm selected drugs that maximized on
target effects and minimized off target effects. This
algorithm also identifies drug combinations of more than two
drugs, which would be unfeasible to predict using a
highthroughput screen. These approaches allow researchers to
take advantage of publicly available drug data to prioritize
combinations for experimental validation.
Moreover, Liu et al. have developed a database of both
successful and unsuccessful drug combinations (DCDB)
. The current version of the database contains 1,363
drug combinations involving 904 individual drugs and
805 targets. The database provides information about
the potential mechanism, drug interactions, indication,
published study, development status, and targets. The
ability to analyze patterns in successful and unsuccessful
combinations with this database will be useful for systems
analysis of drug combinations and rational experimental
design. As an example, Xu et al. constructed a network of
successful drug interactions using data from the Drug
Combination Database . Analysis revealed that drugs
with similar therapeutic effects tended to cluster together
in the network and targets of hub drugs were often
membrane or membrane-associated proteins. They used these
observations to develop a statistical approach to predict
new drug combinations.
Cheng et al. constructed a global human kinome
interaction map by integrating kinase-substrate interactions,
kinase-drug interactions, protein-protein interactions, and
atomic resolution three-dimensional structural
proteinprotein interactions . Their analysis of the topological
features of these networks revealed an enrichment of hubs
as drug targets. While targeting hubs can be beneficial due
to cascading effects, their analysis revealed that targeting
hubs also increases risk of adverse drug reactions and drug
resistance through feedback and crosstalk.
Network based approaches have also been used to
study drug-drug interactions. As an example, Cheng and
Zhao developed a comprehensive drug-drug interaction
network incorporating 6946 interactions of 721
approved drugs using data from DrugBank . They then
calculated drug-drug pair similarities using four features:
phenotypic similarity, therapeutic similarity, chemical
structure similarity, and genomic similarity. They applied
five machine learning-based models to the dataset to
predict drug-drug interaction, with the overall hypothesis that
drugs with similar chemical structure, target proteins,
adverse drug reactions, and therapeutic purposes have high
probability of drug-drug interaction. They tested the
model on antipsychotic drug-drug interactions and found
literature support for predictions of drug-drug interactions
involving weight gain and P450 inhibition. This approach
demonstrates the power of harnessing network-based
drug-drug interactions to reveal new information on
adverse drug effects and provide additional filtering rules for
drug combination studies.
Integrating functional genomics and computational
methods for identifying drug combinations
Large scale knockdown screens using RNA interference
(RNAi) can also be used to identify potential drug
combinations . RNAi screens can identify genes that lead
to sensitivity or resistance to a drug of interest. As an
example, an RNAi screen conducted by Berns et al.
showed that knockdown of PTEN decreased sensitivity
to trastuzumab in BT-474 breast cancer cells .
Follow-up studies showed that assessment of both the
loss of PTEN expression and activating mutations in
PIK3CA could predict the risk for HER2 amplified
tumor progression. Drugs reducing PI3K signaling may
therefore increase response to trastuzumab. In another
study, Prahallad et al. used an RNAi screen to identify
kinases whose knockdown synergized with BRAF
(V600E) inhibition in colon tumors . Follow-up
experiments demonstrated synergy between cetuximab
(EGFR inhibitor) and vemurafenib (BRAF inhibitor).
The rational combination of cetuximab and
vemurafenib is currently being evaluated in clinical trials.
Pritchard et al. used RNAi signatures of eight cell
death genes to determine the mechanism of drug
combination effects in lymphoma cells . Single drugs
were classified based on their similarity to the RNAi
signatures of well-characterized compounds with known
mechanisms. They then generated signatures for drug
combinations to see if the signature was more similar to
results from one of the drugs alone, an average of the
two, or a unique signature. Results showed that the
combination signature was usually a weighted
composite of single drug effects where one drug potentiated the
mechanism of the other or the two drugs produced an
additive effect. Interestingly, they observed that
applying larger pools of drugs to tumors reduced the genetic
heterogeneity, which could be prohibitive in selection of
personalized drug treatments for patients based on
Spreafico et al. identified non-canonical Wnt pathway
mediated resistance to MEK1/2 inhibitor Selumetinib in
colorectal cancer cells by integrating gene set
enrichment analysis and synthetic lethality screens . Using
cyclosporine A (CsA) as a drug to inhibit non-canonical
Wnt pathway, they validated that the combination of
Selumetinib and CsA has synergistic anti-proliferative effects
both in vitro and in vivo patient-derived xenografts. This
rational combination is now being translated into a Phase
I clinical study (ClinicalTrials.gov ID: NCT02188264).
This illustrates the utility of integrating functional
genomics screens with bioinformatics to identify and translate
drug combinations into clinical study.
High-throughput drug combination screens
While the number of FDA approved drugs makes
exhaustive drug screens unfeasible (Figure 1), there have
been efforts to reduce the search space in an unbiased
manner. Tan et al. used pools of ten drugs in 384-well
plates to study all possible pairs of 1,000 compounds in
the minimum number of wells possible in order to find
drugs combinations that synergistically prevent HIV
replication . Synergistic wells from the primary screen
are then deconvolved into possible synergistic pairs for a
secondary screen. Results revealed an enrichment of
anti-inflammatory drugs in the synergistic pairs.
Roller et al. conducted a functional chemical genetic
screen of 300 drug combinations in nine melanoma cell
lines and identified pairs of compounds that synergistic
increase cytotoxicity . Interestingly, the synergistic
cytotoxicities identified did not correlate with the known
oncogene RAS and RAF mutational status of the
melanoma cell lines. From this screen, they identified sorafenib
(a multi-kinase inhibitor) and diclofenac (a non-steroidal
anti-inflammatory drug) to be the most robust drug
combination that has synergistic effects across the
melanoma cell lines. By using this functional chemical
genetic screen, the authors uncovered novel interactions
between signaling inhibitors that would not be predicted
based on current understanding of the signaling networks.
Their results suggest that the underlying signaling
networks controlling drug responses can vary substantially
based on unidentified elements of cell genotype. In
another study, Griner et al. conducted a large scale screen of
multiple concentrations of 500 compounds with ibrutinib
in activated B-cell-like diffuse large lymphoma cells .
They discovered many compounds that interacted
favorably with ibrutinib, including inhibitors of PI3K signaling,
the Bcl2-family, and the B-cell receptor pathway.
One of the major ongoing initiatives at the
Developmental Therapeutics Program of the National Cancer
Institute, U.S. National Institute of Health, is the
largescale high-throughput drug combination screening of
100 FDA approved drugs. The first set of screening
results generated 5,000 possible drug combinations in the
NCI-60 cancer cell lines panel . The goal of this
project is to identify novel drug combinations that are more
active than the single agents alone. As all the drugs
tested have been FDA approved, any drug combinations
identified from this screen may rapidly translate into the
clinic. As the NCI-60 cell lines panel have been fully
characterized by various omics technologies, the
release of this drug combination matrix to the public
could facilitate the development of novel computational
methods to integrate, predict, and mine the interactions
between molecular markers and drug combinations.
While substantial intratumoral heterogeneity has been
detected in cancer patients using next generation
sequencing technologies , current drug combination
prediction methods have primarily focused on targeting
the predominant tumor subpopulation. To study the
effect of different tumor subpopulations on treatment
efficacy, Zhao et al. developed a multi-objective linear
optimization algorithm to select optimal drug
combinations for heterogeneous tumors by maximizing efficacy
and minimizing toxicity . Their goal was to determine
the best drug combinations to minimize all
subpopulations. They experimentally validated the algorithms
prediction of two-drug combinations with three-component
heterogeneous tumors created using RNA interference
. They then expanded the model to simulate more
complex tumors and greater numbers of drugs . Their
results revealed that intratumor heterogeneity influences
the prediction of effective drug combinations. Different
predictions are made depending on if all tumor
subpopulations are considered or just the predominant
subpopulation. This approach represents a step forward of predicting
drug combinations to tackle tumor heterogeneity in the
era of precision oncology.
Given the experimental costs of exhaustively testing
drug combinations, computational models of signaling
networks will be especially useful in pre-clinical
screening of combinations of compounds. Model simulations
reveal non-intuitive effects of drug combinations .
Due to the size and complexity of cancer signaling,
modeling strategies accommodating reconstruction of larger
networks while still being compatible with quantitative
Figure 3 Future strategy for drug combination predictions with parallel integration of computational modeling, preclinical testing, and
clinical trials. A) Future combinatorial drug discovery approaches will benefit from tighter integration of gene signatures and phenotypic screen
data with computational models, tuning the models to specific cancer cell-lines. Model simulations enable prediction of effective drug combinations
for preclinical validation. Preclinical data can then be used to further refine computational models. B) For clinical application, patient gene signatures
can be clustered with gene expression signatures from previously modeled cell lines. Similarity scores can then be computed to find the most similar
model to the patients tumor for selection of the appropriate drug combination.
systems analysis tools will be especially useful . Global
network models allow for more comprehensive and
unbiased discovery of therapeutic targets than experimental
approaches based on a priori selection of important
pathways . One of these methods is the recently described
normalized-Hill modeling approach, which even with
T cell large granular
Table 2 Summaries of reviewed systems approaches for identifying drug combinations
Computational models of cell signaling networks
Mixed integer linear
Systems analysis of
Literature support 
Literature support 
Literature support 
Combined inhibition of MEK and PI3K optimally
decreased cell viability.
the ratio of PTEN to activated PI3K predicts RTK
ErbB3 inhibition inhibits the ErbB-PI3K network
more potently than current therapies.
Combined inhibition of c-MYC and ERBB2 improved
treatment for trastuzumab resistant breast cancer.
Sphingosine kinase 1 and NFKB are essential for survival
of leukemic T cell large granular lymphocytes.
MK2 and MEK are co-regulators of ERK and
EGF induced IKK inhibition.
Ras had the greatest influence on hypertrophy and
correlation between node degree and influence is low.
Identified consistent synergistic and antagonistic motifs.
Synergy is more prevalent in motifs with negative feedback
between the target and an upstream effector or mutual
inhibition between parallel pathways.
Maladaptive and adaptive hypertrophy features were in
separate modules in the simplified hypertrophy network
map generated by k-means clustering of ligands and
Reconstructed classic T cell signaling network using
multiparameter single-cell data and Bayesian network inference.
PI3K inhibition enhanced docetaxel-induced cytotoxicity
mTor inhibition induced glucocorticoid sensitivity
by decreasing MCL1
The combination of bosutinib and gefitinib has synergistic
effects in EGFR mutant non-small cell lung cancer
Developed an algorithm using a training set of drug sensitivities
with known targets to predict responses to new drugs
Target Inhibition inference using Maximization and Minimization
Averaging (TIMMA). Improved computational cost and accuracy
of the above TIM approach. Predicted kinase pairs that could
be inhibited to prevent cancer survival.
Performed phenotypic screen using an optimal set of 32
kinase inhibitors. They used an elastic net regulatization
algorithm to deconvolute the polypharmacology and
identify key kinases regulating cell migration.
Created drug and disease functional networks based on
genomic profiles and interactome data. Drug combinations
are predicted by identifying drugs whose targets are
enriched in the disease network.
Built a network of drug-target interactions from DrugBank.
Given an input gene set, the algorithm selects drug combinations
that maximize on target effects and minimize off target effects
Drugs with similar therapeutic effects cluster together in a
network of successful drug combinations produced using
the Drug Combination Database . Network observations
Table 2 Summaries of reviewed systems approaches for identifying drug combinations (Continued)
Drug-drug interactions Drug-drug interaction
Integration of functional genomics
and computational methods
8-gene RNAi signature
High-throughput drug combination screens
were used to develop a statistical approach
for predicting drug combinations (DCPred)
Applied five machine learning models to a data set
of drug-drug pair similarities including 721 approved
drugs to predict drug-drug interactions.
PTEN downregulation with active PI3K signaling
induce trastuzumab resistance
EGFR inhibition synergizes with BRAF(V600E) inhibition
Drug combination signatures were usually a weighted
composite of single drug effects
The combination of Selumetinib (MEK1/2 inhibitor) and
CsA (Wnt inhibitor) has synergistic anti-proliferative effects
Used pools of 10 drugs in 384-well plates to study all
possibly pairs of 1000 compounds in the minimum
number of wells possible
Sorafenib (a multi-kinase inhibitor) and diclofenac
(NSAID) had synergistic effects across all nine
tested melanoma cell lines.
Screen of 500 compounds with ibrutinib revealed
favorable combinations with inhibitors of PI3K signaling,
the Bcl2 family, and B-cell receptor pathway
Screen of 5,000 combinations of FDA-approved
drugs in the NCI-60 cancer cell line panel.
Intatumor heterogeneity influences the prediction
of effective drug combinations.
Literature support 
minimal biochemical data from literature enables a global
view of quantitative functional relationships between every
node in a network . This method was used to identify
the most important nodes in a integrative network model
of cardiac hypertrophy with 106 nodes and 193 reactions
and default parameter values . Many of these cardiac
hypertrophy signaling pathways also play important roles
in tumor growth. While model predictions need to be
experimentally validated, models can substantially improve
hypothesis generation and experimental planning.
In addition to larger network reconstructions, future
modeling efforts will benefit from tighter integration of
high-throughput sequencing, proteomics, and phenotypic
screen data (Figure 3A). This will enable tuning of a model
to an individual patients tumor, which would be beneficial
for use in personalized medicine. Comprehensive double
knockdown model simulations would enable personalized
prediction of drug combinations for patients. The ability
to readily adapt models to different situations is important
in cancer research since the molecular networks are not
fixed within a particular cancer type . Patients that
share the same mutation and tumor type can have
different responses to a drug . Genetic background, cell
lineage, and exogenous signals can influence the network
behavior . Efficacy data identified from in vitro and
in vivo experiments would then be used for model
refinement so more informed predictions of drug combinations
can be made in future studies.
Predicted drug combinations should be validated in
cancer cell lines and in relevant in vivo human disease
models such as patient-derived tumor xenografts .
These models, however, typically overestimate the
clinical benefit due to factors such as tumor heterogeneity,
differences in tumor microenvironment, and inaccurate
estimates of drug exposure . Therefore, it is
important to have a high threshold when choosing effective
combinations, ignoring modest inhibitions of tumor
growth in favor of combinations promoting cancer cell
death and tumor regression.
Another opportunity for improved design of
combination therapies is through quantitative systems
pharmacology approaches integrating cell signaling network models
with pharmacokinetic-pharmacodynamic (PK/PD) models
. Quantitative systems pharmacology uses multi-scale
data to better understand and ultimately predict how
drugs affect cellular networks and human pathophysiology
. Mechanistic models of cell signaling networks are
linked to PK/PD models of physiological processes at the
level of tissues and organisms. These models will enable
patient-specific prediction of therapeutic and toxic drug
responses and drug resistance mechanisms, improve
translation of in vitro discoveries to patients, and enhance
discovery of pharmacodynamic biomarkers.
Additionally, future work will benefit from parallel
integration of computational modeling, preclinical testing,
and clinical trials, where data from each approach can be
used for refinement of the other. Computational models
tuned to specific cancer cell lines using bioinformatics and
experimental data could be perturbed to make predictions
of effective drug combinations to validate in preclinical
models (Figure 3A). For clinical application, similarity
scores between patients and previously modeled cell lines
could be calculated using statistical clustering (Figure 3B).
Drug combinations predicted using the most similar
model could then be applied in the clinic.
It is becoming increasingly apparent that drug
combinations will be essential for improving therapies for complex
diseases such as cancer . The signaling pathways
controlling these systemic diseases are highly
interconnected, with cross-talk, redundancy, and feedback,
making single-target therapies much less effective .
While combination therapies have the potential to
prevent the development of resistance seen in many single
drug therapies, it is prohibitively expensive to
experimentally test every potential combination, especially when
considering combinations of more than two drugs. Here,
we highlighted a variety of systems biology applications
for advancing the prediction of effective drug
combinations, as summarized in Table 2. These methods include
computational modeling, gene signature analysis,
functional genomics, and high-throughput drug combination
screening. Utilization and integration of these systems
biology approaches hold great promise in speeding up the
development of clinically relevant drug combinations.
Karen Ryall is a post-doctoral research fellow at the University of Colorado School
of Medicine. Her research interests include systems biology, computational
modeling, and bioinformatics.
Aik Choon Tan is an Associate Professor of Bioinformatics at the University of
Colorado School of Medicine. His research interests include translational
bioinformatics and cancer systems biology.
This work is partly supported by the National Institutes of Health under Ruth L.
Kirschstein National Research Service Award T32CA17468 (K.A.R.), the National
Institutes of Health P50CA058187, Cancer League of Colorado, the Department
of Defense Award W81XWH-11-1-0527 and the David F. and Margaret T. Grohne
Family Foundation. Its contents are solely the responsibility of the authors and
do not necessarily represent the official views of the funders.
1. Pammolli F , Magazzini L , Riccaboni M. The productivity crisis in pharmaceutical R& D. Nat Rev Drug Discov . 2011 ; 10 ( 6 ): 428 - 38 .
2. Hopkins AL . Network pharmacology: the next paradigm in drug discovery . Nat Chem Biol . 2008 ; 4 ( 11 ): 682 - 90 .
3. Sams-Dodd F. Target-based drug discovery: is something wrong? Drug Discov Today . 2005 ; 10 ( 2 ): 139 - 47 .
4. Heineke J , Molkentin JD . Regulation of cardiac hypertrophy by intracellular signalling pathways . Nat Rev Mol Cell Biol . 2006 ; 7 ( 8 ): 589 - 600 .
5. Breitkreutz D , Hlatky L , Rietman E , Tuszynski JA. Molecular signaling network complexity is correlated with cancer patient survivability . Proc Natl Acad Sci U S A . 2012 ; 109 ( 23 ): 9209 - 12 .
6. Kitano H. Cancer as a robust system: implications for anticancer therapy . Nat Rev Cancer . 2004 ; 4 ( 3 ): 227 - 35 .
7. Albert R , Jeong H , Barabasi AL. Error and attack tolerance of complex networks . Nature . 2000 ; 406 ( 6794 ): 378 - 82 .
8. Lusis AJ , Weiss JN . Cardiovascular networks: systems-based approaches to cardiovascular disease . Circulation . 2010 ; 121 ( 1 ): 157 - 70 .
9. Berger SI , Iyengar R. Network analyses in systems pharmacology . Bioinformatics . 2009 ; 25 ( 19 ): 2466 - 72 .
10. Glickman MS , Sawyers CL . Converting cancer therapies into cures: lessons from infectious diseases . Cell . 2012 ; 148 ( 6 ): 1089 - 98 .
11. Gorre ME , Mohammed M , Ellwood K , Hsu N , Paquette R , Rao N , et al. Clinical resistance to STI-571 cancer therapy caused by BCR-ABL gene mutation or amplification . Science . 2001 ; 293 ( 5531 ): 876 - 80 .
12. Engelman JA , Zejnullahu K , Mitsudomi T , Song Y , Hyland C , Park JO , et al. MET amplification leads to gefitinib resistance in lung cancer by activating ERBB3 signaling . Science . 2007 ; 316 ( 5827 ): 1039 - 43 .
13. Sharma SV , Lee DY , Li B , Quinlan MP , Takahashi F , Maheswaran S , et al. A chromatin-mediated reversible drug-tolerant state in cancer cell subpopulations . Cell . 2010 ; 141 ( 1 ): 69 - 80 .
14. Weinstein IB , Joe A. Oncogene addiction . Cancer Res . 2008 ; 68 ( 9 ): 3077 - 80 .
15. Sun X , Vilar S , Tatonetti NP. High-throughput methods for combinatorial drug discovery . Sci Transl Med . 2013 ; 5 ( 205 ): 205rv1 .
16. Tang J , Aittokallio T. Network pharmacology strategies toward multi-target anticancer therapies: from computational models to experimental design principles . Curr Pharm Des . 2014 ; 20 ( 1 ): 23 - 36 .
17. Fitzgerald JB , Schoeberl B , Nielsen UB , Sorger PK . Systems biology and combination therapy in the quest for clinical efficacy . Nat Chem Biol . 2006 ; 2 ( 9 ): 458 - 66 .
18. Chou T , Talalay P. Analysis of combined drug effects: a new look at a very old problem . Trends Pharmacol Sci . 1983 ; 4 : 450 - 4 .
19. Keith CT , Borisy AA , Stockwell BR . Multicomponent therapeutics for networked systems . Nat Rev Drug Discov . 2005 ; 4 ( 1 ): 71 - 8 .
20. Goldoni M , Johansson C. A mathematical approach to study combined effects of toxicants in vitro: evaluation of the Bliss independence criterion and the Loewe additivity model . Toxicol In Vitro. 2007 ; 21 ( 5 ): 759 - 69 .
21. Kreeger PK , Lauffenburger DA . Cancer systems biology: a network modeling perspective . Carcinogenesis . 2010 ; 31 ( 1 ): 2 - 8 .
22. Klipp E , Liebermeister W. Mathematical modeling of intracellular signaling pathways . BMC Neurosci . 2006 ;7 Suppl 1: S10 .
23. Schoeberl B , Eichler-Jonsson C , Gilles ED , Mller G. Computational modeling of the dynamics of the MAP kinase cascade activated by surface and internalized EGF receptors . Nat Biotechnol . 2002 ; 20 ( 4 ): 370 - 5 .
24. Chen WW , Schoeberl B , Jasper PJ , Niepel M , Nielsen UB , Lauffenburger DA , et al. Input-output behavior of ErbB signaling pathways as revealed by a mass action model trained against dynamic data . Mol Syst Biol . 2009 ; 5 ( 239 ): 239 .
25. Iadevaia S , Lu Y , Morales FC , Mills GB , Ram PT . Identification of optimal drug combinations targeting cellular networks: integrating phospho-proteomics and computational network analysis . Cancer Res . 2010 ; 70 ( 17 ): 6704 - 14 .
26. Faratian D , Goltsov A , Lebedeva G , Sorokin A , Moodie S , Mullen P , et al. Systems biology reveals new strategies for personalizing cancer medicine and confirms the role of PTEN in resistance to trastuzumab . Cancer Res . 2009 ; 69 ( 16 ): 6713 - 20 .
27. Schoeberl B , Pace EA , Fitzgerald JB , Harms BD , Xu L , Nie L , et al. Therapeutically targeting ErbB3: a key node in ligand-induced activation of the ErbB receptor-PI3K axis . Sci Signal . 2009 ; 2 ( 77 ): ra31 .
28. Papin JA , Palsson BO . Topological analysis of mass-balanced signaling networks: a framework to obtain network properties including crosstalk . J Theor Biol . 2004 ; 227 ( 2 ): 283 - 97 .
29. Morris MK , Saez-Rodriguez J , Sorger PK , Lauffenburger DA. Logic-based models for the analysis of cell signaling networks . Biochemistry . 2010 ; 49 ( 15 ): 3216 - 24 .
30. Sahin O , Frhlich H , Lbke C , Korf U , Burmester S , Majety M , et al. Modeling ERBB receptor-regulated G1/S transition to find novel targets for de novo trastuzumab resistance . BMC Syst Biol . 2009 ; 3 : 1 .
31. Zhang R , Shah MV , Yang J , Shah MV , Yang J , Nyland SB , et al. Network model of survival signaling in large granular lymphocyte leukemia . Proc Natl Acad Sci U S A . 2008 ; 105 ( 42 ): 16308 - 13 .
32. Aldridge BB , Saez-Rodriguez J , Muhlich JL , Sorger PK , Lauffenburger DA . Fuzzy logic analysis of kinase pathway crosstalk in TNF/EGF/insulin-induced signaling . PLoS Comput Biol . 2009 ; 5 ( 4 ): e1000340 .
33. Kraeutler MJ , Soltis AR , Saucerman JJ . Modeling cardiac beta-adrenergic signaling with normalized-Hill differential equations: comparison with a biochemical model . BMC Syst Biol . 2010 ; 4 ( 1 ): 157 .
34. Saucerman JJ , Brunton LL , Michailova AP , McCulloch AD . Modeling beta-adrenergic control of cardiac myocyte contractility in silico . J Biol Chem . 2003 ; 278 ( 48 ): 47997 - 8003 .
35. Ryall KA , Holland DO , Delaney KA , Kraeutler MJ , Parker AJ , Saucerman JJ . Network reconstruction and systems analysis of cardiac myocyte hypertrophy signaling . J Biol Chem . 2012 ; 287 ( 50 ): 42259 - 68 .
36. Molkentin JD , Dorn GW . Cytoplasmic signaling pathways that regulate cardiac hypertrophy . Annu Rev Physiol . 2001 ; 63 : 391 - 426 .
37. Kestler HA , Wawra C , Kracher B , Khl M. Network modeling of signal transduction: establishing the global view . Bioessays . 2008 ; 30 ( 11 - 12 ): 1110 - 25 .
38. Alon U. Biological networks: the tinkerer as an engineer . Science . 2003 ; 301 ( 5641 ): 1866 - 7 .
39. Alon U. Network motifs: theory and experimental approaches . Nat Rev Genet . 2007 ; 8 ( 6 ): 450 - 61 .
40. Pelaez N , Carthew RW . Biological robustness and the role of microRNAs: a network perspective . Curr Top Dev Biol . 2012 ; 99 : 237 - 55 .
41. Yin N , Ma W , Pei J , Ouyang Q , Tang C , Lai L. Synergistic and antagonistic drug combinations depend on network topology . PLoS One . 2014 ; 9 ( 4 ): e93960 .
42. Zhang Y , Smolen P , Baxter DA , Byrne JH . Computational analyses of synergism in small molecular network motifs . PLoS Comput Biol . 2014 ; 10 ( 3 ): e1003524 .
43. Shannon P , Markiel A , Ozier O , Baliga NS , Wang JT , Ramage D , et al. Cytoscape: a software environment for integrated models of biomolecular interaction networks . Genome Res . 2003 ; 13 ( 11 ): 2498 - 504 .
44. Ferro A , Giugno R , Pigola G , Pulvirenti A , Skripin D , Bader GD , et al. NetMatch: a Cytoscape plugin for searching biological networks . Bioinformatics . 2007 ; 23 ( 7 ): 910 - 2 .
45. Ryall KA , Bezzerides VJ , Rosenzweig A , Saucerman JJ . Phenotypic screen quantifying differential regulation of cardiac myocyte hypertrophy identifies CITED4 regulation of myocyte elongation . J Mol Cell Cardiol . 2014 ; 72 : 74 - 84 .
46. Yedidia JS , Freeman WT , Weiss Y. Understanding belief propagation and its generalizations . In: Exploring Artificial Intelligence in the New Millennium . 2002 . p. 239 - 69 .
47. Molinelli EJ , Korkut A , Wang W , Miller ML , Gauthier NP , Jing X , et al. Perturbation biology: inferring signaling networks in cellular systems . PLoS Comput Biol . 2013 ; 9 ( 12 ): e1003290 .
48. Sachs K , Perez O , Pe'er D , Lauffenburger DA , Nolan GP . Causal proteinsignaling networks derived from multiparameter single-cell data . Science . 2005 ; 308 ( 5721 ): 523 - 9 .
49. Lamb J , Crawford ED , Peck D , Modell JW , Blat IC , Wrobel MJ , et al. The connectivity map: using gene-expression signatures to connect small molecules , genes, and disease. Science . 2006 ; 313 ( 5795 ): 1929 - 35 .
50. Qu XA , Rajpal DK . Applications of connectivity map in drug discovery and development . Drug Discov Today . 2012 ; 17 ( 23 - 24 ): 1289 - 98 .
51. Hieronymus H , Lamb J , Ross KN , Peng XP , Clement C , Rodina A , et al. Gene expression signature-based chemical genomic prediction identifies a novel class of HSP90 pathway modulators . Cancer Cell . 2006 ; 10 ( 4 ): 321 - 30 .
52. Wei G , Twomey D , Lamb J , Schlis K , Agarwal J , Stam RW , et al. Gene expression-based chemical genomics identifies rapamycin as a modulator of MCL1 and glucocorticoid resistance . Cancer Cell . 2006 ; 10 ( 4 ): 331 - 42 .
53. Edgar R , Domrachev M , Lash AE . Gene Expression Omnibus: NCBI gene expression and hybridization array data repository . Nucleic Acids Res . 2002 ; 30 ( 1 ): 207 - 10 .
54. Sirota M , Dudley JT , Kim J , Chiang AP , Morgan A , Sweet-Cordero A , et al. Discovery and preclinical validation of drug indications using compendia of public gene expression data . Sci Trans Med . 2011 ; 3 ( 96 ): 96ra77 .
55. Riedel RF , Porrello A , Pontzer E , Chenette EJ , Hsu DS , Balakumaran B , et al. A genomic approach to identify molecular pathways associated with chemotherapy resistance . Mol Cancer Ther . 2012 ; 11 ( 5 ): 1214 - 5 .
56. Kim J , Yoo M , Kang J , Tan AC. K-Map: connecting kinases with therapeutics for drug repurposing and development . Hum Genomics . 2013 ; 7 ( 1 ): 20 .
57. Kim J , Vasu VT , Mishra R , Singleton KR , Yoo M , Leach SM , et al. Bioinformatics-driven discovery of rational combination for overcoming EGFR-mutant lung cancer resistance to EGFR therapy . Bioinformatics . 2014 ; 30 ( 17 ): 2393 - 8 .
58. Pal R , Berlow N. A kinase inhibition map approach for tumor sensitivity prediction and combination therapy design for targeted drugs . Pac Symp Biocomput . 2012 ; 351 - 62 .
59. Berlow N , Davis LE , Cantor EL , Sguin B , Keller C , Pal R. A new approach for prediction of tumor sensitivity to targeted drugs based on functional data . BMC Bioinformatics . 2013 ; 14 ( 1 ): 239 .
60. Tang J , Karhinen L , Xu T , Szwajda A , Yadav B , Wennerberg K , et al. Target inhibition networks: predicting selective combinations of druggable targets to block cancer survival pathways . PLoS Comput Biol . 2013 ; 9 ( 9 ): e1003226 .
61. Gujral TS , Peshkin L , Kirschner MW . Exploiting polypharmacology for drug target deconvolution . Proc Natl Acad Sci U S A . 2014 ; 111 ( 13 ): 5048 - 53 .
62. Huang L , Li F , Sheng J , Xia X , Ma J , Zhan M , et al. DrugComboRanker: drug combination discovery based on target network analysis . Bioinformatics . 2014 ; 30 ( 12 ): i228 - 36 .
63. Pang K , Wan Y , Choi WT , Donehower LA , Sun J , Pant D , et al. Combinatorial therapy discovery using mixed integer linear programming . Bioinformatics . 2014 ; 30 ( 10 ): 1456 - 63 .
64. Liu Y , Hu B , Fu C , Chen X. DCDB : drug combination database . Bioinformatics . 2010 ; 26 ( 4 ): 587 - 8 .
65. Xu KJ , Song J , Zhao XM . The drug cocktail network . BMC Syst Biol . 2012 ;6 Suppl 1: S5 .
66. Cheng F , Zhao Z. Machine learning-based prediction of drug-drug interactions by integrating drug phenotypic, therapeutic, chemical, and genomic properties . J Am Med Inf Assoc . 2014 ; 21 (e2): e278 - 86 .
67. Iorns E , Lord CJ , Turner N , Ashworth A. Utilizing RNA interference to enhance cancer drug discovery . Nat Rev Drug Discov . 2007 ; 6 ( 7 ): 556 - 68 .
68. Berns K , Horlings HM , Hennessy BT , Madiredjo M , Hijmans EM , Beelen K , et al. A functional genetic approach identifies the PI3K pathway as a major determinant of trastuzumab resistance in breast cancer . Cancer Cell . 2007 ; 12 ( 4 ): 395 - 402 .
69. Prahallad A , Sun C , Huang S , Nicolantonio FD , Salazar R , Zecchin D , et al. Unresponsiveness of colon cancer to BRAF(V600E) inhibition through feedback activation of EGFR . Nature. 2012 ; 483 ( 7387 ): 100 - 3 .
70. Pritchard JR , Bruno PM , Gilbert LA , Capron KL , Lauffenburger DA , Hemann MT . Defining principles of combination drug mechanisms of action . Proc Natl Acad Sci U S A . 2012 ; 110 ( 2 ): E170 - 9 .
71. Spreafico A , Tentler JJ , Pitts TM , Tan AC , Gregory MA , Arcaroli JJ , et al. Rational combination of a MEK inhibitor, selumetinib, and the Wnt/calcium pathway modulator, cyclosporin A, in preclinical models of colorectal cancer . Clin Cancer Res . 2013 ; 19 ( 15 ): 4149 - 62 .
72. Tan X , Hu L , Luquette LJ , Gao G , Liu Y , Qu H , et al. Systematic identification of synergistic drug pairs targeting HIV . Nat Biotechnol . 2012 ; 30 ( 11 ): 1125 - 30 .
73. Roller D , Axelrod M , Capaldo B , Jensen K , Mackey A , Michael J. Synthetic lethal screening with small molecule inhibitors provides a pathway to rational combination therapies for melanoma . Mol Cancer Ther . 2012 ; 11 ( 11 ): 2505 - 15 .
74. Matthews Griner LA , Guha R , Shinn P , Young RM , Keller JM , Liu D , et al. High-throughput combinatorial screening identifies drugs that cooperate with ibrutinib to kill activated B-cell-like diffuse large B-cell lymphoma cells . Proc Natl Acad Sci U S A . 2014 ; 111 ( 6 ): 2349 - 54 .
75. Holbeck S , Collins JM , Doroshow JH . NCI-60 combination screening matrix of approved anticancer drugs . Eur J Cancer . 2012 ; 48 (Suppl:6): 11 .
76. Marusyk A , Almendro V , Polyak K. Intra-tumour heterogeneity: a looking glass for cancer ? Nat Rev Cancer . 2012 ; 12 ( 5 ): 323 - 34 .
77. Zhao B , Hemann MT , Lauffenburger DA . Intratumor heterogeneity alters most effective drugs in designed combinations . Proc Natl Acad Sci U S A . 2014 ; 111 ( 29 ): 10773 - 8 .
78. Zhao B , Pritchard JR , Lauffenburger DA , Hemann MT . Addressing genetic tumor heterogeneity through computationally predictive combination therapy . Cancer Discov . 2014 ; 4 ( 2 ): 166 - 74 .
79. Pe'er D , Hacohen N. Principles and strategies for developing network models in cancer . Cell . 2011 ; 144 ( 6 ): 864 - 73 .
80. Sharma SV , Haber DA , Settleman J. Cell line-based platforms to evaluate the therapeutic efficacy of candidate anticancer agents . Nat Rev Cancer . 2010 ; 10 ( 4 ): 241 - 53 .
81. Tentler JJ , Tan AC , Weekes CD , Jimeno A , Leong S , Pitts TM , et al. Patientderived tumour xenografts as models for oncology drug development . Nat Rev Clin Oncol . 2012 ; 9 ( 6 ): 338 - 50 .
82. Lieu CH , Tan AC , Leong S , Diamond JR , Eckhardt SG . From bench to bedside: lessons learned in translating preclinical studies in cancer drug development . J Natl Cancer Inst . 2013 ; 105 ( 19 ): 1441 - 56 .
83. Iyengar R , Zhao S , Chung SW , Mager DE , Gallo JM . Merging systems biology with pharmacodynamics . Sci Transl Med . 2012 ; 4 ( 126 ): 126ps7 .
84. Sorger PK , Allerheiligen SRB . Quantitative and systems pharmacology in the post-genomic Era: New approaches to discovering drugs and understanding therapeutic mechanisms. An NIH white paper by the QSP workshop group . 2011 . p. 0 - 47 .
85. Feala JD , Cortes J , Duxbury PM , Piermarocchi C , McCulloch AD , Paternostro G . Systems approaches and algorithms for discovery of combinatorial therapies . Wiley Interdiscip Rev Syst Biol Med . 2010 ; 2 ( 2 ): 181 - 93 .