Close
The page header's logo
About
FAQ
Home
Collections
Login
USC Login
Register
0
Selected 
Invert selection
Deselect all
Deselect all
 Click here to refresh results
 Click here to refresh results
USC
/
Digital Library
/
University of Southern California Dissertations and Theses
/
Computational tools for drug repurposing
(USC Thesis Other) 

Computational tools for drug repurposing

doctype icon
play button
PDF
 Download
 Share
 Open document
 Flip pages
 More
 Download a page range
 Download transcript
Copy asset link
Request this asset
Transcript (if available)
Content

COMPUTATIONAL TOOLS FOR DRUG REPURPOSING

by

Belinda B. Garana


A Dissertation Presented to the
FACULTY OF THE USC GRADUATE SCHOOL
UNIVERSITY OF SOUTHERN CALIFORNIA
In Partial Fulfillment of the
Requirements for the Degree
DOCTOR OF PHILOSOPHY
(CHEMICAL ENGINEERING)



August 2023






Copyright 2023          Belinda B. Garana

ii
Table of Contents
List of Figures ..................................................................................................................... iv
Abbreviations ...................................................................................................................... vi
Abstract .............................................................................................................................. vii
Chapter 1: Introduction ....................................................................................................... 1
Chapter 2: Drug mechanism enrichment analysis improves prioritization of  
therapeutics for repurposing .............................................................................................. 4
Background ..................................................................................................................... 4
Methods........................................................................................................................... 9
Drug mechanism enrichment analysis (DMEA) .................................................................................... 9
Simulation study of DMEA .................................................................................................................. 10
CMap L1000 query .............................................................................................................................. 11
CMap PRISM query ............................................................................................................................ 12
Weighted gene voting (WGV) ............................................................................................................. 13
DMEA using WGV molecular classification scores ............................................................................. 15
Simulation study of DMEA using WGV molecular classification scores ............................................. 16
DrugEnrichr and Drugmonizome ........................................................................................................ 17
Cell culture .......................................................................................................................................... 18
Cell viability experiments ..................................................................................................................... 18
Implementation ............................................................................................................. 19
DMEA: web application ....................................................................................................................... 19
Drug rank lists ................................................................................................................................. 19
Gene signatures.............................................................................................................................. 21
DMEA: R package ............................................................................................................................... 22
Results .......................................................................................................................... 23
DMEA identifies an enriched drug MOA in simulated data ................................................................. 23
DMEA identifies similar MOAs based on gene expression connectivity scores ................................. 24
DMEA identifies selectively toxic MOAs based on cell viability connectivity scores .......................... 29
DMEA identifies selectively toxic MOAs based on molecular signatures ........................................... 32
DMEA identifies potential senescence-inducing and senolytic drug MOAs for primary human
mammary epithelial cells ..................................................................................................................... 37
Discussion ..................................................................................................................... 40
Chapter 3: Mechanism of action landscape reveals shared sensitivities in large-scale
drug screens ..................................................................................................................... 45
Background ................................................................................................................... 45
Methods......................................................................................................................... 47
Mechanism of action landscape (MOAL) ............................................................................................ 47
MOAL: simulations .............................................................................................................................. 49
MOAL: CEVIChE ................................................................................................................................. 50
MOAL: weighted gene voting (WGV) .................................................................................................. 50
UMAP analysis .................................................................................................................................... 51
Results .......................................................................................................................... 52
MOA landscape detects coordinated enrichment of drug sets ........................................................... 52
MOA landscape removes more than 90% of drug MOA pairs from further consideration ................. 53

iii
MOA enrichment is correlated across large-scale drug screens ........................................................ 55
MOA relationships can depend on tissue type ................................................................................... 60
MOA landscape can evaluate network-level performance of drug sensitivity prediction algorithms .. 63
Discussion ..................................................................................................................... 65
References........................................................................................................................ 70
Appendices ....................................................................................................................... 80
Appendix A: Drug mechanism enrichment analysis improves prioritization of
therapeutics for repurposing ......................................................................................... 80
Appendix B: Mechanism of action landscape reveals shared sensitivities in large-
scale drug screens ........................................................................................................ 86

 

iv
List of Figures
Figure 1. DMEA is more flexible and statistically rigorous than other approaches to
evaluate drug MOA............................................................................................................. 7
Figure 2. Overview of drug mechanism enrichment analysis ......................................... 10
Figure 3. Sensitivity analysis of DMEA using synthetic data .......................................... 23
Figure 4. DMEA identifies similar MOAs based on gene expression connectivity scores
.......................................................................................................................................... 26
Figure 5. DMEA identifies selectively toxic MOAs based on cell viability connectivity
scores................................................................................................................................ 31
Figure 6. DMEA identifies selectively toxic MOAs based on external gene expression
signatures of intrinsic EGFR inhibitor resistance and acquired RAF inhibitor  
resistance, respectively .................................................................................................... 35
Figure 7. DMEA identifies potential senescence-inducing and senolytic drug MOAs  
for primary HMECs ........................................................................................................... 39
Figure 8. The MOA landscape approach to identify shared drug sensitivities ............... 49
Figure 9. MOA landscape detects coordinated enrichment of drug sets ........................ 53
Figure 10. MOA landscape removes more than 90% of drug MOA pairs from further
consideration .................................................................................................................... 55
Figure 11. MOA enrichment is correlated across large-scale drug screens................... 57
Figure 12. MOA landscape identifies drug pairings not obvious from UMAP  
visualization ...................................................................................................................... 59
Figure 13. MOA relationships can depend on tissue type .............................................. 62
Figure 14. MOA landscape can evaluate network-level performance of drug  
sensitivity prediction algorithms ....................................................................................... 64
Supplemental Figure 1. DMEA generates MOA rankings that are improved over  
single-drug rankings and similar to or better than other tools’ rankings of drug MOA ... 80
Supplemental Figure 2. DMEA identifies similar MOAs based on gene expression
connectivity scores ........................................................................................................... 81
Supplemental Figure 3. Overview of drug mechanism enrichment analysis using  
WGV molecular classification scores ............................................................................... 82

v
Supplemental Figure 4. Sensitivity analysis of DMEA pipeline with WGV using  
synthetic data .................................................................................................................... 83
Supplemental Figure 5. DMEA identifies enrichment of the EGFR inhibitor MOA and
other MOAs using external signatures of intrinsic EGFR inhibitor resistance ................ 85
Supplemental Figure 6. MOA relationships can depend on tissue type ......................... 87

 

vi
Abbreviations
AUC: Area under the curve
CCLE: Cancer cell line encyclopedia
CMap: Connectivity map
CSV: Comma-separated values
DMEA: Drug mechanism enrichment analysis
DSEA: Drug-set enrichment analysis
ES: Enrichment score
FDR: False discovery rate
GCT: Gene cluster text
GSEA: Gene set enrichment analysis
HGNC: HUGO Gene Nomenclature Committee
HMEC: Human mammary epithelial cells
HMGCR: 3-hydroxy-3-methylglutaryl-CoA reductase
HUGO: Human Genome Organization
LC-MS: Liquid chromatography-mass spectrometry
MOA: Mechanism of action
MOAL: Mechanism of action landscape
NES: Normalized enrichment score
NSCLC: Non-small cell lung cancer
PCL: Perturbagen class
PRISM: Profiling relative inhibition simultaneously in mixtures
WGV: Weighted gene voting  

vii
Abstract
Most patients are ineligible for targeted therapies and left to choose from other
treatments which can have deadly side effects. Considering this, my work has focused
on repurposing existing targeted therapies for new disease applications (i.e., drug
repurposing).

Though many drug repurposing algorithms exist, they typically only evaluate
individual drugs which can be prone to off-target effects and output a long list of drugs,
making it difficult for researchers to prioritize drug targets. To remove the assumption that
each drug works by its known mechanism(s), I developed software to group drugs by
shared mechanism-of-action (MOA). The utility of this approach is that drug MOAs are
only identified as candidates if most of the drugs annotated with that MOA are ranked as
strong candidates. By aggregating results across many drugs to evaluate drug MOAs
instead of just individual drugs, researchers can be more confident in their treatment
approach.

Altogether, my computational tools improve prioritization of drugs for repurposing
and can even identify therapeutics which may mitigate resistance. My software does not
require the user to possess programming knowledge and is freely available online for the
public. While experimental conditions can affect drug toxicity, the computational tools I
have developed will help researchers prioritize drugs to evaluate further.


1
Chapter 1: Introduction
There is a pressing need to identify targeted therapeutics for diseases. Targeted
therapeutics are an attractive treatment approach because they are designed to target
vulnerabilities unique to diseased cells whereas other treatments such as chemotherapy
can have deadly side effects due to their effects on normal tissues. Unfortunately, most
patients are ineligible for targeted therapies (e.g., more than eighty five percent of
cancer patients) and only about half of eligible cancer patients are responsive to these
treatments (Haslam et al., 2021). Considering this, my objective has been to better
prioritize targeted therapies.

My approach has been to leverage published data on existing drugs to repurpose
them for new disease applications (i.e., drug repurposing) since developing a new drug
can take about ten years on average and cost more than one billion dollars. Though many
drug repurposing algorithms have been developed previously (Parca et al., 2019; Wei et
al., 2019; Nguyen et al., 2017; Gao et al., 2021; Napolitano et al., 2018; Fang et al., 2021;
Rivas-Barragan et al., 2020; Domingo-Fernández et al., 2022; Szalai et al., 2019;
Napolitano et al., 2015; Chan et al., 2019; Emon et al., 2020; Rampášek et al., 2019; He
et al., 2023; Jia et al., 2016; Han et al., 2021; Chen et al., 2022; Wu et al., 2022), most
only evaluate individual drugs which can be prone to off-target effects and output a long
list of drugs with limited information about their molecular targets, making it difficult for
researchers to prioritize drug targets. In fact, several drugs have been found to kill cells
via off-target mechanisms (Lin et al., 2019). To avoid relying on the assumption that each
drug works by its known mechanism(s), I developed a method called Drug Mechanism

2
Enrichment Analysis to group drugs by common mechanism-of-action (MOA). With my
software, drug MOAs are only identified as candidates if most of the drugs annotated with
that MOA are ranked as strong candidates. By evaluating drug MOAs instead of just
individual drugs, we can be more confident in our treatment approach.

While there are other algorithms which evaluate drug mechanisms-of-action
(MOAs) (Subramanian et al., 2017; Kuleshov et al., 2016, 2019; E. Y. Chen et al., 2013;
Kropiwnicki et al., 2021; Huang et al., 2018; Napolitano et al., 2015), they have several
limitations. Firstly, they only query select public databases (e.g., the Connectivity Map
L1000 transcriptional database). Additionally, they have limited statistical rigor (e.g., lack
p or q values or do not use permutation-based metrics). Thirdly, they only accept one type
of unranked input list (e.g., gene symbols or drug names). Lastly, they do not provide
MOA-specific plots. To address these needs, I developed my software to be compatible
with any dataset or other drug repurposing algorithm, use permutation-based statistics,
and generate plots of both the overall and MOA-specific results.

Though I have shown with many different datasets that my approach improves
prioritization of targeted therapies for repurposing, resistance is a major barrier to
successful treatment. For example, even if a treatment is initially successful, cancers
almost always become resistant to single-drug treatments (Szakács et al., 2006). As
such, it would be advantageous to identify combination treatments which can mitigate
resistance. However, there is limited data available for combination treatments because
they are exponentially more expensive to test compared to individual drugs.

3

To prioritize drug MOAs for further testing in combination, I developed a method
termed the MOA Landscape. With this approach, I evaluate drugs and their MOAs based
on similar toxicity in published drug screens. Though this framework does not evaluate
synergy, it can reduce the number of drug MOA pairs by more than 90% for further
evaluation. Furthermore, while it can detect known drug combinations, it also detects less
studied combination treatments which may not be obvious from other dimensionality
reduction approaches and tissue-dependent relationships. Additionally, this software can
be used to evaluate the network-level performance of drug sensitivity prediction
algorithms.

Altogether, my work improves prioritization of therapeutics for repurposing and can
even identify treatments which may target drug-resistant cells. My software and results
are freely available online to the public and designed such that programming knowledge
is not required of the user. While drug toxicity can vary in different experimental
conditions, my results show that the computational tools I have developed will help
researchers prioritize drugs for further testing.

 

4
Chapter 2: Drug mechanism enrichment analysis improves
prioritization of therapeutics for repurposing
Background
Identifying effective therapeutics for diseases remains a pressing challenge. Recent
efforts in large-scale ‘omic profiling (Haoxin Li et al., 2019; Nusinow et al., 2020; Ghandi
et al., 2019; Frejno et al., 2020), pharmacological and genetic loss-of-function screening
(Corsello et al., 2020; Behan et al., 2019; Meyers et al., 2017), and drug perturbation
profiling (Subramanian et al., 2017) have generated a wealth of molecular data
characterizing large numbers of cell lines and their responses to perturbations. Many
computational approaches have been developed to leverage these molecular data for
drug sensitivity predictions and/or drug repurposing (Parca et al., 2019; Wei et al., 2019;
Nguyen et al., 2017; Gao et al., 2021; Napolitano et al., 2018; Fang et al., 2021; Rivas-
Barragan et al., 2020; Domingo-Fernández et al., 2022; Szalai et al., 2019; Napolitano
et al., 2015; Chan et al., 2019; Emon et al., 2020; Rampášek et al., 2019; He et al.,
2023; Jia et al., 2016; Han et al., 2021; Chen et al., 2022; Wu et al., 2022), and these
efforts have successfully identified drugs for a wide variety of diseases, including HIV
(Fang et al., 2021), osteoporosis (Brum et al., 2015), diabetes (Zhang et al., 2015), and
cancer (Gonçalves et al., 2020; Tsoi et al., 2018; Viswanathan et al., 2017). Despite
these successes, many patients remain ineligible for targeted therapies, including over
80% of cancer patients (Haslam et al., 2021). Furthermore, only about half of eligible
cancer patients are responsive to targeted therapy, emphasizing the need for improved
drug discovery and repurposing methods.


5
One common drawback of many drug repurposing tools is that they output a long list
of candidate drugs with limited information about how the top candidates are related.
For example, the gene2drug algorithm (Napolitano et al., 2018) returns a ranked list of  >
1,300 drugs without any information about molecular targets or pathways of these
drugs. This complicates efforts to prioritize drugs on the list for validation because
researchers must consider many drugs targeting different molecular pathways with the
caveat that some targeted therapies may not actually inhibit their intended target (Lin et
al., 2019). Therefore, given a list of candidate drugs, we reasoned that grouping drugs
with similar mechanisms of action (MOAs) into a “set” and then statistically evaluating
the enrichment of the drug set in the list would increase on-target signal and reduce off-
target effects compared to evaluating drugs on an individual basis. Here, MOA refers to
both the biological pathway targeted and the direction of action of each drug (e.g.,
“EGFR inhibitor”). Our approach, called drug mechanism enrichment analysis (DMEA),
is an adaptation of the popular gene set enrichment analysis (GSEA) algorithm
(Subramanian et al., 2005) in which drugs, rather than genes, are grouped into sets
based on annotated MOAs. Each drug set is then statistically evaluated against a
background of all other drug sets. If multiple drugs which share a common MOA are all
highly ranked candidates, then this indicates that the identified MOA is more likely to be
a true on-target sensitivity.

Notable alternatives to our approach for analyzing enriched MOAs in drug lists
include the Connectivity Map (CMap) L1000 Query (Subramanian et al., 2017),
DrugEnrichr (Kuleshov et al., 2016, 2019; E. Y. Chen et al., 2013), Drugmonizome

6
(Kropiwnicki et al., 2021), DrugPattern (Huang et al., 2018), and drug set enrichment
analysis (DSEA) (Napolitano et al., 2015). However, these tools have several key
limitations (Fig. 1) including that they: (1) can only query preselected public data sets
(e.g., CMap’s L1000 transcriptional database); (2) have limited statistical rigor (e.g., lack
of p values with CMap; lack of multiple hypothesis correction with DSEA; lack of
permutation-based metrics with DrugEnrichr, Drugmonizome, and DrugPattern); (3)
accept only one type of unranked input list (i.e., gene symbols for CMap L1000 Query;
drug names for DrugEnrichr, Drugmonizome, DrugPattern, and DSEA); and (4) do not
generate plots of MOA-specific results. In addition, we note that DSEA queries gene
sets (e.g., gene ontology terms like “cellular protein localization”) rather than drug MOAs
(e.g., “HDAC inhibitor”). To address these shortcomings, we sought to make DMEA
compatible with any data set or drug repurposing algorithm, maintain the statistical rigor
of GSEA, and generate plots of both the overall results (i.e., volcano plot of all MOA
normalized enrichment scores) and MOA-specific results (i.e., mountain plots).


7

Figure 1. DMEA is more flexible and statistically rigorous than other approaches to
evaluate drug MOA
The Venn diagram compares our method, DMEA, with the Connectivity Map (CMap)
L1000 query of gene expression signatures (Subramanian et al., 2017) and the
DrugEnrichr (Kuleshov et al., 2016, 2019; E. Y. Chen et al., 2013) and Drugmonizome
methods (Kropiwnicki et al., 2021).

Furthermore, to address a lack of web-accessible tools to predict selectively toxic
drugs based on an input gene signature, we included a feature to pair DMEA with a
simple molecular classification method (i.e., weighted gene voting, WGV (Golub et al.,

8
1999); see Methods). Although the CMap L1000 Query also accepts input gene
symbols to rank drug MOA (Subramanian et al., 2017), it is limited to use of the CMap
L1000 expression database and cannot accept more than 150 input gene symbols or
consider gene ranks. Similar to the CMap L1000 Query, the gene2drug web tool also
ranks individual drugs based on gene inputs using a gene set-based analysis of the
original CMap gene expression database (e.g., gene ontology terms) without
considering input gene ranks (Napolitano et al., 2018). In addition to sharing these
limitations of the CMap L1000 Query, gene2drug does not evaluate drugs in terms of
MOA. On the other hand, the CMap PRISM Query does rank drugs based on selective
toxicity, but it only accepts cell line names as input features, restricting its applicability to
cell lines present in the CMap database, and does not evaluate drug MOA. In contrast,
DMEA can accept an input gene signature with any number of ranked genes and
requires their directionality to evaluate drug MOA based on selective toxicity.

In summary, DMEA can help researchers better prioritize potential drug treatments
by aggregating results across many drugs which share MOAs to identify global trends.
By quantifying the coordinated enrichment of drugs annotated with the same MOA and
normalizing scores across a large background of drug MOAs, DMEA can improve on-
target prioritization of candidates for drug repurposing. DMEA is publicly available as a
web application or an R package at https://belindabgarana.github.io/DMEA.


9
Methods
Drug mechanism enrichment analysis (DMEA)
DMEA tests whether drugs known to share a MOA are enriched in a rank-
ordered drug list. DMEA can be applied to any rank-ordered list of drugs with
annotations for at least two MOAs. For a drug MOA to be evaluated, at least six drugs
(or the minimum number of drugs per MOA set by the user) must be annotated with that
MOA and each drug must be ranked by a nonzero numeric value. DMEA uses the same
algorithm as GSEA (Subramanian et al., 2005) but applies it to sets of drugs, rather
than genes, to identify drug MOAs which are overrepresented at either end of the input
rank-ordered drug list (further detail below). If a drug MOA is positively enriched, then
drugs annotated with that MOA are overrepresented at the top of the list. Conversely, if
a drug MOA is negatively enriched, then drugs which share that MOA annotation are
overrepresented at the bottom of the list.

Specifically, for each MOA, DMEA calculates an enrichment score (ES) as the
maximum deviation from zero of a running-sum, weighted Kolmogorov–Smirnov-like
statistic. The p value is estimated using an empirical permutation test wherein drugs are
randomly assigned MOA labels in 1,000 independent permutations to calculate a
distribution of null enrichment scores (ES_null); the p value is then calculated as the
percentage of same-signed ES_null equal to or greater than the ES divided by the
percentage of same-signed ES_null. The normalized enrichment score (NES) is then
calculated by dividing the ES by the mean of the same-signed portion of the ES_null
distribution. Finally, the q value or false discovery rate (FDR) is calculated as the

10
percentage of same-signed NES in the null distribution (i.e., NES_null) with NES equal
or greater to the observed NES divided by the percentage of same-signed NES equal or
greater. We use a significance threshold of p < 0.05 and FDR < 0.25 by default per the
recommendation for GSEA, but this FDR cutoff can be customized by the user. Given a
rank-ordered drug list, DMEA generates (1) enrichment results for all tested drug MOAs;
(2) a volcano plot summarizing the NES and -log10(p value) for all tested drug MOAs;
and (3) mountain plot(s) for individual drug MOA(s) which pass the given FDR cutoff
(Fig. 2).


Figure 2. Overview of drug mechanism enrichment analysis
DMEA is an adaptation of GSEA which analyzes a rank-ordered drug list to identify drug
MOAs that are overrepresented at either end of the input drug list. Given a rank-ordered
drug list where drugs have been annotated with known MOAs, DMEA runs an enrichment
analysis for each individual MOA. After calculating p values and FDR q values, DMEA
outputs (1) enrichment results for all tested drug MOAs; (2) a volcano plot summarizing
the NES and -log10(p value) for all tested drug MOAs; and (3) mountain plot(s) for
individual drug MOA(s) which pass the FDR cutoff.

Simulation study of DMEA
To evaluate the sensitivity of DMEA, we first simulated a rank-ordered drug list
by randomly assigning values from a normal distribution (μ  =  0, σ  =  0.5) for 1,351 drugs
with MOA annotations in the PRISM drug screen. Next, a number of drugs, X, were

11
randomly sampled as a synthetic drug set and their rank values were selected from a
shifted normal distribution (μ  =  Y, σ  =  0.5); the size of the synthetic drug set, X, was
varied from 5 to 50 drugs, and the perturbation value Y was varied from -1 to  +1. This
rank-ordered drug list was then analyzed by DMEA to determine the enrichment of the
synthetic drug set relative to the known drug MOA sets provided by the PRISM drug
screen. For each pair of X and Y values, the simulation was repeated 50 times to
assess reproducibility (i.e., the ability of DMEA to consistently detect a true difference
between the synthetic drug set and the background drug sets determined by MOA
annotations from PRISM).

CMap L1000 query
The Connectivity Map (CMap) web portal (https://clue.io) (Subramanian et al.,
2017) allows users to query the L1000 gene expression database using 10 to 150 input
up- and down-regulated gene IDs. The output is a normalized connectivity score which
indicates the similarity between the query and differentially expressed gene sets
induced by drug treatments. A positive score indicates similarity between the query and
the perturbagen signature, whereas a negative score indicates dissimilarity. Specifically,
we used the “query.gct” file from the zip file output by the CMap L1000 Query (found
within their “arfs/TAG” folder), including the MOA annotations provided in the file. Since
this file includes results for all cell lines in the L1000 database as well as information for
quality control, we removed any scores indicated to be low quality and averaged scores
across cell lines for each drug with MOA annotations. Here, we used example data sets
from the CMap web portal, including: (1) GSE32547, HUVEC cells treated with the

12
HMGCR inhibitor pitavastatin (1 μM, 4 h) or DMSO (Maejima et al., 2014); (2)
GSE35230, A375 melanoma clones treated with the MEK inhibitor GSK212 (30 nM,
24 h) or DMSO (Greger et al., 2012); (3) GSE14003, JEKO1 cells treated with the
proteasome inhibitor bortezomib (10 h) or untreated (Wang et al., 2009); (4) GSE28896,
human CD34
+
cells treated with the glucocorticoid agonist dexamethasone (24 h) or
untreated (Narla et al., 2011); and (5) GSE33643, A2058 cells treated with the
PI3K/MTOR inhibitor BEZ235 (3 doses at 24 h) or DMSO (Brachmann et al., 2012a).
We also used the up- and down-regulated biomarkers from a proteomic signature of
senescence in primary human mammary epithelial cells (HMECs) (Delfarah et al.,
2021). To compare DMEA’s results to CMap’s MOA enrichment results, we used the
“gsea_result.gct” file found within the “gsea/TAG/arfs/NORM_CS” folder. We compared
the results for chemical perturbagens combined across all cell lines, specifically
“cell_iname = -666” with pert_type = “TRT_CP”, for either set_type = “PCL” or
set_type = “MOA_CLASS”.

CMap PRISM query
The Connectivity Map (CMap) web portal (https://clue.io) (Subramanian et al.,
2017) also allows users to query PRISM viability data (Corsello et al., 2020) for 3 to 489
input cell line IDs classified as having “UP” or “DOWN” phenotypes. The query outputs
a normalized connectivity score ranking the drugs based on their toxicity towards the
“UP” versus “DOWN” cell lines. A positive score indicates toxicity towards “UP” cell
lines, whereas a negative score indicates toxicity towards “DOWN” cell lines or lack of
toxicity towards “UP” cell lines. In particular, we used the “ncs.gct” file from the zip file

13
output by the CMap PRISM Query, including the MOA annotations provided in the file.
Again, we only considered drugs with MOA annotations. Here, we used examples
provided by the CMap web portal including: (1) cell lines with the EGFR activating
mutation p.E746_A750del (i.e., “UP” cell lines: NCIH1650, PC14, and HCC827); (2) cell
lines with high expression of PDGFRA (i.e., “UP” cell lines: 42MGBA, A204, A2780,
G292CLONEA141B1, G402, GB1, HS618T, KNS42, LMSU, MG63, MON, NCIH1703,
SBC5, SKNAS, SNU685, SW579, U118MG, U251MG, and YH13); and (3) cell lines
sensitive to the HMGCR inhibitor lovastatin (i.e., “UP” cell lines: HUH28, SNU1079,
MG63, LOXIMVI, MDAMB231, SF295, SNU1105, YKG1, ACHN, HCT15, SNU423,
SNU886, CALU1, HCC4006, HCC44, HCC827, NCIH1915, NCIH661, NCIH838, PC14,
RERFLCMS, SQ1, SW1573, KYSE150, A2780, COV434, JHOM1, MCAS, KP3,
SW1990, MSTO211H, YD15, HS944T, MDAMB435S, MELJUSO, A204, HT1080,
RH30, LMSU, FTC238, YD8, 5637, and AGS).

Weighted gene voting (WGV)
To calculate a molecular classification score for cell lines based on external
molecular signatures, we used weighted gene voting (WGV) (Golub et al., 1999). The
WGV score is the dot product between an external gene signature of interest and
normalized RNAseq expression values for 327 adherent cancer cell lines from the
Cancer Cell Line Encyclopedia (CCLE, version 19Q4) (Ghandi et al., 2019). In other
words, the WGV score for each cell line is the sum across all genes available in both
the input gene signature and CCLE RNAseq data set, where each gene’s value is the
product between their gene signature ranking and CCLE RNAseq normalized

14
expression value. This WGV score ranks each cell line based on the similarity of its
gene expression to that of the input gene signature, such that cell lines with expression
more similar to the positive phenotype of the gene signature are more positively ranked
and cell lines with expression more similar to the negative phenotype of the gene
signature are more negatively ranked. In this study, we analyzed four independent
transcriptomic signatures, three derived from data sets for intrinsic resistance to EGFR
inhibitors and one derived from a data set for acquired resistance to a RAF inhibitor. For
each transcriptomic data set, we used the R package limma (Ritchie et al., 2015) to
perform an eBayes statistical analysis for differential expression comparing sensitive
and resistant samples. Then, the top 500 genes based on |log2(fold-change)| with q
value  <  0.05 were used for WGV (with the log2(fold-change) being the gene “weight” or
rank value).

For gene signatures of EGFR inhibitor sensitivity, we used data sets GSE12790
(Hoeflich et al., 2009), GSE31625 (Balko et al., 2006), and Coldren et al. (Coldren et al.,
2006). In GSE12790, transcriptomic profiles were provided for breast cancer cell lines
classified as either sensitive (EC50  <  1 µM: HDQ-P1, CAL85-1, and HCC1806) or
resistant to erlotinib (EC50  >  10 µM: CAL-51, CAL-120, MDA-MB-231, BT-20,
HCC1569, EFM-192A, HCC1954, MDA-MB-453, BT474, HCC1428, T47D, ZR-75-1,
KPL-1, BT-483, MDA-MB-415, HCC1500, CAMA-1, and MCF7). For GSE31625, we
used 17 transcriptomic profiles from 3 non-small cell lung cancer cell lines sensitive
(H1650, H3255, and PC-9) and 12 profiles of 2 cell lines resistant to erlotinib (A549 and
UKY-29). Finally, based on classifications from Coldren et al., we used CCLE RNAseq

15
profiles of 5 non-small cell lung cancer cell lines sensitive (NCIH1650, HCC95,
NCIH1975, NCIH1648, and NCIH2126) and 7 cell lines resistant to gefitinib (NCIH520,
NCIH460, NCIH1299, HCC44, A549, NCIH1703, and HCC15). For a gene signature of
RAF inhibitor sensitivity, we used data set GSE66539 with paired biopsy samples of
melanoma from 3 patients before vemurafenib treatment and after emergence of
resistance to vemurafenib (Kemper et al., 2016).

DMEA using WGV molecular classification scores
To identify drug MOAs with selective toxicity towards cells represented by an
input gene signature, DMEA can be used in combination with a molecular classification
method such as WGV, correlations, and large public databases for gene expression and
drug screens. To do this, we accessed the Cancer Cell Line Encyclopedia version 19Q4
for RNAseq data and calculated WGV scores for 327 adherent cancer cell lines using
external molecular signatures. To avoid overfitting, we did not include WGV scores from
any cell lines that had been used to generate the external molecular signature. Next, we
calculated the Pearson correlation between the WGV scores and PRISM drug
sensitivity scores (i.e., area under the curve (AUC) values for cell viability as a function
of drug concentration) for each drug (Corsello et al., 2020) using data from the most
recent PRISM screen available (e.g., HTS002, MTS005, MTS006, and MTS010). Lastly,
drugs were ranked by the Pearson correlation coefficient, and the rank-ordered drug list
was analyzed by DMEA using the MOA annotations provided in the PRISM data set.


16
Simulation study of DMEA using WGV molecular classification scores
For 200 synthetic cell lines, we sampled drug sensitivity scores for 1,351 drugs
with MOA annotations in the PRISM drug screen from a bimodal mixture of two normal
distributions (μ1  =  0.83, σ1  =  0.08 and μ2  =  1.31, σ2  =  0.08) with the lower distribution
containing 72% of all drugs. This distribution was chosen to reflect the distribution of the
PRISM drug sensitivity data (i.e., AUC) (Corsello et al., 2020). Additionally, we
simulated gene expression for each cell line by sampling from a normal distribution with
a mean (μ) of 0 and standard deviation (σ) of 0.5. This distribution was chosen to reflect
the distribution of the normalized CCLE RNAseq data (Ghandi et al., 2019).

To introduce a synthetic association between gene expression and drug
sensitivity, we randomly sampled a synthetic gene set of 25 genes and a synthetic drug
set of 10 drugs. Next, expression values for the synthetic gene set and sensitivity
scores for the synthetic drug set were each sampled from a shifted distribution, where
the magnitude of the shift for each synthetic cell line is determined by a perturbation
value ranging from 0 (no perturbation) to 0.1. For example, for a perturbation value of
0.1, the mean gene expression for the 25 perturbed genes in cell line 1 was μ  =  -0.1,
and the mean sequentially increased by 0.001 for cell lines 2–200; similarly, the mean
drug sensitivity of cell line 1 to the 10 perturbed drugs was shifted by -0.1, and this shift
value sequentially increased by 0.001 for cell lines 2–200. This created a gradient of
perturbations in the 200 cell lines, meaning cell line 1 had the largest negative
perturbation and cell line 200 had the largest positive perturbation. Then, we calculated
WGV scores for each cell line by taking the dot product of the expression values of the

17
synthetic gene set and the difference in average gene expression between the top and
bottom 10 percent of cell lines (i.e., gene weights from cell lines 181–200 which had the
highest mean expression versus cell lines 1–20 which had the lowest mean
expression). Afterwards, we calculated the Pearson correlation between the WGV and
drug sensitivity scores for each of the 1,351 drugs in the synthetic data set. Finally,
drugs were ranked by their Pearson correlation coefficient and DMEA was performed to
measure the enrichment of the synthetic drug set relative to the background drug sets
which were determined by the MOA annotations in the PRISM drug screen. To assess
reproducibility, this entire process was repeated 50 times.

DrugEnrichr and Drugmonizome
To benchmark DMEA against DrugEnrichr and Drugmonizome, we input the top
50 or 100 drugs from each analysis detailed above into DrugEnrichr and
Drugmonizome. Though DMEA evaluates the full list of ranked drugs from each
analysis, DrugEnrichr and Drugmonizome are designed to accept just one input drug
set, so we chose to use the top 50 or 100 drugs as input drug sets to match the orders
of magnitude of the example inputs provided on their websites. The top drugs were
positively ranked in analyses of drug rank lists from the CMap Query tools or negatively
ranked in analyses where drugs were ranked based on Pearson correlations between
their PRISM drug sensitivity scores (i.e., AUC) and CCLE WGV scores. We recorded
their rankings of MOA sets from the PRISM drug repurposing hub based on each tool’s
default ranking when viewing their results in tabulated form on the web.


18
Cell culture
HMEC cells were purchased from Thermo Scientific and cultured in M87A
medium (50% MM4 medium and 50% MCDB170 supplemented with 5 ng/ml EGF,
300 ng/ml hydrocortisone, 7.5 µg/ml insulin, 35 µg/ml BPE, 2.5 µg/ml transferrin, 5 µM
isoproterenol, 50 µM ethanolamine, 50 µM o-phosphoethanolamine, 0.25% FBS, 5 nM
triiodothyronine, 0.5 nM estradiol, 0.5 ng/ml cholera toxin, 0.1 nM oxytocin, 1% anti-anti,
and no AlbuMax) in atmospheric oxygen. Glucose and glutamine-free DMEM was
purchased from Corning (90-113-PB), Ham’s F12 was purchased from US Biological
(N8542-12), and MCD170 medium was purchased from Caisson Labs (MBL04).
Glucose and glutamine were added to the media at the appropriate concentration for
each media type. At each passage, cells were lifted with TrypLE at 80–90% confluency
and seeded at a density of 2.3  ×  10
3
/cm
2
.

Cell viability experiments
Proliferating HMECs (PD  ~  12) were seeded at a concentration of 2.1  ×  10
3
/cm
2
or
9.5  ×  10
3
/cm
2
for DMSO and triapine treatment, respectively. The following day, cells
were treated with DMSO (vehicle control) or 2 µM triapine for 3 days. The cells were
counted and then treated with either DMSO (vehicle control), dacomitinib, AZD8931, or
navitoclax at 100 nM or 500 nM for 3 days. Cell viability and live cell number were
measured with trypan blue assay using a TC20 automated cell counter (Bio-Rad).
Chemical inhibitors were from Sigma (triapine) or MedChemExpress (dacomitinib,
AZD8931, and navitoclax).


19
Implementation
DMEA: web application
Using our web application, researchers can either input a drug rank list or gene
signature to identify enriched drug MOA without any programming knowledge required.
With both input types, the outputs contain: (1) the processed input (after any averaging
across duplicate input features if applicable), (2) MOA annotations used (either provided
by the user or the default PRISM drug annotations provided on our GitHub repository
(MOA_gmt_file_n6_no_special_chars.gmt), (3) the results of our DMEA analysis, (4)
any drug sets which were removed because they were not represented by at least 6
drugs or the minimum set by the user, (5) any unannotated drugs which could not be
matched into drug sets, (6) a volcano plot with an overview of the normalized
enrichment scores and -log10(p values) for all evaluated drug MOA with significant
enrichments highlighted in red, and (7) mountain plots for each significantly enriched
drug MOA so that users can visually confirm that most of the drugs for these MOA were
ranked as strong candidates. All the analyses of non-simulated data sets in this study
can easily be replicated as examples from a drop-down menu on our web application.
These example inputs are available on our GitHub repository (Examples) and described
below as well as on our website’s “How to Use” page
(https://belindabgarana.github.io/DMEA/howtouse.html), which also includes a video
tutorial.

Drug rank lists

20
The input drug rank list can either be a CSV file with headers or direct output
from a CMap Query tool (GCT file). If inputting a drug rank list as a CSV file, the user
can either provide MOA annotations in the third column separated by the “|” character in
each row or rely on our default MOA annotations derived from the PRISM drug screen
database; the first column must contain drug names and the second column must
contain signed (i.e., nonzero) drug ranks. If there are multiple ranks provided for each
drug, the user must select the checkbox option to average results across drugs. In the
advanced settings, there is also an option to convert drug synonyms if no MOA
annotations are provided, which is enabled by default. For transparency, the CSV file
used to convert drug synonyms is uploaded on our GitHub repository
(PRISM_drug_synonyms.csv) and contains PRISM drug names, PubChem CIDs, and
synonyms found on PubChem (https://pubchem.ncbi.nlm.nih.gov). The conversion of
drug synonyms feature allows more input drugs to be matched into PRISM drug sets
and even allows input of unique PubChem CIDs as drug names to avoid issues with
different naming systems. Drugs which are not matched into drug sets (i.e.,
unannotated drugs) even after drug synonym conversion is performed are still
considered as background drugs in our analysis and also output as a CSV file to allow
users to repeat the analysis after either converting their drug names manually (e.g., with
unique PubChem CIDs using the PubChem search tool
at https://pubchem.ncbi.nlm.nih.gov/) or adding any known MOA annotations for these
drugs. If a CMap L1000 drug rank list is input, the positively enriched drug MOA may
induce similar expression to the “UP” phenotype genes and vice-versa. If a CMap
PRISM drug rank list is input, then the positively enriched drug MOA may be selectively

21
toxic to the “UP” phenotype cell lines and vice-versa. If a custom drug rank list is input,
then the positively enriched drug MOA are overrepresented at the top of the input drug
rank list and vice-versa.

Gene signatures
The input gene signature should be formatted as a CSV file with headers where
the gene symbols are in the first column and the gene ranks are in the second column.
If there are multiple ranks provided for each gene, the user must select the checkbox
option to average results across genes. In the advanced settings, there is an option to
use either HUGO Gene Nomenclature Committee (HGNC)-approved gene symbols to
analyze your input gene signature (Seal et al., 2023) or gene symbols used in the CCLE
19Q4 release. Our software also outputs a CSV file containing any genes which were
not matched to the CCLE 19Q4 data set so that users can repeat their analysis after
either using the HGNC multi-symbol checker web tool
(https://www.genenames.org/tools/multi-symbol-checker/) to convert their gene symbols
to approved symbols or searching for the CCLE 19Q4 versions of their gene symbols in
the CSV file on our GitHub repository which includes approved symbols, aliases,
previous symbols, and HGNC IDs for each gene symbol from the CCLE 19Q4 release
(CCLE_gene_symbols_20230404.csv). Using the input gene signature, DMEA runs
weighted gene voting (WGV) to rank 372 adherent cancer cell lines in the CCLE 19Q4
RNAseq data set based on their expression of the gene signature and then correlates
their WGV scores with their drug sensitivity scores (i.e., AUC) across 1,351 drugs with
MOA annotations in the PRISM drug screen. For transparency in this process, our

22
software outputs the WGV scores, correlation results, and scatter plots with the
correlations for each drug in addition to the other outputs mentioned above. With this
analysis, negatively enriched drug MOA may be selectively toxic to samples with higher
expression of positively ranked genes in the input gene signature and vice-versa.

DMEA: R package
Our R package is available on GitHub to allow many automated DMEA analyses
to be run in sequence and querying of any database beyond just the CCLE RNAseq
19Q4 release and the PRISM drug screen. In brief, if analyzing an input drug rank list,
users can run the “drugSEA” function; if analyzing an input gene signature, users can
run the “DMEA” function. These functions accept the same inputs except formatted as
data frames in the R software environment as well as custom MOA annotations and, in
the case of the DMEA function, also expression and drug sensitivity data frames
representing the same biological samples (e.g., cell lines). The same outputs are
available for each input type as described above for the web application. Documentation
is available on our GitHub repository (https://github.com/BelindaBGarana/DMEA),
including a vignette and man pages with installation instructions and runnable
examples, and also on our website’s “How to Use” page
(https://belindabgarana.github.io/DMEA/howtouse.html).


23
Results
DMEA identifies an enriched drug MOA in simulated data
To evaluate the ability of DMEA to identify the enrichment of drug sets, we tested
it on a normally distributed, synthetic ranked list of 1,351 drugs (see Methods). For a
randomly sampled set of drugs ranging in size from 5 to 50 drugs, we shifted these
drugs’ rankings by a perturbation value ranging from -1 to 1. Next, we ran DMEA using
the full rank list of drugs to assess enrichment of the synthetic drug set. This process
was repeated 50 times for each synthetic drug set size, after which the average
normalized enrichment score (NES) and the percentage of replicates with significant
enrichment of the synthetic drug set were visualized as heatmaps (Fig. 3).


Figure 3. Sensitivity analysis of DMEA using synthetic data
Synthetic rank-ordered drug lists were generated with varying perturbations (y-axis) of
different drug set sizes (x-axis), then analyzed by DMEA (see Methods). For each
combination of drug set size and perturbation value, 50 replicates were
performed. A) Heatmap showing the average DMEA NES for the perturbed drug
set. B) Heatmap showing the percent of DMEA replicates with FDR q value  <  0.25 for the
perturbed drug set.

24

As expected, we observed no false discovery except for very small drug sets
(i.e., when evaluating a set of 5 drugs we observed 2% of replicates were falsely
enriched). As the magnitude of the perturbation was increased or decreased, the
average NES of the synthetic drug set increased or decreased, respectively. Likewise,
the percentage of replicates passing the significance threshold of p  <  0.05 and
FDR  <  0.25 increased as the magnitude of the perturbation increased. These results
demonstrate that DMEA can successfully identify an enriched set of drugs in simulated
data.

DMEA identifies similar MOAs based on gene expression connectivity scores
Next, we sought to test whether DMEA could identify enriched drug MOAs in
rank-ordered drug lists generated by the Connectivity Map (CMap), (Subramanian et al.,
2017) a popular tool for drug discovery that contains  > 1 million gene expression
signatures measured using L1000, a reduced representation transcriptomic profiling
method. Specifically, we analyzed example data sets from the CMap L1000 Query tool
to identify perturbagen signatures that are similar or dissimilar to an input gene set.
First, we used a gene expression signature from HUVEC cells treated with pitavastatin,
an inhibitor of 3-hydroxy-3-methylglutaryl-CoA reductase (HMGCR) (Maejima et al.,
2014), to rank 3,868 drugs based on the similarity of their effects on gene expression.
Because pitavastatin itself was found in the list of 3,868 drugs, one might have
expected it to be the top-ranked, most similar drug produced by this analysis, but in fact
it ranked 24th out of 3,868 drugs (Appendix A: Fig. S1A). In contrast, analysis of the

25
rank-ordered list of drugs using DMEA identified the HMGCR MOA as the only
significant similar MOA (Fig. 4A). This demonstrates that analysis of MOAs by DMEA
can generate clearer and more statistically significant results than analysis of individual
drugs in results from the CMap L1000 Query.


26

Figure 4. DMEA identifies similar MOAs based on gene expression connectivity
scores
Rank-ordered drug lists were generated by querying the CMap L1000 gene expression
perturbation signatures and then analyzed by DMEA. A) HUVEC cells treated with the
HMGCR inhibitor pitavastatin (Maejima et al., 2014), B) A375 melanoma clones treated
with the MEK inhibitor GSK212 (Greger et al., 2012), and C) JEKO1 cells treated with the
proteasome inhibitor bortezomib (Wang et al., 2009). Volcano plots summarizing the NES
and −log10(p value) for all tested drug MOAs and mountain plots of the expected MOAs

27
are shown. Red text indicates MOAs with p value  <  0.05 and FDR  <  0.25. For each
mountain plot, the inhibitors with the most positive connectivity scores are highlighted.

Next, we tested a gene expression signature from A375 melanoma cells treated
with the MEK inhibitor GSK212 (Greger et al., 2012). Again, DMEA correctly identified
that MEK inhibitors were the most similar MOA in the rank-ordered list of drugs (Fig. 4B;
Appendix A: Fig. S1A). In this case, comparison to the single-drug rankings was not
possible because the L1000 database does not contain the true drug treatment,
GSK212 (Appendix A: Fig. S1A). Subsequently, we analyzed a gene expression
signature from JEKO1 mantle cell lymphoma cells treated with the proteasome inhibitor
bortezomib (Wang et al., 2009). DMEA again accurately found that proteasome
inhibitors were the most similar MOA (Fig. 4C) and DMEA’s MOA-level ranking (#1) was
improved upon the single-drug ranking of the true drug treatment, bortezomib
(#14/3,868) (Appendix A: Fig. S1A). Finally, we used DMEA to test data sets from
human CD34
+
cells treated with the glucocorticoid agonist dexamethasone (Narla et al.,
2011) and A2058 melanoma cells treated with the PI3K/mTOR inhibitor BEZ235
(Brachmann et al., 2012a). In both cases, DMEA correctly identified the expected MOA
as significantly enriched (glucocorticoid receptor agonist and PI3K/mTOR inhibitor
MOAs, respectively) (Appendix A: Fig. S2) and DMEA’s MOA-level rankings were
improved upon CMap’s individual drug rankings of the true drug treatments (Appendix
A: Fig. S1A). Taken together, these results show that DMEA can correctly identify
enriched MOAs in rank-ordered lists of drugs generated by the CMap L1000 Query and
that the MOA-level rankings of the true drug treatments are improved compared to the
single-drug rankings.


28
Next, we compared DMEA’s MOA-level results to those of the CMap L1000
Query (found in an output sub-folder generated by CMap L1000 Query called
“gsea/TAG/arfs/NORM_CS”). Like our DMEA results, CMap’s MOA-level rankings
revealed the expected MOA as the top-ranked MOA in all cases except for
glucocorticoid receptor agonists and PI3K inhibitors which were not found in the L1000
output (Appendix A: Fig. S1A). We also compared our DMEA results to the CMap
L1000’s perturbagen classes (PCLs), which are derived from MOA sets but exclude
drugs which do not fit the overall trend of the MOA (Subramanian et al., 2017). Again,
CMap’s PCL rankings were similar to that of DMEA (Appendix A: Fig. S1A). Thus,
DMEA and the CMap L1000 generate similar MOA-level rankings. However, in contrast
to DMEA, the CMap L1000 MOA-based analysis has less statistical rigor (i.e.,
no p values provided by CMap) and does not generate any plots of the overall and
MOA-specific results (e.g., volcano or mountain plots).

We also compared our results with those of DrugEnrichr and Drugmonizome
(Appendix A: Fig. S1). Since DrugEnrichr and Drugmonizome are only designed to
evaluate one input drug set, whereas DMEA and the CMap L1000 Query evaluate the
full list of drugs, we input the top 50 and 100 positively ranked drug names into
DrugEnrichr and Drugmonizome. DrugEnrichr and Drugmonizome were only able to
evaluate the expected drug MOA for three out of five CMap L1000 examples (HMGCR
inhibitor pitavastatin, mTOR/PI3K inhibitor BEZ235, and glucocorticoid receptor agonist
dexamethasone). Drugmonizome performed similarly in ranking the expected MOA to
both DMEA and the CMap L1000 Query tool. DrugEnrichr also performed similarly with

29
the mTOR/PI3K inhibitor example but had lower rankings for the glucocorticoid receptor
agonist and HMGCR inhibitor examples (Appendix A: Fig. S1A). However, it is
important to consider that we cannot make a fair comparison since DrugEnrichr and
Drugmonizome were only evaluating drug MOA based on a small fraction of the drug
lists compared to DMEA and the CMap L1000 Query.

DMEA identifies selectively toxic MOAs based on cell viability connectivity scores
To evaluate if DMEA can identify enriched MOAs in a different type of rank-
ordered drug list, we used the CMap PRISM Query tool to query data from the PRISM
drug repurposing database (Corsello et al., 2020). Given an input list of cell line names,
the CMap PRISM Query generates a list of  ~1,200 drugs ranked by normalized
connectivity scores which represent the predicted sensitivity of the input cell lines to
each drug. The higher the normalized connectivity score, the more toxic the drug is
predicted to be for the input cell lines. Again, we analyzed example data sets from the
CMap PRISM Query tool to test DMEA, including: (1) cell lines with the activating EGFR
mutation p.E746_A750del (n  =  3), (2) cell lines with high expression
of PDGFRA (n  =  19), and (3) cell lines with sensitivity to the HMGCR inhibitor lovastatin
(n  =  43). As hypothesized, DMEA identified EGFR inhibitors (Fig. 5A), PDGFR inhibitors
(Fig. 5B), and HMGCR inhibitors (Fig. 5C), respectively, as significantly positively
enriched in these rank-ordered drug lists. DMEA also improved upon the rank of the
true drug sensitivity for the HMGCR inhibitor lovastatin (#1 with DMEA’s MOA-level
rankings versus #3 in the single-drug rankings); DrugEnrichr and Drugmonizome also
improved on the rank of the true drug sensitivity, though they evaluated fewer drug

30
MOA sets because their analyses only considered the top 50 or 100 positively ranked
drugs instead of the full list of evaluated drugs (Appendix A: Fig. S1B). Altogether, these
examples demonstrate that DMEA can identify enriched MOA in rank-ordered lists of
drugs generated by CMap Query of the PRISM drug screen and that the MOA-level
ranking of the true drug sensitivity is higher than that of the single-drug ranking.


31

Figure 5. DMEA identifies selectively toxic MOAs based on cell viability
connectivity scores
Rank-ordered drug lists were generated by querying the PRISM database with input cell
line sets characterized by A) the activating EGFR mutation p.E746_A750del, B) high
expression of PDGFRA, and C) sensitivity to the HMGCR inhibitor lovastatin. Volcano
plots summarizing the NES and −log10(p value) for all tested drug MOAs and mountain
plots of the expected MOAs are shown. Red text indicates MOAs with p value  <  0.05 and

32
FDR  <  0.25. For each mountain plot, the inhibitors with the most positive connectivity
scores are highlighted.

DMEA identifies selectively toxic MOAs based on molecular signatures
To offer a web-accessible method to identify selectively toxic drug MOAs based
on an input molecular signature (i.e., up- or down-regulated genes that characterize a
disease or cell type), we paired DMEA with a simple molecular classification method,
namely weighted gene voting (WGV) (Golub et al., 1999). Specifically, we: (1) used
WGV to classify adherent cancer cell lines in the CCLE database based on similarity to
the input gene signature; (2) correlated WGV scores with drug sensitivity scores (i.e.,
AUC) for each of 1,351 drugs in the PRISM database; and (3) ranked drugs by the
correlation coefficient of WGV scores and drug AUC values (Appendix A: Fig. S3; see
Methods).

To test this approach, we first performed a simulation study. Specifically, we
simulated gene expression and drug sensitivity scores for 200 cell lines by randomly
sampling values from distributions that reflected the CCLE RNAseq and PRISM drug
sensitivity data. Next, to create a synthetic association between gene expression and
drug sensitivity, we perturbed a subset of the gene expression data and drug sensitivity
scores. We then ran 50 replicates to determine if DMEA could consistently identify
enrichment of the synthetic drug set in this simulated data. To visualize the results, we
plotted the average normalized enrichment score (NES) and the percent of replicates
which pass the significance threshold of FDR  <  0.25 as heatmaps (Appendix A: Fig. S4).
Importantly, when there was no perturbation in drug sensitivity (AUC), the tested drug

33
set was never significantly enriched (0% of replicates) regardless of the size of the
perturbation in RNA expression. This demonstrates that DMEA is not prone to false
positive results using this WGV-based approach. In addition, increasing the perturbation
in either RNA expression or drug sensitivity led to increased enrichment scores (i.e.,
average NES) and increased significance (i.e., higher percentage of significant
replicates). These results illustrate that DMEA can successfully identify associations
between gene expression and drug sensitivity with high reproducibility in simulated
data.

Next, we tested whether DMEA could successfully identify drug MOAs with
selective toxicity using published transcriptomic signatures of drug resistance. First, we
tested three different signatures of intrinsic resistance to EGFR inhibitors: (1) non-small
cell lung cancer (NSCLC) cell lines treated with erlotinib (GSE31625) (Balko et al.,
2006); (2) breast cancer cell lines treated with erlotinib (GSE12790) (Hoeflich et al.,
2009); and (3) NSCLC cell lines treated with gefitinib (Coldren et al.) (Coldren et al.,
2006). Notably, there was little overlap in the genes used for WGV (GSE12790 and
Coldren et al.share zero genes, GSE12790 and GSE31625 share 15 genes, and
GSE31625 and Coldren et al.share 19 genes; Fig. 6B). Despite the lack of overlap in
input gene signatures, all three DMEA analyses correctly identified EGFR inhibitors as
the top toxic drug MOA for the EGFR inhibitor-sensitive cancer cell lines (Fig. 6A,B,
Appendix A: Fig. S5). Again, DMEA’s MOA-level rankings were improved compared to
the single-drug rankings (#1 for EGFR inhibitors in all cases versus #16 for erlotinib
based on GSE31625, #10 for erlotinib based on GSE12790, and #13 for gefitinib based

34
on Coldren et al.); DrugEnrichr and Drugmonizome’s MOA rankings were also improved
over the individual drug rankings, though they evaluated fewer drug MOA sets since
they are only designed to accept one input drug set instead of the complete drug rank
list (Appendix A: Fig. S1B). In addition, DMEA revealed consistent results across all
three input gene signatures for drug MOAs identified as potentially toxic for EGFR
inhibitor-resistant cancer cell lines, including HMGCR and MDM inhibitors (Fig. 6B,
Appendix A: Fig. S5). These results support that DMEA can identify selectively toxic
drug MOAs given a molecular signature of intrinsic drug resistance and that DMEA’s
MOA-level rankings improve upon single-drug rankings of toxicity.

35

Figure 6. DMEA identifies selectively toxic MOAs based on external gene
expression signatures of intrinsic EGFR inhibitor resistance and acquired RAF
inhibitor resistance, respectively
Using gene expression signatures of intrinsic resistance to EGFR inhibition and acquired
resistance to RAF inhibition, we calculated WGV molecular classification scores for 327
adherent cancer cell lines in the CCLE database. For each signature, the WGV scores
were correlated with drug sensitivity scores (i.e., AUC) for 1,351 drugs from the PRISM
database. Drugs were then ranked by Pearson correlation coefficient, and DMEA was

36
performed to identify selectively toxic MOAs. A) DMEA analysis of GSE12790 (Hoeflich
et al., 2009) transcriptomic signature of intrinsic resistance to EGFR inhibitor erlotinib,
including a volcano plot of NES versus -log10(p value) for MOA evaluated where red text
indicates MOAs with p value  <  0.05 and FDR  <  0.25 and a mountain plot showing that
DMEA identified the EGFR inhibitor MOA as negatively enriched. The most negatively
correlated EGFR inhibitors are labeled along with their correlation
coefficients. B) Comparison of three transcriptomic signatures for intrinsic resistance to
EGFR inhibition analyzed using DMEA, including a Venn diagram showing the number of
shared genes among the signatures and a dot plot illustrating the consistency of MOA
enrichment across DMEA’s analyses. C) DMEA analysis of GSE66539 (Kemper et al.,
2016) transcriptomic signature of acquired resistance to RAF inhibitor vemurafenib,
including a volcano plot of NES versus -log10(p value) for MOA evaluated where red text
indicates MOAs with p value  <  0.05 and FDR  <  0.25 and a mountain plot showing that
DMEA identified the RAF inhibitor MOA as negatively enriched. The most negatively
correlated RAF inhibitors are labeled along with their correlation coefficients.

Next, we tested whether DMEA could identify selectively toxic drug MOAs given
a transcriptomic signature of acquired drug resistance. Specifically, we analyzed
RNAseq data from patient biopsies of BRAF-mutant melanoma before treatment with
the BRAF inhibitor vemurafenib and after the development of acquired resistance
(Kemper et al., 2016). Again, DMEA correctly identified the RAF inhibitor MOA as the
top toxic drug MOA for the samples collected prior to BRAF inhibitor treatment (Fig. 6C),
and the ranking of RAF inhibitors at the MOA-level (#1) was improved compared to the
ranking of vemurafenib alone (#35); this ranking was also improved when inputting the
top 50 negatively ranked drugs into DrugEnrichr or Drugmonizome (Appendix A: Fig.
S1B). Additionally, inhibitors of HDAC, EGFR, CHK, and SYK were identified as
possibly beneficial for tumors with acquired resistance to RAF inhibition. Conversely,
DMEA identified that drugs inhibiting MEK, MDM, nerve growth factor receptor, and
HMGCR may be toxic towards tumors which are sensitive to RAF inhibitors (Fig. 6C).
These results demonstrate that DMEA can amplify on-target signal to identify acquired

37
resistance in tumors and other drug MOAs which may be beneficial based on patient
biopsies.

DMEA identifies potential senescence-inducing and senolytic drug MOAs for primary
human mammary epithelial cells
Lastly, we sought to demonstrate how DMEA can be used as a discovery tool. As
an example, we analyzed our recently published proteomic signature of replicative
senescence, a stress-induced irreversible growth arrest associated with aging, in
primary human mammary epithelial cells (HMECs) (Delfarah et al., 2021). To highlight
the versatility of DMEA to either identify similar or selectively toxic drug MOAs, we
analyzed the same molecular signature using either the CMap L1000 Query or our
WGV-based approach to rank drugs, respectively (Fig. 7A). First, we performed a CMap
L1000 Query using the gene names for the up- and down-regulated proteins to predict
drug MOAs that could induce senescence in HMECs. Using the CMap results, DMEA
revealed positive enrichment for MOAs including proteasome, HDAC, HMGCR, and
MDM inhibitors (Fig. 7B), suggesting that treatment with drugs from these MOAs may
induce senescence in primary HMECs. Among the MOAs with significant negative
enrichments were Na/K-ATPase inhibitors and matrix metalloprotease inhibitors,
suggesting that these drug MOAs might antagonize senescence in primary HMECs.


38


39
Figure 7. DMEA identifies potential senescence-inducing and senolytic drug MOAs
for primary HMECs
A) Schematic detailing how the proteomic signature of replicative senescence in primary
HMECs (Delfarah et al., 2021) was used to identify either senescence-inducing or
senolytic drug MOAs. B) DMEA results for senescence-inducing drug MOAs. (Left)
Volcano plot of NES versus -log10(p value) for drug MOAs from DMEA. Red text indicates
MOAs with p value  <  0.05 and FDR  <  0.25. (Right) Mountain plot showing the positive
enrichment of the proteasome inhibitor MOA in the rank-ordered drug list of CMap L1000
connectivity scores. The proteasome inhibitors with the most positive connectivity scores
are highlighted. C) DMEA results for senolytic drug MOAs. (Left) Volcano plot of NES
versus -log10(p value) for drug MOAs from DMEA. Red text indicates MOAs
with p value  <  0.05 and FDR  <  0.25. (Right) Mountain plot showing the positive
enrichment of the EGFR inhibitor MOA in the rank-ordered drug list of correlation
coefficients. The EGFR inhibitors with the most positive correlation coefficients are
highlighted. D) The EGFR inhibitors dacomitinib and AZD8931 and the senolytic
compound navitoclax exhibited senolytic activity in HMECs. Proliferating HMECs
(PD  ~  12) were treated with DMSO or 2 μM triapine for 3 days to induce proliferating or
senescent phenotypes, respectively, as in our previous work (Delfarah et al., 2019).
Proliferating and senescent HMECs were then treated with DMSO (negative control),
100 nM/500 nM dacomitinib, 100 nM / 500 nM AZD8931, or 100 / 500 nM navitoclax for
3 days, after which cell viability and live cell number were measured by trypan blue
staining. The live cell number was normalized to the number of live cells present at the
time of drug treatment. * and ** represent p  <  0.05 and 0.01, respectively, compared to
the senescent DMSO control calculated by Student’s t-test.

Second, we analyzed the same proteomic signature of senescence with our own
WGV-based method to identify selectively toxic MOAs based on the input molecular
signature (Appendix A: Fig. S3). In contrast to analysis of the CMap L1000 Query which
identified MOAs which induce similar gene expression, this DMEA pipeline instead
predicted drug MOAs that may be toxic to cells with similar gene expression profiles as
senescent HMEC (i.e., senolytic MOAs). Using this approach, we found that the EGFR
and MEK inhibitor MOAs were significantly positively enriched in the rank-ordered drug
list (Fig. 7C). This suggests that compounds from these MOAs may be senolytic in
HMECs. Among the negatively enriched drug MOAs which may be more toxic to non-

40
senescent, proliferating HMECs, we found MDM inhibitors, bromodomain inhibitors, and
other MOAs.

Next, we experimentally tested the hypothesis that EGFR inhibitors exhibit
selective toxicity against senescent HMECs using the same cells with which we
generated the proteomic signature of senescence (Delfarah et al., 2019). We treated
proliferating and senescent HMECs with DMSO (negative control), the EGFR inhibitor
dacomitinib, the EGFR inhibitor AZD8931, or the known senolytic compound navitoclax
(positive control) at 100 and 500 nM. After 3 days of treatment, dacomitinib, AZD8931,
and navitoclax significantly reduced the viability of senescent but not proliferating
HMECs at both concentrations (Fig. 7D). Additionally, all three drugs reduced the total
number of viable senescent cells. Notably, low concentrations of dacomitinib and
navitoclax (100 nM) were selectively toxic to senescent cells without affecting the
growth rate of non-senescent, proliferating HMECs. In contrast, although 100 nM of
AZD8931 was selectively toxic to senescent HMEC, AZD8931 also reduced the growth
rate of proliferating, non-senescent HMECs. These experimental results support that the
EGFR inhibitors dacomitinib and AZD8931 are novel senolytic compounds in HMECs,
validating a hypothesis generated by DMEA.

Discussion
Here, we introduce drug mechanism enrichment analysis (DMEA), a user-friendly
bioinformatic method to better prioritize drug candidates for repurposing by grouping
drugs based on shared MOA. Similar to how GSEA enhances biological interpretation of

41
transcriptomic data (Subramanian et al., 2005), DMEA improves drug repurposing by
aggregating information across many drugs with a common MOA instead of considering
each drug independently. We have demonstrated the power and sensitivity of DMEA
first with simulated data (Fig. 3; Appendix A: Fig. S4) and then with real examples
including gene expression connectivity scores (Fig. 4), cell viability connectivity scores
(Fig. 5), and weighted gene voting molecular classification scores (Figs. 6, 7). In all
cases, DMEA ranked the true drug MOA sensitivity or similarity higher than the original
ranking of the single-drug agent and performed similarly to or better than existing tools
for evaluating drug MOA (Appendix A: Fig. S1). In addition, DMEA improves upon
existing tools for analyzing enriched MOA in drug lists in terms of flexibility, statistical
rigor, and visual outputs (Fig. 1). This demonstrates that DMEA helps better prioritize
drug treatments by improving the on-target identification of candidate drugs.

Importantly, our results demonstrate the ability of DMEA to analyze a variety of input
rank-ordered drug lists from different drug repurposing algorithms to identify enrichment
of diverse MOAs (e.g., kinase inhibitors, proteasome inhibitors, metabolic pathway
inhibitors). In these validation cases, DMEA not only identified the expected drug MOAs
(e.g., EGFR inhibitor MOA given a signature of EGFR inhibitor resistance) but also
MOAs which may exhibit toxicity against tumor cells resistant to the input signature of
interest. One interesting example is that DMEA identified HMGCR inhibitors as
potentially toxic to cancer cells with intrinsic resistance to EGFR inhibitors (Fig. 6A,B;
Appendix A: Fig. S5). Indeed, this finding is supported by published work demonstrating
that HMGCR inhibitors can overcome resistance to EGFR inhibitors in NSCLC cells by

42
inhibiting AKT (Hwang et al., 2014; J. Chen et al., 2013). In addition, DMEA also
identified that EGFR inhibitor-resistant cells may be sensitive to MDM inhibitors
(Fig. 6A,B; Appendix A: Fig. S5), a finding that is supported by published work showing
that MDM2 mediates resistance to EGFR inhibitors in mouse models of NSCLC (Shen
et al., 2020). Furthermore, our analysis suggested that melanomas sensitive to BRAF
inhibitors may also be sensitive to MEK inhibitors (Fig. 6C), an observation that is
supported by clinical trials showing that combination treatment with BRAF and MEK
inhibitors is more effective than inhibition of BRAF alone in BRAF-mutant melanoma
patients (Long et al., 2014). Finally, for melanomas with acquired resistance to RAF
inhibitors, DMEA identified CHK inhibitors and SYK inhibitors as potentially beneficial. In
fact, both CHK1 and SYK kinases have been identified as drug targets for melanomas
resistant to RAF inhibitors (Hwang et al., 2018; Abecunas et al., 2023). Collectively,
these results support that DMEA can even identify drug mechanisms beneficial for
combination treatments and drug-resistant cancers.

To demonstrate the power of DMEA for biological discovery, we analyzed our
recently published proteomic signature of replicative senescence in primary HMECs
(Delfarah et al., 2021) (Fig. 7). To illustrate the difference between CMap and DMEA’s
interpretation of an input gene signature, we used: 1) the CMap L1000 Query followed
by DMEA to identify similar (e.g., senescence-inducing) drug MOAs and 2) DMEA with
WGV molecular classification scores to identify selectively toxic (e.g., senolytic) drug
MOAs. Both senescence-inducing and senolytic compounds have great therapeutic
promise in aging (Xu et al., 2018; Baar et al., 2017; Baker et al., 2011, 2016; Cai et al.,

43
2020), cancer (Dörr et al., 2013; Guerrero et al., 2019), and other diseases (Hickson et
al., 2019). Among the potential senescence-inducing drug MOAs we identified were
proteasome, HDAC, HMGCR, and MDM inhibitors (Fig. 7B). Indeed, experimental
evidence has shown that proteasome inhibitors induce senescence in primary
fibroblasts (Chondrogianni and Gonos, 2004; Torres et al., 2006) and that HDAC
inhibitors can induce senescence in cancer cells (Lorenz et al., 2011; Venkatesh et al.,
2015; Vargas et al., 2014). For potential senolytic MOAs, we identified EGFR and MEK
inhibitors (Fig. 7C) and subsequently experimentally validated the senolytic activity of
the EGFR inhibitors dacomitinib and AZD8931 in a drug-induced model of HMEC
senescence (Fig. 7D). To our knowledge, this is the first demonstration that EGFR
inhibitors can exhibit senolytic activity. Although we did not test whether MEK inhibitors
would exhibit senolytic activity in primary HMECs, it has been shown in Ras-expressing
cells that MEK inhibitors selectively kill senescent cells (Kochetkova et al., 2017). Taken
together, our results indicate that DMEA is a powerful tool for repurposing drug MOAs
based on selectivity against or similarity to a given molecularly characterized cell state.

Despite the success of our validation examples, we note that DMEA is limited by our
knowledge of drug MOA. For many targeted therapeutics, the putative MOA may be
incorrect (Lin et al., 2019). Nevertheless, DMEA mitigates the risk of false positives by
evaluating groups of drugs which share a MOA rather than relying on results from
individual drugs alone. Thus, even if some drugs are misannotated, DMEA may still
correctly identify enriched MOAs by aggregating information across multiple drugs,
rather than considering each drug independently. However, improved annotation of drug

44
MOAs and targets, potentially through newer approaches including metabolomics
(Holbrook ‐Smith et al., 2022) and proteomics (Savitski et al., 2014), will improve the
power of DMEA.

In summary, drug mechanism enrichment analysis (DMEA) improves prioritization of
drugs for repurposing by grouping drugs that share mechanisms of action (MOAs).
DMEA can thus be used to further process rank-ordered lists of drugs from drug
repurposing algorithms to sharpen on-target signal. To provide an easily accessible tool
for drug repurposing, we also added the option to pair DMEA with WGV molecular
classification as well as public databases of transcriptomic profiles (e.g., L1000, CCLE)
and drug screens (e.g., PRISM). With this feature, DMEA can interpret an input gene
signature to identify drug mechanisms which exhibit selective toxicity towards cell states
(e.g., cancer, senescence). Furthermore, our results support that DMEA has potential to
aid in the discovery of therapeutics for combination treatments or drug-resistant
cancers. DMEA is publicly available to use either as a web application or an R package
at https://belindabgarana.github.io/DMEA.
 

45
Chapter 3: Mechanism of action landscape reveals shared
sensitivities in large-scale drug screens
Background
While targeted therapies have the potential for efficacy with minimal effect on
normal tissues, it is estimated that only seven percent of cancer patients respond to
genome-targeted therapies, representing about half of patients deemed eligible based on
expression of the target gene (Haslam et al., 2021). Furthermore, even if a patient’s
cancer is responsive to the treatment, resistance is a major obstacle to cancer remission
and occurs almost universally. Consequently, a combination of multiple drugs is often
used, though multidrug resistance can still be observed (Szakács et al., 2006).
Additionally, disparities in health care coverage can leave patients without access to new
technologies such as next generation sequencing which can better inform treatment
strategies. Considering this leaves over 90% of cancer patients to face treatment options
with potentially fatal side effects, we need a better way to prioritize combination therapies
and treatments for resistant cells for further evaluation.

 Recent developments in large-scale pharmacological screenings across
molecularly characterized cancer cell lines enable the systematic analysis of drug
sensitivity in cancer. However, there have been limited experimental studies on drug
combinations because they are exponentially more expensive than single-drug studies.
Many computational frameworks have been previously developed in an effort to reduce
costs by identifying potential drug combinations (Feala et al., 2010; Sun et al., 2016; Pang
et al., 2014; Goswami et al., 2015; Yang et al., 2015; Li et al., 2015; Wildenhain et al.,

46
2015; Bansal et al., 2014; Preuer et al., 2018; Wooten and Albert, 2021; Flobak et al.,
2017; Sun et al., 2020; Sidorov et al., 2019; Hongyang Li et al., 2019). However, these
methods have only studied combinations for a limited number of drugs (on the order of
10-100). Additionally, these other methods only evaluate individual drug pairings, rather
than pairings of drug mechanisms of action (MOAs). As we have shown previously
(Garana et al., 2023), evaluating many drugs grouped by their MOA annotations can
better prioritize drug candidates compared to just evaluating individual drugs. In this
study, we introduce the mechanism of action landscape (MOAL) tool which can analyze
published large-scale drug screens to prioritize MOA pairings to further evaluate in
combination. In particular, we analyzed the PRISM drug screen with 1,351 drugs
annotated with mechanisms of action (MOAs) screened in 320 adherent cancer cell lines
and the CTRPv2 spanning 481 drugs screened in 406 adherent cancer cell lines.

Though our approach does not examine the physical or chemical similarities of
drugs to determine their potential value in combination (i.e., synergy), it can leverage
published drug screens to identify: (1) drug MOAs which may be toxic to samples known
to be resistant to another drug MOA and (2) drug MOAs which tend to be toxic to the
same samples. By grouping drugs by shared MOA and statistically evaluating these MOA
enrichments against a large background of drugs, MOAL removes most potential drug
MOA combinations from consideration to greatly reduce the list of candidates to evaluate
in combination (e.g., just 131 MOA pairings or less than four percent pass our significance
threshold when starting with more than 1.8 million drug pairings). While further
experimental evaluation of the value of MOAL for identifying drug combinations is

47
warranted, this approach has potential to fast-track valuable drug combinations out of a
long candidate list.

Our results show consistency in shared MOA sensitivities across different large-
scale drug screens and added value over UMAP data visualization. Nevertheless, our
findings reveal that results can vary depending on tissue type. We also demonstrate that
this method can be used to evaluate the ability of drug sensitivity prediction algorithms to
recover known relationships identified in large-scale drug screens. The results from this
study are available for query online and the MOA Landscape tool is available as an R
package so that researchers can apply this algorithm to future drug screen datasets.
These resources are publicly available at https://belindabgarana.github.io/MOAL.

Methods
Mechanism of action landscape (MOAL)
To identify drug mechanisms of action (MOAs) which may benefit cancer cell lines
sensitive or resistant to another drug MOA, we applied the drug mechanism enrichment
analysis (DMEA) method (Garana et al., 2023) we developed previously to evaluate
relationships between MOA pairs in a drug screen dataset (Fig. 8). First, we ran Pearson
correlations between the sensitivity scores (AUC) for each drug-drug pairing. Next, for
each drug, we ran DMEA by ranking all other drugs based on these Pearson correlation
estimates to calculate the enrichment of drug MOAs. After this step, drugs are labeled as
potentially misannotated if their enrichment for their own MOA is not greater than or equal
to the average positive enrichment (normalized enrichment score, NES >= 1) or a user-

48
defined minimum NES threshold. Optionally, drugs labeled as potentially misannotated
can be removed from downstream analyses as a quality control measure. Though our
figures feature results without drugs failing this quality control step unless otherwise
specified, we have completed our analyses both ways and made all these results
available on our web application and GitHub repository for transparency.

For each MOA evaluated, we ran DMEA again, this time ranking the drugs based
on the NES scores for the drug MOA at hand. Finally, we perform quality control measures
to identify significant drug MOA pairings. These quality control measures include
removing drug MOA pairs from consideration unless: (1) each MOA in the pairing is
significantly positively self-enriched (e.g., EGFR inhibitors should be identified as
potentially beneficial for cancers sensitive to EGFR inhibitors), (2) the MOA pairing is
significant both ways (e.g., EGFR inhibitors should be enriched for HMGCR inhibitors and
vice-versa), and (3) the MOA pairing relationship should be in the same direction both
ways (e.g., EGFR inhibitors should be negatively enriched for HMGCR inhibitors if
HMGCR inhibitors are negatively enriched for EGFR inhibitors). Significance thresholds
are p < 0.05 and q < 0.25 unless changed by the user.


49

Figure 8. The MOA landscape approach to identify shared drug sensitivities
First, drug sensitivity (AUC) is correlated between each pair of drugs in the given dataset
(e.g., the PRISM drug screen of 1,351 drugs annotated with 85 drug mechanisms of
action). Second, based on Pearson correlations with each drug, we calculate enrichment
of each drug mechanism of action (MOA). After this step, we identify drugs which are not
positively enriched for their own MOA which could optionally be removed from
downstream analyses based on a user-defined threshold. Third, we rank drugs based on
enrichment of each drug MOA to calculate enrichment of drug MOAs again to ultimately
obtain MOA vs. MOA enrichments. Finally, quality control is performed to reduce the
number of significant MOA pairings for further consideration.

MOAL: simulations
To verify if the MOAL method can detect true relationships between drug sets, we
first sampled drug sensitivity scores for 60 drugs grouped into 10 sets (each containing 6
drugs) screened against 25 cell lines (representing the sample size of our analysis for
cancers of the central nervous system) such that the sensitivity scores were
representative of the PRISM drug screen AUC values for 1,351 drugs tested in 320
adherent cancer cell lines (a bimodal mixture of two normal distributions where μ1=0.83,
σ1=0.08 and μ2=1.31, σ2=0.08 and 72% of drugs were in the lower (i.e., more toxic)
distribution). Next, to introduce a synthetic association between two drug sets, we
sampled the drug sensitivity scores for these two drug sets from shifted distributions
where the magnitude of the shift for each cell line is determined by a perturbation value

50
ranging from -0.2 to 0.2. For example, for a perturbation value of -0.2, the mean drug
sensitivity of cell line 1 to the set of 6 perturbed drugs was shifted by +0.2, and this shift
value sequentially decreased such that the mean drug sensitivity of cell line 25 to the
perturbed drug set was shifted by -0.2. Lastly, we processed this synthetic drug sensitivity
data set using the MOAL method and stored the enrichment results between the two
perturbed drug sets. This process was repeated 25 times to test reproducibility.

MOAL: CEVIChE
The CEVIChE’s predictions of sensitivity to drugs in the CTRPv2 drug screen were
downloaded from https://saezlab.shinyapps.io/ceviche/ (Szalai et al., 2019). To reduce
the size of the dataset such that there is only one sensitivity score per drug for each cell
line, we only considered viability predictions for six-hour, 10 uM drug treatments because
this combination of treatment time and concentration covered the greatest number of cell
lines. This drug matrix (containing 455 drugs screened in 67 cancer cell lines) was used
as an input to MOAL to evaluate how well the algorithm recovers shared drug MOA
sensitivities from the CTRPv2 drug screen.

MOAL: weighted gene voting (WGV)
To demonstrate evaluation of a drug sensitivity prediction algorithm which can
easily be applied to the PRISM drug screen, we used the simple molecular classification
method weighted gene voting (WGV). For each drug, we identified the top 500
differentially expressed genes (ranked by the absolute value of the log2(fold-change) and
filtered for p & q < 0.05) between the top and bottom 10% of adherent cancer cell lines

51
ranked by the sensitivity metric AUC (Area Under the Curve) using the eBayes function
from the R package limma (Liu et al., 2015). Next, we performed WGV for each drug by
using the log2(fold-change) of the top 500 differentially expressed genes as weights. In
WGV, each cell line is ranked by the sum across all 500 genes of the product of the gene
weight and the RNAseq normalized expression value, where higher scores suggest the
cell line has more similar expression to the sensitive versus resistant cell lines. To do this,
we accessed CCLE RNAseq data (v19Q4) (Ghandi et al., 2019) for 320 adherent cell
lines which were also studied in the PRISM drug screen. After repeating this procedure
for 374 drugs which had differential gene expression as described above, we used the
resulting drug matrix as an input to MOAL to evaluate how well this WGV-based algorithm
recovers MOA-level drug sensitivity relationships in the PRISM drug screen.

UMAP analysis
UMAP data for the PRISM drug sensitivity scores (AUC) was downloaded from the
PRISM website (https://www.theprismlab.org) and is available to view as a figure in their
journal article (Corsello et al., 2020). The average x- and y-positions were calculated for
each provided MOA annotation and used to calculate the distance between each MOA
pair. The Pearson correlation was evaluated between this average UMAP distance and
our average MOA landscape NES score across all 78 non-self, quality controlled MOA
pairings which were given as MOA annotations for the UMAP.


52
Results
MOA landscape detects coordinated enrichment of drug sets
To test if our method can detect real coordinated shifts in sensitivity between
pairs of drug sets, we first applied our method to synthetic data sets. These synthetic
data sets were designed to replicate real drug sensitivity distributions observed in the
PRISM drug screen (see Methods). By shifting the sensitivity scores (i.e., AUC values)
of drugs belonging to artificial sets, we tested if our method detects real relationships
between two drug sets.

Importantly, no false discovery was observed in our simulations (i.e., 0% of
simulations pass quality control if the drug sets’ sensitivity scores were not shifted; Fig.
9B). Also, 100% of simulations pass quality control if both perturbed drug sets were
shifted by at least two standard deviations. Additionally, the sign of the normalized
enrichment scores between the drug sets was positive if the shifts were in the same
direction (i.e., drugs were more toxic in both sets or vice-versa) or negative if the shifts
were in opposite direction (i.e., drugs were more toxic in one set and less toxic in
another set) as expected (Fig. 9A). These results show that the MOA landscape method
can detect real relationships between drug sets.


53

Figure 9. MOA landscape detects coordinated enrichment of drug sets  
Synthetic data was generated to be representative of the PRISM drug screen’s sensitivity
scores (AUC). Two drug sets were shifted by varying amounts (x- and y-axes). For each
combination of drug set shifts, 25 simulations were performed. A) Tile plot showing the
average normalized enrichment score (NES) between the two drug sets. B) Tile plot
showing the percent of simulations passing quality control (QC) thresholds.

MOA landscape removes more than 90% of drug MOA pairs from further consideration
To evaluate our method, we first determined if we could identify self-enrichments
across 85 drug mechanisms of action (MOAs) in the PRISM dataset (e.g., that cell lines
sensitive to EGFR inhibitors would be identified as sensitive to EGFR inhibitors; Fig. 10A).
Since cell behavior can depend on environmental factors, we also evaluated MOA
relationships in adherent cancer cell lines stratified by tissue type (Fig. 13; Fig. S6). With
our default method of not removing potentially misannotated drugs to avoid bias, only
about 58% of drug MOAs were significantly positively self-enriched as expected when
considering all 320 adherent cancer cell lines and this percentage drops with the number
of cell lines considered. If the user decides to remove potentially misannotated drugs as
a form of quality control, most MOAs are self-enriched as expected even when only a few

54
cell lines are considered (more than 89% of MOAs for all 320 adherent cancer cell lines,
or more than 76% of MOAs regardless of the number of adherent cancer cell lines
considered in various tissue-based subsets).  

However, regardless of if potentially misannotated drugs are removed from
downstream analyses or not, less than nine percent of MOA pairings pass the final quality
control measures for further consideration as potential combination treatments (Fig. 10B).
If potentially misannotated drugs are not removed from downstream analyses, more
MOAs are evaluated for MOA-MOA enrichments and less than four percent of MOA
pairings pass the final quality control measures. Additionally, the normalized enrichment
scores (NES) for MOA pairings are strongly correlated between our analyses with or
without inclusion of potentially misannotated drugs (r = 0.95, p < 1E-323; Fig. 10C); this
correlation is even stronger if only considering MOA pairings which pass quality control
with at least one of these methods (r = 0.99, p = 2.1E-232), though more MOA pairings
pass quality control when potentially misannotated drugs are removed (Fig. 10D).


55

Figure 10. MOA landscape removes more than 90% of drug MOA pairs from further
consideration  
A) Drug mechanisms of action (MOAs) were considered as passing quality control if they
were significantly positively self-enriched (i.e., NES > 0, p < 0.05, and FDR < 0.25). B)
Drug MOA pairs were considered as passing quality control if: (1) both MOAs pass quality
control individually, (2) both MOAs were significantly enriched for each other, and (3) both
MOAs were enriched for each other in the same direction (see Methods). C) Normalized
enrichment scores (NES) for MOA pairings are correlated whether potentially
misannotated drugs are removed from downstream analyses or not. D) Venn diagram
showing that most of the same MOA pairings pass quality control whether potentially
misannotated drugs are removed from downstream analyses or not.  

MOA enrichment is correlated across large-scale drug screens
To corroborate the shared sensitivities identified from our analysis of the PRISM
drug screen, we also analyzed drug sensitivity scores from the CTRPv2 drug screen.
However, the CTRPv2 analysis only evaluated 14 drug MOAs compared to PRISM’s 85
MOAs because it only considers 481 drugs compared to PRISM’s 1,351 drugs. Despite
the CTRPv2 dataset being smaller than PRISM, 100% of the drug MOAs in the CTRPv2

56
dataset were self-enriched as expected. Since the PRISM dataset covers more MOAs
(85), it also identified more shared drug sensitivities (i.e., MOA-MOA enrichments) but the
overlapping sensitivities were in the same direction and strongly correlated (Fig. 11A; r =
0.76, p = 4.2E-19); this correlation is slightly stronger if drugs identified as potentially
misannotated are removed before the final enrichment analyses (Fig. 11B; r = 0.78, p =
8.3E-12). All the significant, overlapping MOA pairings were for sensitive cancer cell lines
(meaning samples sensitive to one MOA may also be sensitive to another). We suspect
this may be due to the cell lines in the PRISM drug screen demonstrating sensitivity to
most of the drugs studied and the relatively few MOAs represented in the CTRPv2
dataset.  


57

Figure 11. MOA enrichment is correlated across large-scale drug screens  
Average normalized enrichment scores (NES) are correlated across MOAL results for the
PRISM and CTRPv2 drug screens (r = 0.76 and p = 4.2E-19 with our default analysis
method shown in panel A, or r = 0.78 and p = 8.3E-12 if potentially misannotated drugs
are removed from downstream analyses as shown in panel B). Orange datapoints
represent MOA pairings which pass quality control thresholds in both datasets. Non-self
MOA pairings with less than 15% overlap are labeled.

Now that we verified our results are consistent across drug screens, we studied
MOA pairings from the PRISM drug screen more closely since they may be worth
prioritizing for evaluation in combination and the size of the PRISM dataset enables
researchers to study about six-fold more drug MOAs than the CTRPv2 dataset.
Specifically, we scored the strength of each MOA pairing passing our quality control
measures by their average normalized enrichment scores (NES). The pairings with the
greatest positive score represent the top potential drug combinations for sensitive cancer

58
cell lines. Conversely, the pairings with the greatest negative score represent the top
potential drug combinations for resistant cancer cell lines.  

It is worth noting that the Broad Institute also published a UMAP visualization which
they called a “drug response landscape,” but our landscape approach differs because our
algorithm groups drugs by common MOA and performs correlations and enrichment
analyses to evaluate the statistical significance of these MOA relationships instead of
solely considering results for individual drugs. Furthermore, we find that most MOA
pairings are determined to be insignificant in our MOA landscape even if they are close
together (i.e., show similarity) in the UMAP visualization and some of our strongest MOA
pairings have an above-average UMAP distance (e.g., bromodomain and HDAC
inhibitors; Fig. 12). With our statistical analysis and quality control methods, the MOA
landscape can eliminate insignificant MOA relationships which represent more than 90%
of MOA pairings evaluated and identify significant MOA pairings which may not be
obvious from distances on a UMAP.


59

Figure 12. MOA landscape identifies drug pairings not obvious from UMAP
visualization  

60
A) Whether or not potentially misannotated drugs are removed, the MOA landscape
(MOAL) statistically narrows down the number of candidate drug pairings to identify both
drug pairings which might be obvious from a UMAP visualization (e.g., MEK & RAF
inhibitors, PLK & tubulin polymerization inhibitors, and EGFR & HMGCR inhibitors,
labeled in the plot) and others which may not be obvious (e.g., bromodomain & HDAC
inhibitors may be toxic towards the same cell lines despite being relatively far apart on
the UMAP visualization and CHK & MDM inhibitors or KIT inhibitors & vitamin D receptor
agonists might be toxic towards different cell lines despite being relatively close on the
UMAP visualization, labeled in the plot). QC: quality control. B) Mountain plots illustrating
that adherent cancer cell lines sensitive to tubulin polymerization, RAF, or HDAC
inhibitors may also be sensitive to PLK, MEK, or bromodomain inhibitors, respectively,
and vice-versa. C) Mountain plots illustrating that adherent cancer cell lines resistant to
EGFR, CHK, or KIT inhibitors may be sensitive to HMGCR inhibitors, MDM inhibitors, or
vitamin D receptor agonists, respectively, and vice-versa.

Additionally, our comprehensive analysis serves as a resource to help researchers
prioritize drug mechanisms to evaluate for toxicity towards cells sensitive or resistant to
another drug mechanism. Using our web application, researchers can filter based on
dataset (i.e., PRISM, CTRPv2, or DEPMAP); tissue type; drug MOA; drug name; and
choose whether to examine results which did or did not include potentially misannotated
drugs. With this query feature, researchers can more closely examine drug-MOA or MOA-
MOA relationships across different cancer cell line subsets and datasets before further
evaluating drug combinations experimentally. Next, we highlight differences in MOA
relationships across different tissue types to illustrate the value of querying these different
results on our web application.

MOA relationships can depend on tissue type
It is well known that the tumor microenvironment can influence cell behavior.
Thus, we investigated if MOA pairs identified across all adherent cancer cell lines are
relevant for individual tissue types (Appendix B: Table S1; Fig. S6). As expected, MOA

61
pairings significant across all adherent cancer cell lines may not necessarily be
significant in individual tissue types. In fact, none of the pairwise drug MOA
relationships we studied are significant across all tissue types in our analysis (Appendix
B: Fig. S6).

When considering the top MOA pairings across tissue types represented by at
least 25 adherent cancer cell lines, it is still apparent that results vary (Fig. 13). Though
several positive enrichments (i.e., co-sensitivities) may be shared, negative enrichments
(i.e., opposite sensitivities) spanning the top 3 MOA pairs are not shared across cancers
of the lung (N = 63), skin (N = 29), or central nervous system (N = 25). As expected, we
found that MEK inhibitors may be toxic towards skin cancers sensitive to RAF inhibitors
and vice-versa. We also found that skin cancers resistant to MEK inhibitors may be
sensitive to PARP inhibitors and vice-versa. Cancers of the central nervous system
resistant to CHK inhibitors may be sensitive to HMGCR inhibitors and vice-vera.
Additionally, lung cancers resistant to EGFR inhibitors may be sensitive to HMGCR
inhibitors and vice-versa. To examine results across more tissue types, please refer to
Appendix B: Fig. S6 or our web application at https://belindabgarana.github.io/MOAL.


62

Figure 13. MOA relationships can depend on tissue type  

63
A) The top 3 positively and top 3 negatively enriched drug MOA pairings for each tissue-
based subset of adherent cancer cell lines reveal that MOA relationships can depend on
tissue type. B) Adherent lung cancer cell lines resistant to EGFR inhibitors may be
sensitive to HMGCR inhibitors and vice-versa. C) Adherent cancer cell lines of the central
nervous system resistant to CHK inhibitors may be sensitive to HMGCR inhibitors and
vice-versa. D) Adherent skin cancer cell lines resistant to MEK inhibitors may be sensitive
to PARP inhibitors and vice-versa.

MOA landscape can evaluate network-level performance of drug sensitivity prediction
algorithms
The MOA landscape not only has potential as a new framework to analyze large-
scale drug screens, but also to evaluate how well drug sensitivity prediction algorithms
can recover MOA relationships observed in these large-scale drug screens. To
demonstrate this capability, we analyzed a well-respected machine learning-based drug
sensitivity prediction algorithm for the CTRPv2 drug screen, known as CEVIChE, in
addition to a simple molecular classification algorithm which can easily be applied to
larger datasets, known as weighted gene voting (WGV). The MOA landscape results
generated from CEVIChE’s predicted drug sensitivity scores are correlated with that of
the CTRPv2 drug screen on which it was based (Fig. 14A; r = 0.65, p = 2.1E-5). The
results from the simpler algorithm WGV were not as strongly correlated with the PRISM
drug screen on which it was based (Fig. 14B; r = 0.5, p = 5.6E-20). However, the WGV
algorithm was able to recover more than twice as many significant MOA relationships
compared to the CEVIChE algorithm, though it did evaluate more than four times as many
MOA pairs because of the difference in the number of drugs between the PRISM and
CTRPv2 drug screens (Fig. 14C-E). Still, a benefit of the CEVIChE algorithm is that it did
not falsely identify any MOAs as significantly enriched whereas the WGV algorithm
portrayed six MOA pairings as significant which did not pass quality control in our analysis

64
the original PRISM drug screen. Nevertheless, neither of these algorithms were able to
recover most of the significant MOA relationships found in these large-scale drug screens.
This analysis illustrates how the MOA landscape can be used to evaluate the global
performance of drug sensitivity prediction algorithms based on relationships between
MOAs instead of just considering individual drugs.


Figure 14. MOA landscape can evaluate network-level performance of drug
sensitivity prediction algorithms  
A & C) CEVIChE predictions of drug sensitivity scores recover most MOA relationships
from the CTRPv2 drug screen, and our MOA landscape analysis of these datasets shows
correlation between normalized enrichment scores (NES) for MOA pairings (r = 0.65, p =
2.1E-5). B & D) Weighted gene voting (WGV) predictions of drug sensitivity scores

65
recover most MOA relationships from the PRISM drug screen and our MOA landscape
analysis of these datasets shows correlation between normalized enrichment scores
(NES) for MOA pairings (r = 0.5, p = 5.6E-20). E) Success, false positive, and false
negative rates for drug sensitivity prediction algorithms CEVIChE and WGV relative to the
MOA pairs evaluated by both the algorithm and the original drug screen, where: (1) the
success rate is the percentage of MOA pairings passing quality control with the prediction
algorithm which also passed quality control in the original drug screen, (2) the false
positive rate is the percentage of MOA pairings passing quality control with the prediction
algorithm which did not pass quality control with the original drug screen, and (3) the false
negative rate is the percentage of MOA pairings passing quality control with the original
drug screen which did not pass quality control with the prediction algorithm.

Discussion
Considering resistance is a major obstacle to successful cancer therapy and
observed almost universally in single-agent treatments, it is crucial to identify drug
combinations which may increase efficacy or target resistant cells. Here, we showed that
our method called the MOA landscape (MOAL) can analyze large-scale drug screens to
prioritize drug MOA pairings for further evaluation as combination treatments or
treatments for resistant cancers. To evaluate our method, we confirmed that it: (1) detects
true coordinated shifts between drug sets based on simulations (Fig. 9), (2) correctly
identifies self-enrichments and removes more than 90% of drug MOA pairings from
further consideration (Fig. 10), and (3) produces results which are consistent between the
PRISM and CTRPv2 drug screens (Fig. 11). Furthermore, we recovered known drug
MOA combinations in addition to previously unreported drug MOA pairings (Figs. 12-13;
Appendix B: Fig. S6). We also demonstrate that it is important to consider that
relationships between drug MOAs may vary based on factors such as tissue type (Fig.
13; Appendix B: Fig. S6).


66
Though we cannot determine if these MOA pairings may be synergistic
combinations without further experimental testing, it is encouraging that some MOA
pairs identified are already established synergistic combinations used in clinical
practice. In particular, MEK and RAF inhibitors are a well-known combination used for
melanoma treatment (Long et al., 2014) and we found that MEK and RAF inhibitors may
be toxic towards the same skin cancers (Fig. 13). Based on our analysis, MEK and RAF
inhibitors may also be beneficial in combination in breast and kidney cancers, which is
supported by a study demonstrating their synergy in models of triple-negative breast
cancer (Nagaria et al., 2017). Furthermore, skin cancers resistant to MEK inhibition may
be sensitive to PARP inhibition and vice-versa (Fig. 13). Indeed, this pairing has been
found to be a synergistic combination in RAS-mutant cancers (Yang et al., 2021; Sun et
al., 2017). Another standard treatment is EGFR inhibitors for EGFR-mutant lung
cancers. As shown in Figure 13, we also recovered a previously discovered
combination where HMGCR inhibitors may be toxic towards lung cancers which are
resistant to EGFR inhibition (Hwang et al., 2014; J. Chen et al., 2013).

Furthermore, our analysis also uncovered potential drug pairings which were not
obvious from a UMAP visualization of the same dataset (Fig. 12). For example,
bromodomain inhibitors were identified as potentially toxic to adherent cancer cell lines
sensitive to HDAC inhibition and vice-versa even though they were relatively far apart
on the PRISM UMAP visualization. Indeed, bromodomain and extra-terminal motif
(BET) inhibitors have been shown to be synergistic in combination with HDAC inhibitors
in many cancers, including AML (Alcitepe et al., 2022; Shao et al., 2017; Fiskus et al.,

67
2014), breast cancer (Borbely et al., 2015; B. Zhang et al., 2020; Huang et al., 2023),
chondrosarcoma (Huan et al., 2020), glioblastoma (Meng et al., 2018; Gusyatiner et al.,
2021), mantle cell lymphoma (Sun et al., 2015), medulloblastoma (Kling et al., 2022),
melanoma (Heinemann et al., 2015), multiple myeloma (Carew et al., 2019),
neuroblastoma (Shahbazi et al., 2016), leukemia (Schäker-Hübner et al., 2021),
lymphoma (Bhadury et al., 2014), B cell lymphoma (Gaudio et al., 2016; Spriano et al.,
2022; Cortiguera et al., 2019), cutaneous T-cell lymphoma (Zhao et al., 2019, 2022),
osteosarcoma (Yu et al., 2022), pancreatic cancer (X. Zhang et al., 2020; Mazur et al.,
2015; He et al., 2020), non-small cell lung cancer (Evanno et al., 2017), small-cell lung
cancer (Liu et al., 2018), rhabdomyosarcoma (Enßle et al., 2018), testicular germ cell
tumors (Hamazaki et al., 1989), urothelial carcinoma (Hölscher et al., 2018), and
Waldenström macroglobulinemia (Matissek et al., 2021). On the other hand, checkpoint
(CHK) inhibitors were identified as potentially toxic to adherent cancer cell lines
resistant to MDM inhibition and vice-versa even though they were relatively close
together on the PRISM UMAP visualization. It has been found that toxicity of MDM
inhibition was improved in renal cancer cell lines when combined with DNA damage and
G2/M checkpoint inhibition (Tsao and Corn, 2010). Another study showed that breast
cancer cell lines resistant to checkpoint inhibition had higher levels of MDM2 (Nieves-
Neira and Pommier, 1999). These examples demonstrate the utility of the MOA
landscape algorithm to identify significant drug MOA pairings which may not be
apparent from a simple UMAP visualization.


68
Despite this algorithm’s ability to recover known drug combinations in tissue-
specific contexts, it has several limitations. Firstly, it is important to distinguish that our
analysis only identifies shared sensitivities (i.e., drug MOAs which tend to be toxic in the
same cell lines if positively enriched, or different cell lines if negatively enriched). This
means that the MOA pairings we have identified may not necessarily be synergistic
combinations. However, across all adherent cancer cell lines, this analysis reduced the
number of MOA pairings for evaluation by over 90%. This addresses a limiting factor in
the evaluation of drug combinations: excessive cost due to the exponential number of
combinations to test given an increasing number of drugs to consider.

Another limitation is that this method requires a dataset with a minimum of two
drug MOAs each represented by at least six drugs studied in at least three biological
samples. Since recent technological advances have made automated drug screens
feasible (Mitchell et al., 2023; Zecha et al., 2023; Garana and Graham, 2022; Holbrook ‐
Smith et al., 2022; Corsello et al., 2020), we foresee this algorithm becoming applicable
to more datasets in the future. Due to limited data availability, our tissue-specific analyses
have relatively few samples so they may not necessarily be representative of all cancers
derived from these tissues and we could not control for other potentially relevant factors
such as sex, media type, lineage subtype (e.g., NSCLC versus SCLC), molecular
subtype, or primary tumor versus metastatic status. Also due to data availability, we were
only able to study drug sensitivity in adherent cancer cell lines, which may not be
representative of drug sensitivity in vivo or in other diseases.


69
Overall, our MOA landscape (MOAL) method has identified several drug
combinations using the PRISM drug screen dataset. This framework can be applied to
any drug screen dataset with mechanism of action (MOA) annotations. This work
represents a significant step forward in the application of computational algorithms to
identify potential drug combinations for sensitive or resistant cancers because it evaluates
many drugs grouped by MOA rather than just relatively few individual drugs. Prioritizing
potential drug combinations with this kind of analysis will be important for designing cost-
effective drug combination experiments and ultimately for addressing the issue of
resistance in cancer and other diseases. The MOA landscape can be accessed as a web
app or R package at https://belindabgarana.github.io/MOAL.

 

70
References
Abecunas,C. et al. (2023) Loss of NF1 in Melanoma Confers Sensitivity to SYK Kinase
Inhibition. Cancer Res, 83, 316–331.
Alcitepe,İ. et al. (2022) HDAC inhibitor Vorinostat and BET inhibitor Plx51107 epigenetic
agents’ combined treatments exert a therapeutic approach upon acute myeloid
leukemia cell model. Med Oncol, 39, 257.
Baar,M.P. et al. (2017) Targeted Apoptosis of Senescent Cells Restores Tissue
Homeostasis in Response to Chemotoxicity and Aging. Cell, 169, 132-147.e16.
Baker,D.J. et al. (2011) Clearance of p16Ink4a-positive senescent cells delays ageing-
associated disorders. Nature, 479, 232–236.
Baker,D.J. et al. (2016) Naturally occurring p16(Ink4a)-positive cells shorten healthy
lifespan. Nature, 530, 184–189.
Balko,J.M. et al. (2006) Gene expression patterns that predict sensitivity to epidermal
growth factor receptor tyrosine kinase inhibitors in lung cancer cell lines and
human lung tumors. BMC Genomics, 7, 289.
Bansal,M. et al. (2014) A community computational challenge to predict the activity of
pairs of compounds. Nat Biotechnol, 32, 1213–1222.
Behan,F.M. et al. (2019) Prioritization of cancer therapeutic targets using CRISPR-Cas9
screens. Nature, 568, 511–516.
Bhadury,J. et al. (2014) BET and HDAC inhibitors induce similar genes and biological
effects and synergize to kill in Myc-induced murine lymphoma. Proc Natl Acad
Sci U S A, 111, E2721-2730.
Borbely,G. et al. (2015) Induction of USP17 by combining BET and HDAC inhibitors in
breast cancer cells. Oncotarget, 6, 33623–33635.
Brachmann,S.M. et al. (2012a) Characterization of the mechanism of action of the pan
class I PI3K inhibitor NVP-BKM120 across a broad range of concentrations. Mol
Cancer Ther, 11, 1747–1757.
Brachmann,S.M. et al. (2012b) Characterization of the mechanism of action of the pan
class I PI3K inhibitor NVP-BKM120 across a broad range of concentrations. Mol
Cancer Ther, 11, 1747–1757.
Brum,A.M. et al. (2015) Connectivity Map-based discovery of parbendazole reveals
targetable human osteogenic pathway. Proc Natl Acad Sci U S A, 112, 12711–
12716.

71
Cai,Y. et al. (2020) Elimination of senescent cells by β-galactosidase-targeted prodrug
attenuates inflammation and restores physical function in aged mice. Cell Res,
30, 574–589.
Carew,J.S. et al. (2019) Rational cotargeting of HDAC6 and BET proteins yields
synergistic antimyeloma activity. Blood Adv, 3, 1318–1329.
Chan,J. et al. (2019) Breaking the paradigm: Dr Insight empowers signature-free,
enhanced drug repurposing. Bioinformatics, 35, 2818–2826.
Chen,E.Y. et al. (2013) Enrichr: interactive and collaborative HTML5 gene list
enrichment analysis tool. BMC Bioinformatics, 14, 128.
Chen,J. et al. (2013) Atorvastatin overcomes gefitinib resistance in KRAS mutant
human non-small cell lung carcinoma cells. Cell Death Dis, 4, e814.
Chen,R. et al. (2022) CPDR: An R Package of Recommending Personalized Drugs for
Cancer Patients by Reversing the Individual’s Disease-Related Signature. Front
Pharmacol, 13, 904909.
Chondrogianni,N. and Gonos,E.S. (2004) Proteasome inhibition induces a senescence-
like phenotype in primary human fibroblasts cultures. Biogerontology, 5, 55–61.
Coldren,C.D. et al. (2006) Baseline Gene Expression Predicts Sensitivity to Gefitinib in
Non–Small Cell Lung Cancer Cell Lines. Mol Cancer Res, 4, 521–528.
Corsello,S.M. et al. (2020) Discovering the anticancer potential of non-oncology drugs
by systematic viability profiling. Nat Cancer, 1, 235–248.
Cortiguera,M.G. et al. (2019) Suppression of BCL6 function by HDAC inhibitor mediated
acetylation and chromatin modification enhances BET inhibitor effects in B-cell
lymphoma cells. Sci Rep, 9, 16495.
Delfarah,A. et al. (2021) Identification of a Proteomic Signature of Senescence in
Primary Human Mammary Epithelial Cells. J. Proteome Res., 20, 5169–5179.
Delfarah,A. et al. (2019) Inhibition of nucleotide synthesis promotes replicative
senescence of human mammary epithelial cells. J. Biol. Chem., 294, 10564–
10578.
Domingo-Fernández,D. et al. (2022) Causal reasoning over knowledge graphs
leveraging drug-perturbed and disease-specific transcriptomic signatures for drug
discovery. PLoS Comput Biol, 18, e1009909.
Dörr,J.R. et al. (2013) Synthetic lethal metabolic targeting of cellular senescence in
cancer therapy. Nature, 501, 421–425.

72
Emon,M.A. et al. (2020) PS4DR: a multimodal workflow for identification and
prioritization of drugs based on pathway signatures. BMC Bioinformatics, 21,
231.
Enßle,J.C. et al. (2018) Co-targeting of BET proteins and HDACs as a novel approach
to trigger apoptosis in rhabdomyosarcoma cells. Cancer Lett, 428, 160–172.
Evanno,E. et al. (2017) Tri-methylation of H3K79 is decreased in TGF-β1-induced
epithelial-to-mesenchymal transition in lung cancer. Clin Epigenetics, 9, 80.
Fang,M. et al. (2021) Drug perturbation gene set enrichment analysis (dpGSEA): a new
transcriptomic drug screening approach. BMC Bioinformatics, 22, 22.
Feala,J.D. et al. (2010) Systems approaches and algorithms for discovery of
combinatorial therapies. WIREs Syst Biol Med, 2, 181–193.
Fiskus,W. et al. (2014) Highly active combination of BRD4 antagonist and histone
deacetylase inhibitor against human acute myelogenous leukemia cells. Mol
Cancer Ther, 13, 1142–1154.
Flobak,Å. et al. (2017) CImbinator: a web-based tool for drug synergy analysis in small-
and large-scale datasets. Bioinformatics, 33, 2410–2412.
Frejno,M. et al. (2020) Proteome activity landscapes of tumor cell lines determine drug
responses. Nat Commun, 11, 3639.
Gao,S. et al. (2021) Modeling drug mechanism of action with large scale gene-
expression profiles using GPAR, an artificial intelligence platform. BMC
Bioinformatics, 22, 17.
Garana,B.B. et al. (2023) Drug mechanism enrichment analysis improves prioritization
of therapeutics for repurposing. BMC Bioinformatics, 24, 215.
Garana,B.B. and Graham,N.A. (2022) Metabolomics paves the way for improved drug
target identification. Molecular Systems Biology, 18.
Gaudio,E. et al. (2016) Bromodomain inhibitor OTX015 (MK-8628) combined with
targeted agents shows strong in vivo antitumor activity in lymphoma. Oncotarget,
7, 58142–58147.
Ghandi,M. et al. (2019) Next-generation characterization of the Cancer Cell Line
Encyclopedia. Nature, 569, 503–508.
Golub,T.R. et al. (1999) Molecular Classification of Cancer: Class Discovery and Class
Prediction by Gene Expression Monitoring. Science, 286, 531–537.
Gonçalves,E. et al. (2020) Drug mechanism-of-action discovery through the integration
of pharmacological and CRISPR screens. Mol Syst Biol, 16, e9405.

73
Goswami,C.P. et al. (2015) A New Drug Combinatory Effect Prediction Algorithm on the
Cancer Cell Based on Gene Expression and Dose–Response Curve. CPT
Pharmacometrics Syst. Pharmacol., 4, 80–90.
Greger,J.G. et al. (2012) Combinations of BRAF, MEK, and PI3K/mTOR inhibitors
overcome acquired resistance to the BRAF inhibitor GSK2118436 dabrafenib,
mediated by NRAS or MEK mutations. Mol Cancer Ther, 11, 909–920.
Guerrero,A. et al. (2019) Cardiac glycosides are broad-spectrum senolytics. Nat Metab,
1, 1074–1088.
Gusyatiner,O. et al. (2021) BET inhibitors repress expression of interferon-stimulated
genes and synergize with HDAC inhibitors in glioblastoma. Neuro Oncol, 23,
1680–1692.
Hamazaki,S. et al. (1989) Thiobarbituric acid-reactive substance formation of rat kidney
brush border membrane vesicles induced by ferric nitrilotriacetate. Arch Biochem
Biophys, 274, 348–354.
Han,X. et al. (2021) SubtypeDrug: a software package for prioritization of candidate
cancer subtype-specific drugs. Bioinformatics, 37, 2491–2493.
Haslam,A. et al. (2021) Updated estimates of eligibility for and response to genome-
targeted oncology drugs among US cancer patients, 2006-2020. Ann Oncol.
He,B. et al. (2023) ASGARD is A Single-cell Guided Pipeline to Aid Repurposing of
Drugs. Nat Commun, 14, 993.
He,S. et al. (2020) Potent Dual BET/HDAC Inhibitors for Efficient Treatment of
Pancreatic Cancer. Angew Chem Int Ed Engl, 59, 3028–3032.
Heinemann,A. et al. (2015) Combining BET and HDAC inhibitors synergistically induces
apoptosis of melanoma and suppresses AKT and YAP signaling. Oncotarget, 6,
21507–21521.
Hickson,L.J. et al. (2019) Senolytics decrease senescent cells in humans: Preliminary
report from a clinical trial of Dasatinib plus Quercetin in individuals with diabetic
kidney disease. EBioMedicine, 47, 446–456.
Hoeflich,K.P. et al. (2009) In vivo antitumor activity of MEK and phosphatidylinositol 3-
kinase inhibitors in basal-like breast cancer models. Clin Cancer Res, 15, 4649–
4664.
Holbrook ‐Smith,D. et al. (2022) High ‐throughput metabolomics predicts drug–target
relationships for eukaryotic proteins. Molecular Systems Biology, 18.

74
Hölscher,A.S. et al. (2018) Combined inhibition of BET proteins and class I HDACs
synergistically induces apoptosis in urothelial carcinoma cell lines. Clin
Epigenetics, 10, 1.
Huan,S. et al. (2020) Combination BET Family Protein and HDAC Inhibition
Synergistically Elicits Chondrosarcoma Cell Apoptosis Through RAD51-Related
DNA Damage Repair. Cancer Manag Res, 12, 4429–4439.
Huang,C. et al. (2018) The DrugPattern tool for drug set enrichment analysis and its
prediction for beneficial effects of oxLDL on type 2 diabetes. Journal of Genetics
and Genomics, 45, 389–397.
Huang,Y. et al. (2023) BET-HDAC Dual Inhibitors for Combinational Treatment of
Breast Cancer and Concurrent Candidiasis. J Med Chem, 66, 1239–1253.
Hwang,B.-J. et al. (2018) Chk1 inhibition as a novel therapeutic strategy in melanoma.
Oncotarget, 9, 30450–30464.
Hwang,K.-E. et al. (2014) Effect of simvastatin on the resistance to EGFR tyrosine
kinase inhibitors in a non-small cell lung cancer with the T790M mutation of
EGFR. Exp Cell Res, 323, 288–296.
Jia,Z. et al. (2016) Cogena, a novel tool for co-expressed gene-set enrichment analysis,
applied to drug repositioning and drug mode of action discovery. BMC Genomics,
17, 414.
Kemper,K. et al. (2016) BRAF(V600E) Kinase Domain Duplication Identified in Therapy-
Refractory Melanoma Patient-Derived Xenografts. Cell Rep, 16, 263–277.
Kling,M.J. et al. (2022) A novel dual epigenetic approach targeting BET proteins and
HDACs in Group 3 (MYC-driven) Medulloblastoma. J Exp Clin Cancer Res, 41,
321.
Kochetkova,E.Y. et al. (2017) Targeted elimination of senescent Ras-transformed cells
by suppression of MEK/ERK pathway. Aging (Albany NY), 9, 2352–2375.
Kropiwnicki,E. et al. (2021) Drugmonizome and Drugmonizome-ML: integration and
abstraction of small molecule attributes for drug enrichment analysis and
machine learning. Database (Oxford), 2021, baab017.
Kuleshov,M.V. et al. (2016) Enrichr: a comprehensive gene set enrichment analysis
web server 2016 update. Nucleic Acids Res, 44, W90-97.
Kuleshov,M.V. et al. (2019) modEnrichr: a suite of gene set enrichment analysis tools
for model organisms. Nucleic Acids Res, 47, W183–W190.
Li,Hongyang et al. (2019) TAIJI: approaching experimental replicates-level accuracy for
drug synergy prediction. Bioinformatics, 35, 2338–2339.

75
Li,Haoxin et al. (2019) The landscape of cancer cell line metabolism. Nat Med, 25, 850–
860.
Li,P. et al. (2015) Large-scale exploration and analysis of drug combinations.
Bioinformatics, 31, 2007–2016.
Lin,A. et al. (2019) Off-target toxicity is a common mechanism of action of cancer drugs
undergoing clinical trials. Sci. Transl. Med., 11, eaaw8412.
Liu,R. et al. (2015) Why weight? Modelling sample and observational level variability
improves power in RNA-seq analyses. Nucleic Acids Res, 43, e97–e97.
Liu,Y. et al. (2018) NK Cells Mediate Synergistic Antitumor Effects of Combined
Inhibition of HDAC6 and BET in a SCLC Preclinical Model. Cancer Res, 78,
3709–3717.
Long,G.V. et al. (2014) Combined BRAF and MEK inhibition versus BRAF inhibition
alone in melanoma. N Engl J Med, 371, 1877–1888.
Lorenz,V. et al. (2011) Sodium butyrate induces cellular senescence in neuroblastoma
and prostate cancer cells. Horm Mol Biol Clin Investig, 7, 265–272.
Maejima,T. et al. (2014) Direct evidence for pitavastatin induced chromatin structure
change in the KLF4 gene in endothelial cells. PLoS One, 9, e96005.
Matissek,S.J. et al. (2021) Epigenetic targeting of Waldenström macroglobulinemia cells
with BET inhibitors synergizes with BCL2 or histone deacetylase inhibition.
Epigenomics, 13, 129–144.
Mazur,P.K. et al. (2015) Combined inhibition of BET family proteins and histone
deacetylases as a potential epigenetics-based therapy for pancreatic ductal
adenocarcinoma. Nat Med, 21, 1163–1171.
Meng,W. et al. (2018) Enhanced efficacy of histone deacetylase inhibitor combined with
bromodomain inhibitor in glioblastoma. J Exp Clin Cancer Res, 37, 241.
Meyers,R.M. et al. (2017) Computational correction of copy number effect improves
specificity of CRISPR–Cas9 essentiality screens in cancer cells. Nat Genet, 49,
1779–1784.
Mitchell,D.C. et al. (2023) A proteome-wide atlas of drug mechanism of action. Nat
Biotechnol.
Nagaria,T.S. et al. (2017) Combined targeting of Raf and Mek synergistically inhibits
tumorigenesis in triple negative breast cancer model systems. Oncotarget, 8,
80804–80819.

76
Napolitano,F. et al. (2015) Drug-set enrichment analysis: a novel tool to investigate drug
mode of action. Bioinformatics, btv536.
Napolitano,F. et al. (2018) gene2drug: a computational tool for pathway-based rational
drug repositioning. Bioinformatics, 34, 1498–1505.
Narla,A. et al. (2011) Dexamethasone and lenalidomide have distinct functional effects
on erythropoiesis. Blood, 118, 2296–2304.
Nguyen,L. et al. (2017) Systematic assessment of multi-gene predictors of pan-cancer
cell line sensitivity to drugs exploiting gene expression data. F1000Res, 5, 2927.
Nieves-Neira,W. and Pommier,Y. (1999) Apoptotic response to camptothecin and 7-
hydroxystaurosporine (UCN-01) in the 8 human breast cancer cell lines of the
NCI Anticancer Drug Screen: multifactorial relationships with topoisomerase I,
protein kinase C, Bcl-2, p53, MDM-2 and caspase pathways. Int J Cancer, 82,
396–404.
Nusinow,D.P. et al. (2020) Quantitative Proteomics of the Cancer Cell Line
Encyclopedia. Cell, 180, 387-402.e16.
Pang,K. et al. (2014) Combinatorial therapy discovery using mixed integer linear
programming. Bioinformatics, 30, 1456–1463.
Parca,L. et al. (2019) Modeling cancer drug response through drug-specific informative
genes. Sci Rep, 9, 15222.
Preuer,K. et al. (2018) DeepSynergy: predicting anti-cancer drug synergy with Deep
Learning. Bioinformatics, 34, 1538–1546.
Rampášek,L. et al. (2019) Dr.VAE: improving drug response prediction via modeling of
drug perturbation effects. Bioinformatics, 35, 3743–3751.
Ritchie,M.E. et al. (2015) limma powers differential expression analyses for RNA-
sequencing and microarray studies. Nucleic Acids Research, 43, e47–e47.
Rivas-Barragan,D. et al. (2020) Drug2ways: Reasoning over causal paths in biological
networks for drug discovery. PLoS Comput Biol, 16, e1008464.
Savitski,M.M. et al. (2014) Tracking cancer drugs in living cells by thermal profiling of
the proteome. Science, 346, 1255784.
Schäker-Hübner,L. et al. (2021) 4-Acyl Pyrrole Capped HDAC Inhibitors: A New
Scaffold for Hybrid Inhibitors of BET Proteins and Histone Deacetylases as
Antileukemia Drug Leads. J Med Chem, 64, 14620–14646.
Seal,R.L. et al. (2023) Genenames.org: the HGNC resources in 2023. Nucleic Acids
Res, 51, D1003–D1009.

77
Shahbazi,J. et al. (2016) The Bromodomain Inhibitor JQ1 and the Histone Deacetylase
Inhibitor Panobinostat Synergistically Reduce N-Myc Expression and Induce
Anticancer Effects. Clin Cancer Res, 22, 2534–2544.
Shao,M. et al. (2017) Structure-based design, synthesis and in vitro antiproliferative
effects studies of novel dual BRD4/HDAC inhibitors. Bioorg Med Chem Lett, 27,
4051–4055.
Shen,H. et al. (2020) S6K1 blockade overcomes acquired resistance to EGFR-TKIs in
non-small cell lung cancer. Oncogene, 39, 7181–7195.
Sidorov,P. et al. (2019) Predicting Synergism of Cancer Drug Combinations Using NCI-
ALMANAC Data. Front Chem, 7, 509.
Spriano,F. et al. (2022) Pharmacologic screen identifies active combinations with BET
inhibitors and LRRK2 as a novel putative target in lymphoma. EJHaem, 3, 764–
774.
Subramanian,A. et al. (2017) A Next Generation Connectivity Map: L1000 Platform and
the First 1,000,000 Profiles. Cell, 171, 1437-1452.e17.
Subramanian,A. et al. (2005) Gene set enrichment analysis: a knowledge-based
approach for interpreting genome-wide expression profiles. Proc Natl Acad Sci U
S A, 102, 15545–15550.
Sun,B. et al. (2015) Synergistic activity of BET protein antagonist-based combinations in
mantle cell lymphoma cells sensitive or resistant to ibrutinib. Blood, 126, 1565–
1574.
Sun,C. et al. (2017) Rational combination therapy with PARP and MEK inhibitors
capitalizes on therapeutic liabilities in RAS mutant cancers. Sci Transl Med, 9,
eaal5148.
Sun,X. et al. (2016) Modeling of signaling crosstalk-mediated drug resistance and its
implications on drug combination. Oncotarget, 7, 63995–64006.
Sun,Z. et al. (2020) DTF: Deep Tensor Factorization for predicting anticancer drug
synergy. Bioinformatics, 36, 4483–4489.
Szakács,G. et al. (2006) Targeting multidrug resistance in cancer. Nat Rev Drug Discov,
5, 219–234.
Szalai,B. et al. (2019) Signatures of cell death and proliferation in perturbation
transcriptomics data—from confounding factor to effective prediction. Nucleic
Acids Research, 47, 10010–10026.
Torres,C. et al. (2006) Proteasome inhibitors shorten replicative life span and induce a
senescent-like phenotype of human fibroblasts. J Cell Physiol, 207, 845–853.

78
Tsao,C.C. and Corn,P.G. (2010) MDM-2 antagonists induce p53-dependent cell cycle
arrest but not cell death in renal cancer cell lines. Cancer Biology & Therapy, 10,
1315–1325.
Tsoi,J. et al. (2018) Multi-stage Differentiation Defines Melanoma Subtypes with
Differential Vulnerability to Drug-Induced Iron-Dependent Oxidative Stress.
Cancer Cell, 33, 890-904.e5.
Vargas,J.E. et al. (2014) Inhibition of HDAC increases the senescence induced by
natural polyphenols in glioma cells. Biochem Cell Biol, 92, 297–304.
Venkatesh,R. et al. (2015) Luotonin-A based quinazolinones cause apoptosis and
senescence via HDAC inhibition and activation of tumor suppressor proteins in
HeLa cells. Eur J Med Chem, 94, 87–101.
Viswanathan,V.S. et al. (2017) Dependency of a therapy-resistant state of cancer cells
on a lipid peroxidase pathway. Nature, 547, 453–457.
Wang,Q. et al. (2009) ERAD inhibitors integrate ER stress with an epigenetic
mechanism to activate BH3-only protein NOXA in cancer cells. Proc Natl Acad
Sci U S A, 106, 2200–2205.
Wei,D. et al. (2019) Comprehensive anticancer drug response prediction based on a
simple cell line-drug complex network model. BMC Bioinformatics, 20, 44.
Wildenhain,J. et al. (2015) Prediction of Synergism from Chemical-Genetic Interactions
by Machine Learning. Cell Syst, 1, 383–395.
Wooten,D.J. and Albert,R. (2021) synergy: a Python library for calculating, analyzing
and visualizing drug combination synergy. Bioinformatics, 37, 1473–1474.
Wu,J. et al. (2022) DRviaSPCN: a software package for drug repurposing in cancer via
a subpathway crosstalk network. Bioinformatics, 38, 4975–4977.
Xu,M. et al. (2018) Senolytics improve physical function and increase lifespan in old
age. Nat Med, 24, 1246–1256.
Yang,B. et al. (2021) MEK Inhibition Remodels the Immune Landscape of Mutant KRAS
Tumors to Overcome Resistance to PARP and Immune Checkpoint Inhibitors.
Cancer Res, 81, 2714–2729.
Yang,J. et al. (2015) DIGRE: Drug-Induced Genomic Residual Effect Model for
Successful Prediction of Multidrug Effects: DIGRE Multidrug Prediction. CPT
Pharmacometrics Syst. Pharmacol., 4, 91–97.
Yu,B. et al. (2022) The synergistic anticancer effect of the bromodomain inhibitor
OTX015 and histone deacetylase 6 inhibitor WT-161 in osteosarcoma. Cancer
Cell Int, 22, 64.

79
Zecha,J. et al. (2023) Decrypting drug actions and protein modifications by dose- and
time-resolved proteomics. Science, 380, 93–101.
Zhang,B. et al. (2020) Class I histone deacetylase inhibition is synthetic lethal with
BRCA1 deficiency in breast cancer cells. Acta Pharm Sin B, 10, 615–627.
Zhang,M. et al. (2015) Drug repositioning for diabetes based on ‘omics’ data mining.
PLoS One, 10, e0126082.
Zhang,X. et al. (2020) Characterization of a dual BET/HDAC inhibitor for treatment of
pancreatic ductal adenocarcinoma. Int J Cancer, 147, 2847–2861.
Zhao,L. et al. (2019) Preclinical Studies Support Combined Inhibition of BET Family
Proteins and Histone Deacetylases as Epigenetic Therapy for Cutaneous T-Cell
Lymphoma. Neoplasia, 21, 82–92.
Zhao,L. et al. (2022) The Robust Tumoricidal Effects of Combined BET/HDAC Inhibition
in Cutaneous T-Cell Lymphoma Can Be Reproduced by ΔNp73 Depletion. J
Invest Dermatol, 142, 3253-3261.e4.


 

80
Appendices
Appendix A: Drug mechanism enrichment analysis improves
prioritization of therapeutics for repurposing

Supplemental Figure 1. DMEA generates MOA rankings that are improved over
single-drug rankings and similar to or better than other tools’ rankings of drug MOA
A) Comparison of DMEA’s MOA rankings to rankings generated by the CMap L1000
Query, DrugEnrichr, and Drugmonizome. DMEA’s rankings are improved compared to
single-drug rankings and similar or improved compared to other tools for evaluating drug
MOAs in all cases. NF: not found. B) DMEA’s MOA rankings also improve upon single-
drug rankings in other cases where CMap cannot provide MOA or PCL analysis and are
similar compared to DrugEnrichr and Drugmonizome which only evaluate a relatively
small input list of drug names.

81


Supplemental Figure 2. DMEA identifies similar MOAs based on gene expression
connectivity scores
Rank-ordered drug lists were generated by querying the CMap L1000 gene expression
perturbation signatures and then analyzed by DMEA. A) Human CD34+ cells treated with
the glucocorticoid agonist dexamethasone (24 h) or untreated (Narla et al., 2011); B)
A2058 cells treated with the PI3K/MTOR inhibitor BEZ235 or DMSO (Brachmann et al.,
2012b). Volcano plots summarizing the NES and -log10(p value) for all tested drug MOAs

82
and mountain plots of the expected MOAs are shown. Red text indicates MOAs with p
value < 0.05 and FDR < 0.25.


Supplemental Figure 3. Overview of drug mechanism enrichment analysis using
WGV molecular classification scores
DMEA can accept an input gene signature when paired with a molecular classification
method such as weighted gene voting (WGV), correlations, and large public databases
of gene expression and drug screens. First, gene weights are calculated as log2(fold-
change) values calculated between two groups of samples from transcriptional data. Only
the x genes with log2(fold-change) with q value < 0.05 are selected for subsequent steps
(maximum of 500 genes). Second, the WGV molecular classifier score is calculated as
the dot product between x gene weights and gene expression values from cancer cell
lines (x genes by n cell lines). Here, we used 327 adherent cancer cell lines from the
CCLE. Third, the Pearson correlation is calculated between the WGV score and drug

83
sensitivity score (i.e., area under the curve (AUC) of cell viability versus drug
concentration from the PRISM database) for each of the 1,351 drugs. An example is
shown for the EGFR inhibitor poziotinib. Fourth, the 1,351 drugs are ranked by Pearson
correlation coefficient, and the ranked list is analyzed using a GSEA-like algorithm to
identify drug mechanisms that are enriched at either end of the ranked list. In total, 85
mechanisms-of-action are tested (Supp. Table 1). An example for EGFR inhibitors is
shown. Each red tick mark represents a drug annotated as an EGFR inhibitor in the
waterfall plot of correlation coefficients.


Supplemental Figure 4. Sensitivity analysis of DMEA pipeline with WGV using
synthetic data
Synthetic gene expression and drug sensitivity data was simulated with varying
perturbations to RNA expression (x-axis) and drug sensitivity score (i.e., AUC, y-axis)
(see Methods). A) Heatmap showing the average DMEA NES from 50 simulations for
varying perturbations of the drug AUC and gene expression data. B) Heatmap showing
the percent of DMEA replicates with FDR q value < 0.25 from 50 simulations for varying
perturbations of the drug AUC and gene expression data.


84


85
Supplemental Figure 5. DMEA identifies enrichment of the EGFR inhibitor MOA and
other MOAs using external signatures of intrinsic EGFR inhibitor resistance
Using gene expression signatures of intrinsic erlotinib (GSE31625 (Balko et al., 2006))
and gefitinib resistance (Coldren et al. (Coldren et al., 2006)), we calculated WGV scores
for 327 adherent cancer cell lines in the CCLE database. For each signature, the WGV
scores were correlated with drug sensitivity data (i.e., AUC) for 1,351 drugs from the
PRISM database. Drugs were then ranked by correlation coefficient, and DMEA was
performed to identify enriched MOAs at either end of the rank-ordered drug list. A)
Volcano plot of NES versus -log10(p value) for DMEA using the GSE31625 signature of
erlotinib resistance. Red text indicates MOAs with p value < 0.05 and FDR < 0.25. B)
Mountain plot showing that DMEA identified the EGFR inhibitor MOA as negatively
enriched in the rank-ordered drug list generated using the GSE31625 signature of
erlotinib resistance. The most negatively correlated EGFR inhibitors are highlighted along
with their correlation coefficients. C) Mountain plot showing that DMEA identified the MDM
inhibitor MOA as positively enriched in the rank-ordered drug list generated using the
GSE31625 signature of erlotinib resistance. The most positively correlated MDM
inhibitors are highlighted along with their correlation coefficients. D) Volcano plot of NES
versus -log10(p value) for DMEA using the Coldren et al. signature of gefitinib resistance.
Red text indicates MOAs with p value < 0.05 and FDR < 0.25. E) Mountain plot showing
that DMEA identified the EGFR inhibitor MOA as negatively enriched in the rank-ordered
drug list generated using the Coldren et al. signature of gefitinib resistance. The most
negatively correlated EGFR inhibitors are highlighted along with their correlation
coefficients. F) Mountain plot showing that DMEA identified the HMGCR inhibitor MOA
as positively enriched in the rank-ordered drug list generated using the Coldren et al.
signature of gefitinib resistance. The most positively correlated HMGCR inhibitors are
highlighted along with their correlation coefficients.

 

86
Appendix B: Mechanism of action landscape reveals shared
sensitivities in large-scale drug screens
Supplemental Table 1. Demographics of adherent cancer cell lines studied for
each tissue type
Tissue Type N
N
RPMI
Media
N DMEM
Media
N Female N Male
Age
(Mean +/-
SD (N))
Any 320 45 11 131 171
53.2 +/-
17.4
(240)
Lung 63 45 11 16 45
55.0 +/-
11.8 (49)
Skin 29 16 9 8 16
48.1 +/-
15.3 (21)
Central
Nervous
System
25 10 7 0 0 NA
Ovary 24 12 6 24 0
54.2 +/-
13.6 (16)
Pancreas 23 15 7 8 14
58.4 +/-
12.6 (18)
Esophagus 20 19 1 6 14
63.9 +/-
10.2 (14)
Breast 17 11 3 17 0
51.4 +/-
12.0 (17)
Urinary Tract 16 3 7 0 0 NA
Liver 16 6 3 1 15
49.9 +/-
10.0 (15)
Upper
Aerodigestive
14 6 4 0 0 NA
Uterus 13 5 1 13 0
66.6 +/-
10.7 (9)
Gastric 11 8 0 3 6
56.3 +/-
15.5 (7)
Kidney 11 4 4 1 9
49.0 +/-
24.5 (6)
Colorectal 10 3 1 0 0 NA
Soft Tissue 8 2 1 5 3
32.0 +/-
24.2 (8)
Bone 6 3 1 3 2
26.8 +/-
27.7 (5)

87
Peripheral
Nervous
System
5 3 1 0 0 NA
Bile Duct 4 4 0 0 2 NA
Thyroid 3 1 2 2 1
65.3 +/-
5.7 (3)


Supplemental Figure 6. MOA relationships can depend on tissue type
Top positively and negatively enriched MOA pairs for all tissue types studied. Tissue
types except for soft tissue (N = 8), bone (N = 6), peripheral nervous system (N = 5), bile
duct (N = 4), and thyroid (N = 3) are represented by more than 10 adherent cancer cell
lines (Supp. Table 1). CNS: central nervous system. CRC: colorectal cancer. PNS:
peripheral nervous system. UA: upper aerodigestive. 
Asset Metadata
Creator Garana, Belinda (author) 
Core Title Computational tools for drug repurposing 
Contributor Electronically uploaded by the author (provenance) 
School Andrew and Erna Viterbi School of Engineering 
Degree Doctor of Philosophy 
Degree Program Chemical Engineering 
Degree Conferral Date 2023-08 
Publication Date 06/12/2023 
Defense Date 06/09/2023 
Publisher University of Southern California (original), University of Southern California. Libraries (digital) 
Tag drug repurposing,enrichment analysis,gene expression analysis,mechanism of action,OAI-PMH Harvest,precision medicine,proteomic analysis,senescence,senolytic,targeted therapeutics 
Format theses (aat) 
Language English
Advisor Graham, Nicholas (committee chair), MacLean, Adam (committee member), Malmstadt, Noah (committee member) 
Creator Email bgarana@usc.edu,garana.belinda@gmail.com 
Permanent Link (DOI) https://doi.org/10.25549/usctheses-oUC113169910 
Unique identifier UC113169910 
Identifier etd-GaranaBeli-11946.pdf (filename) 
Legacy Identifier etd-GaranaBeli-11946 
Document Type Dissertation 
Format theses (aat) 
Rights Garana, Belinda 
Internet Media Type application/pdf 
Type texts
Source 20230612-usctheses-batch-1054 (batch), University of Southern California (contributing entity), University of Southern California Dissertations and Theses (collection) 
Access Conditions The author retains rights to his/her dissertation, thesis or other graduate work according to U.S. copyright law.  Electronic access is being provided by the USC Libraries in agreement with the author, as the original true and official version of the work, but does not grant the reader permission to use the work if the desired use is covered by copyright.  It is the author, as rights holder, who must provide use permission if such use is covered by copyright. 
Repository Name University of Southern California Digital Library
Repository Location USC Digital Library, University of Southern California, University Park Campus MC 2810, 3434 South Grand Avenue, 2nd Floor, Los Angeles, California 90089-2810, USA
Repository Email uscdl@usc.edu
Abstract (if available)
Abstract Most patients are ineligible for targeted therapies and left to choose from other treatments which can have deadly side effects. Considering this, my work has focused on repurposing existing targeted therapies for new disease applications (i.e., drug repurposing).

Though many drug repurposing algorithms exist, they typically only evaluate individual drugs which can be prone to off-target effects and output a long list of drugs, making it difficult for researchers to prioritize drug targets. To remove the assumption that each drug works by its known mechanism(s), I developed software to group drugs by shared mechanism-of-action (MOA). The utility of this approach is that drug MOAs are only identified as candidates if most of the drugs annotated with that MOA are ranked as strong candidates. By aggregating results across many drugs to evaluate drug MOAs instead of just individual drugs, researchers can be more confident in their treatment approach.

Altogether, my computational tools improve prioritization of drugs for repurposing and can even identify therapeutics which may mitigate resistance. My software does not require the user to possess programming knowledge and is freely available online for the public. While experimental conditions can affect drug toxicity, the computational tools I have developed will help researchers prioritize drugs to evaluate further. 
Tags
drug repurposing
enrichment analysis
gene expression analysis
mechanism of action
precision medicine
proteomic analysis
senescence
senolytic
targeted therapeutics
Linked assets
University of Southern California Dissertations and Theses
doctype icon
University of Southern California Dissertations and Theses 
Action button