A Transcriptome-Based Precision Oncology Platform for Patient–Therapy Alignment in a Diverse Set of Treatment-Resistant Malignancies

A major goal of precision cancer medicine (PCM) is to improve clinical outcomes by leveraging the molecular-level properties of a tumor—as encoded by mutational, gene expression, epigenetic modification, and proteomic profiles—to accurately predict sensitivity to candidate therapeutic agents. Application of PCM principles may help generate responder-enriched cohorts for clinical trials when predictions are conserved across a substantial fraction of patients (1, 2), and even help prioritize personalized treatments on an individual patient basis.

Systematic application of the current PCM paradigm is largely predicated on two complementary approaches. The first one (oncogene addiction) aims to identify targeted therapies based on the presence of activating genetic alterations in druggable oncoproteins (3); the second (immunotherapy) is based on the discovery that specific tumor-initiated immunosuppressive programs can be abrogated by pharmacologic targeting of immune checkpoints or by sensitizing the immune system to tumor antigens (4).

Despite the remarkable clinical success of these approaches within specific tumor subtypes (5), many tumors may lack actionable genetic alterations, fail to respond to therapy, or develop drug resistance, suggesting an acute need for complementary approaches targeting nononcogene tumor dependencies (6). In particular, despite its critical role in tumor subtype stratification, use of transcriptome-based approaches in precision medicine has lagged.

We and others have shown that, within each tumor histology, cancer cells adopt only a relatively limited, discrete, and remarkably stable repertoire of molecularly distinct transcriptional states (7). These states are mechanistically controlled by tightly autoregulated tumor-checkpoint modules (TCM; ref. 8), comprising a small, highly conserved set of master regulator (MR) proteins responsible for canalizing the effect of mutations in their upstream pathways (7–9). About 30% of MR proteins were shown to represent tumor-essential, nononcogene dependencies, either individually (8) or in combination (10, 11). Indeed, genetic (10–12) or pharmacologic (13–15) inhibition of MR proteins has been shown effective in inverting the activity of TCM MRs, resulting in the abrogation of tumor viability in vitro and in vivo. TCM inversion can be effectively achieved due to their hyperconnected, heavily autoregulated nature, supporting their ability to behave as homeostatic on/off control modules (ref. 12; bioRxiv 2020.10.27.357269). As such, individual MRs—and the TCMs they comprise—represent an actionable class of nononcogene dependencies.

In this study, we test two approaches to leverage these conceptual advances for therapeutic purposes, by targeting either individual, pharmacologically actionable candidate MRs with a high-affinity inhibitor (OncoTarget), or the entire TCM (OncoTreat). These methods rely on the ability to accurately measure the transcriptional activity of the regulatory proteins that maintain tumor cell state, by RNA-sequencing (RNA-seq) profile analysis using the virtual proteomics by enriched regulon (VIPER) algorithm (16), which has been shown to compare favorably with antibody-based protein measurements (17).

Specifically, OncoTarget uses VIPER to identify the most aberrantly activated proteins for which a high-affinity inhibitor is already available, thus representing a straightforward, mutation-agnostic extension of the oncogene addiction paradigm. Indeed, aberrant protein activity can result not only from activating mutations in the encoding gene but also from mutations and signals in upstream pathways (7), as well as from autocrine, paracrine, and endocrine signals (18). In contrast, OncoTreat leverages large-scale perturbational profiles—i.e., RNA-seq profiles representing the cell's response to drug treatments—in MR-matched cell lines selected as high-fidelity (cognate) models of a patient tumor. These profiles support direct experimental assessment of context-specific drug mechanism of action (MoA) and TCM-activity inversion in drug-treated versus vehicle control–treated tumor cells (13). Compared with other related efforts, such as SynergySeq (19), which predict drug sensitivity based on the greatest divergence between drug perturbation and tumor transcriptomic signatures, our methods capitalize on robust network-based assessment of functional protein activity and are built on a broad mechanistic framework incorporating insights on tumor-specific essentiality of MR proteins. Further, we leverage a careful selection of cognate models to assess context-specific drug MoA, and, most critically, we provide in vivo validation.

In this first-in-class application of MR-based tools to predict drug sensitivity, we designed a tumor- and mutation-agnostic noninterventional clinical study to enroll patients with advanced malignancies, across multiple histologies, who had progressed on several lines of therapy (the N of 1 study at Columbia University, IRB-AAAN7562). We present a series of results from the treatment of early passages of the first seven patient-derived xenograft (PDX) models established, which were used to assess the overall preclinical efficacy of OncoTarget- and OncoTreat-predicted drugs. For the latter, we also performed pharmacodynamic assays to assess in vivo recapitulation of TCM-activity inversion, as predicted from in vitro perturbations. Our results demonstrate that OncoTarget and OncoTreat are highly predictive of anti­tumor drug activity in vivo, supporting further development of these clinically actionable tests—both of which are New York and California Department of Health approved and Clinical Laboratory Improvement Amendments (CLIA) compliant [Columbia University Irving Medical Center (CUIMC) pathology department]—for biomarker-driven PCM clinical trials. We discuss the strengths, limitations, and practical challenges encountered in implementing and validating these tools in a clinical context.

To assess OncoTarget and OncoTreat's ability to predict drug sensitivity in tumors from pretreatment RNA-seq profile analysis, we designed an innovative, proof-of-concept clinical study with a preclinical endpoint. The N of 1 study enrolled patients with advanced malignancies refractory or intolerant to standard-of-care treatment, including several rare tumors (Supplementary Table S1). Due to trial design, selected participants had generally progressed on most if not all standard-of-care therapies by time of enrollment, and lacked actionable genetic alterations suggesting potential targeted drug efficacy.

Clinically indicated biopsies or tumor resections were performed at the request of the treating oncologist; consent was required to allow a portion of the fresh tumor tissue from biopsies with ≥70% cellularity to be processed for RNA-seq profiling and transplantation into immunodeficient mice. Here, we report the results of 35 distinct drug arms—including 21 OncoTarget-predicted and 22 OncoTreat-predicted drugs, 8 of which were predicted by both methods—in the first seven, consecutively established PDX models that could be expanded for in vivo drug testing. Given its research nature, the protocol did not require patients to be treated with predicted drugs. Instead, subsequent treatment was chosen by the patient's oncologist using currently available approaches for cancers that have progressed on standard treatments, e.g., mutational profiling, off-label drug use, and referral to therapeutic clinical trials when eligible.

Note that the study was neither designed nor sufficiently powered to assess the efficacy of each individual predicted drug. Rather, our goal was to validate the overall ability of the two methodologies (tests) to predict drugs that elicit in vivo response.

The seven PDX models included three basal-like breast cancers (BLBC; BC-32398, BC-97359, and BC-50291), a pancreatic ductal carcinoma (PDA; PAC-05647), a colon adenocarcinoma (COAD; CAR-23659), a KIT^WT/PDGFR^WT gastrointestinal stromal tumor (GIST) harboring both a KRAS^G12D mutation and a germline SDHB deletion (GIST-81050), and a recurrent WHO grade II anaplastic meningioma (MEN; CNS-16474). Clinical characteristics of the seven patients are summarized (Fig. 1A), including extensive prior therapies and the results of targeted genomic sequencing, when performed at the discretion of the treating oncologist. For two patients (CNS-16474 and PAC-05647), mutational profiling was not available because the oncologist did not expect these tumor types to harbor clinically actionable mutations. Notably, four of the seven subjects had progressed on at least three lines of treatment before enrollment and all seven represented aggressive, drug-resistant tumors.

Due to the constraints associated with the expansion and use of early PDX passages for preclinical therapeutic studies, we applied several practical selection criteria to further prioritize predicted drugs for in vivo validation. Specifically: (a) Only drugs classified as antineoplastic agents were considered, (b) drugs were eliminated if the patient had been previously treated with them, (c) when perturbation profiles were available from suitable models—i.e., for GIST, meningioma, and breast cancer—OncoTreat-predicted drugs were selected over OncoTarget-predicted ones with comparable P values, to ensure an overall balanced number of tested predictions from the two methods, (d) drugs predicted for the patient but not for the corresponding PDX model were eliminated, and (e) when multiple inhibitors sharing the same canonical mechanism were identified—e.g., multiple HDAC or topoisomerase II inhibitors—only the top FDA-approved and/or the most statistically significant one was selected.

The 35 unique patient-specific drugs prioritized for in vivo validation, as well as prediction rationale and dosing schedule, are summarized in Supplementary Table S2 (see Supplementary Table S3 for further details on drug prediction and selection in each of the seven models).

We have developed two complementary, RNA-based assays to transform VIPER-based (16) tumor sample-specific protein activity profiles into actionable drug response predictions (see schematics in Fig. 1B; Supplementary Fig. S1A–S1C).

The first one (OncoTarget) identifies aberrantly activated proteins (threshold P ≤ 10⁻⁵, Bonferroni corrected, as measured by VIPER) for which a high-affinity inhibitor is currently available (see Methods). To assess target actionability, we analyzed DrugBank (RRID:SCR_002700), the SelleckChem database (RRID:SCR_003823), published literature, and public information from drug development pipelines, resulting in a curated list of 180 proteins representing validated, high-affinity targets of clinically relevant small-molecule compounds (Supplementary Table S4). These include proteins that are rarely if ever mutated in cancer, such as topoisomerases, chromatin remodeling enzymes, and proteins aberrantly activated by autocrine, paracrine, or endocrine signals.

The second one (OncoTreat) leverages a large compendium of RNA-seq profiles, generated to represent the transcriptional response of cell lines to a comprehensive repertoire of antineoplastic agents. These allow the identification of compounds capable of inverting the transcriptional activity of a TCM module (TCM-inverters)—as defined by the 25 most activated (25↑) and 25 most inactivated (25↓) candidate MR proteins in a patient tumor—at a conservative statistical significance threshold [P ≤ 10⁻⁵, Bonferroni corrected, by one-tailed analytic-rank based enrichment analysis (aREA); Fig. 1B; Supplementary Fig. S1C; refs. 13, 16; see Methods]. As such, rather than using a priori knowledge, OncoTreat predicts TCM-inverter drugs based on drug MoA assessed de novo from experimental perturbational assays—i.e., based on the differential activity of 2,556 regulatory proteins in drug-treated versus vehicle control–treated cells.

The number of candidate MRs in a TCM was selected based on the average number of MRs necessary to integrate the effect of genetic alterations in their upstream pathways, as assessed in ref. 8. Specifically, we had reported that, for the vast majority of tumor subtypes in The Cancer Genome Atlas (TCGA), there is rapid saturation of mutational events in pathways upstream of the first 1 to 100 candidate MRs of each individual tumor, with 50 MRs sufficient to account for 80% of all mutations in all but 5 of 20 tumor types analyzed (COAD, HNSC, SKCM, STAD, and OV). However, we also showed that drug prediction is extremely robust and reproducible for any TCM size ranging from n = 10 to 200 activated and inactivated MRs, with mean Spearman pairwise correlation ranging from r_s = 0.89 in GIST-81050 to r_s = 0.98 in BC-32398 (Supplementary Fig. S2A–S2E).

Drug MoA and TCM-inversion potential were assessed based on perturbational profiles in selected, high-fidelity (cognate) cell lines. These were identified based on their ability to recapitulate the TCM-activity signature—i.e., the activity of the top and bottom 25 most differentially active MRs—of the greatest fraction of a histology-matched patient cohort. Thus, cognate cell lines are meant to represent a biological surrogate of the tumor of interest only in terms of providing an optimal in vitro context to assess tumor-relevant drug mechanisms of action, thus maximizing statistical power. Specifically, assessing MR protein activity decrease or increase following drug perturbation is best accomplished in cells where these are already significantly activated or inactivated, respectively. Likewise, to account for potential tumor drift effects in PDX passaging, we assessed whether the PDX models used in the study also recapitulated their corresponding patients’ TCM-activity signatures. As such, we assess model fidelity to a human tumor based on TCM MR enrichment (i.e., top 25↑+25↓) of the human tumor in differentially active and inactive proteins in the model, respectively (OncoMatch analysis, refs. 13, 20; bioRxiv 2019.677435; threshold Bonferroni P ≤ 10⁻¹⁰; see Methods).

To illustrate the selection process, in Supplementary Fig. S3A we show the MR-based fidelity of the top 12 breast cancer cell lines identified as candidate cognate models of BLBC tumors in TCGA, as annotated by PAM50 classification (21, 22). Cell lines were selected from a total of 97 profiled breast cancer cell lines, unbiased to receptor status or PAM50 classification, as part of a comprehensive repository that included both the Cancer Cell Line Encyclopedia (CCLE; ref. 23) and the Genentech Cell Line Screening Initiative (gCSI; ref. 24). BT20 emerged among the top five candidates based on the number of patient tumors (78 of 173) whose TCM-activity signature it recapitulated (Bonferroni P ≤ 10⁻¹⁰, by OncoMatch). Because cognate cell lines are identified strictly based on TCM-activity recapitulation, we do not expect them to necessarily recapitulate other phenotypic or even transcriptomic characteristics of their matched patients. For instance, although none of the top 12 cell lines identified for patients with BLBC would be classified as luminal, several of them were claudin-low or mesenchymal-like (e.g., MDAMB231, SUM159PT, BT549, and CAL120), unclassifiable (e.g., BT20, HCC1395), or even ER-negative, HER2-amplified/enriched (e.g., JIMT1, HCC1954) by transcriptome-based classifications (25, 26). This is consistent with the fact that patients with BLBC in TCGA presented conserved activity of the most differentially active MRs, regardless of claudin or HER2 status (7).

In Supplementary Fig. S3B, we show the fidelity of four cognate cell lines that were used to generate drug perturbation profiles to support OncoTreat analyses for five study patients. Although ASPC1 was also identified as a high-fidelity model for the pancreatic tumor in the study, ASPC1-based perturbation profiles were not completed in time for drug prediction and in vivo validation. Nonetheless, for completeness, we are sharing the ASPC1 drug perturbation data as part of this study. As shown, BT20 represents a high-fidelity model for tumor BC-32398 [normalized enrichment score (NES), 14.5, P = 10⁻⁴⁸] and BC-97359 (NES 8.0, P = 10⁻¹⁵), but not for BC-50291 (NES −3.9, P = 1); both GIST cell lines, GIST430 and GISTT1, were identified as high-fidelity models for GIST-81050 (P < 10⁻⁴⁰) despite not harboring the patient SDHB^Del/KRAS^G12D alterations, but rather canonical KIT mutations; finally, the meningioma cell line IOMM was borderline for CNS-16474 (NES 3.3, P = 0.0005).

On average, ∼350 drugs were profiled in each cognate cell line model, including 138 FDA-approved antineoplastics, about 170 late-stage experimental drugs in phase II and III oncology clinical trials, as well as a variable number of additional compounds from diverse libraries with cell line-specific EC₅₀ ≤2 μmol/L (Supplementary Table S5; see Methods). Cells were harvested at 24 hours following perturbation with each compound at two sublethal concentrations—the 48-hour EC₂₀ and one tenth of this concentration, as determined by 10-point dose–response curves—as well as at 6 hours (in selected cell lines) to assess early MR activity changes. We used the highest sublethal drug concentration to reveal signatures optimally reflective of the drug's MoA rather than nonspecific effector proteins in stress or death pathways (27, 28). To avoid testing drugs at nonphysiologically relevant concentrations, we also capped concentrations at their C_Max, defined as the peak serum concentration for the drug's maximum tolerated dose, from published pharmacokinetic studies in humans, when available.

Multiplexed, low-depth (1–2M reads) RNA-seq profiles were generated for each compound using the PLATE-seq technology (29), which supports low-cost, fully automated pooled library generation from 96 or 384-well plates. Drug MoA—representing the drug-mediated differential activity of all regulatory proteins—was then assessed by VIPER analysis of drug-treated versus vehicle-treated (DMSO) samples. In Fig. 2A, as a representative example, we show the hierarchical clustering of VIPER-assessed drug MoA in BT20. Interestingly, although several MoA-related drugs were identified within the same clusters (e.g., topotecan and irinotecan in cluster 1, sorafenib and vandetanib in cluster 5, and dasatinib, ponatinib, and nilotinib in cluster 6), several drugs induced similar differential protein activity profiles, despite having distinct high-affinity targets (e.g., thioguanine, a guanine analogue; sorafenib, a multikinase inhibitor; and vorinostat, a pan-HDAC inhibitor, in cluster 5). This suggests that cellular networks in BT20 cells may effectively canalize drug MoA toward a small number of relatively distinct cellular responses. Additionally, as shown, several drugs produced highly conserved MoA at multiple concentrations (e.g., topotecan at 0.076 and 0.76 μmol/L, crizotinib at 0.72 and 7.2 μmol/L, flubendazole at 1.44 and 14.5 μmol/L, among several such examples, suggesting high reproducibility of these assays, as also previously shown; refs. 13, 29).

Drug MoA profiles from the appropriate cognate cell line(s) were used to generate OncoTreat predictions for patient and PDX tumors. All regulatory proteins represented in the drug MoA were ranked from most inhibited to most activated, thus providing an optimal reporter assay to assess TCM inversion. TCM-inversion assessment in vitro for the 22 OncoTreat-predicted drugs that were prioritized for in vivo validation is shown in Fig. 2B.

Several groups have described clonal drift that occurs with sequential passages in PDX models (30). As a result, to minimize drift, we performed all therapeutic studies in the earliest feasible passage, P1–P5. Additionally, prior to commencing therapeutic testing, we assessed both model fidelity and patient/PDX conservation of drug prediction, as proposed in ref. 31. Following successful engraftment of tumors (P0 passage), we performed RNA-seq and subsequent VIPER, OncoTarget, and OncoTreat analyses to determine (a) the MR-based fidelity of the PDX tumor to the patient tumor and (b) conservation of drug predictions. Drugs predicted from patient sample analysis, but not predicted by analysis of the PDX, were used only in the therapeutic study if alternative options were not available.

Six of the seven PDX models, GIST-81050, BC-32398, CAR-23659, BC-97359, CNS-16474, and PAC-05647, met the predefined fidelity threshold (Bonferroni P ≤ 10⁻¹⁰, by OncoMatch), with NES ranging from 13.8 to 17.3 (Fig. 3A; Supplementary Fig. S4). In fact, in GIST-81050 and CAR-23659 there was almost perfect fidelity to patient TCM activity, whereas in BC-32398, BC-97359, CNS-16474, and PAC-05647 there were a handful of MRs with different activity rank between patient tumor and PDX. In BC-50291, however, there was no appreciable fidelity (P = 0.89). Consequently, in BC-50291 we tested drug predictions for the patient tumor that were not conserved in the PDX.

A total of 35 individual drugs predicted by the analysis were evaluated in individual PDX therapeutic arms—including 21 OncoTarget-predicted and 22 OncoTreat-predicted, 8 of which were predicted by both methods. In a few cases, the same drug was predicted and tested in more than one model. After the expansion of PDX models for the therapeutic study, animals were enrolled once tumor volumes reached ∼100 mm³. Models were treated, and tumor volume measurements were recorded up to a fixed end-of-study time point with a median of 29 days (range, 21–30 days). A waterfall plot of the response of each individual mouse, grouped by PDX model, to prioritized and negative control drugs is shown in Fig. 3B and C and comparisons summarized in Fig. 3D and E; responses are summarized in Supplementary Table S6. Overall, across all predictions, disease control—comprising stable disease (SD), partial response (PR), or complete response (CR)—was observed for 30 of 35 predictions (85.7%).

We evaluated overall treatment response, at the end of treatment, for drugs that were OncoTarget-predicted, OncoTreat-predicted, and predicted by both methods. Compared with vehicle control, there were significant differences in disease control rate (DCR)—defined as the rate of SD + PR + CR—and objective response rate (ORR)—defined as the rate of PR + CR—for both OncoTarget-predicted (P_DCR < 10⁻³, P_ORR = 0.03, by Fisher exact test) and OncoTreat-predicted (P_DCR < 10⁻³, P_ORR = 0.01; Fig. 3B–E) drugs. Drugs predicted by both had response rates comparable to drugs predicted by just one methodology (DCR = 76%, n = 19/25; ORR = 20%, n = 5/25). Due to unanticipated toxicity possibly related to study treatment (e.g., tumor ulceration in a breast cancer PDX), 3 mice were unevaluable for response in the OncoTarget cohort (BC-97359, n = 2 for MK-2206 arm, n = 1 for panobinostat arm), whereas 1 mouse was unevaluable in the OncoTreat cohort (BC-97359, n = 1 for irinotecan arm), and 3 mice were unevaluable for response in the Both cohort (GIST-81050, n = 2 for daunorubicin arm, n = 1 for topotecan arm). Although cautioning that the study was neither designed nor sufficiently powered to evaluate the efficacy of individual drug predictions, for completeness, we provide growth curves for each drug (Supplementary Fig. S5A–S5C).

As a set of appropriate negative controls, we evaluated response to randomly selected antineoplastic drugs that were not statistically significant by either OncoTarget or OncoTreat (when available) analysis (i.e., P = 1 from both methods). We note that these are not true random controls but rather potential false-negative predictions meant to assess the algorithms’ predictive power. Four PDX models (GIST-81050, CAR-23659, PAC-05647, and BC-50291) were thus treated with 13 negative control drugs (Supplementary Table S2) and vehicle control. To avoid bias due to differences in tumor growth rates, these drugs were tested in conjunction with an additional vehicle control arm in each PDX, thus providing a tumor growth-independent reference. Disease progression was observed across all models treated with negative control drugs and statistically indistinguishable from concurrent vehicle control (Fisher exact P = 0.6; Supplementary Table S6).

Cumulative Kaplan–Meier analysis of animals in the OncoTarget, OncoTreat, and OncoTreat + OncoTarget cohorts was performed (Fig. 4A). The analysis demonstrates highly statistically significant improvement in disease control using agents predicted by either methodology, compared with vehicle control (P < 10⁻⁴, two-tailed log-rank test). In contrast, there was no statistically significant difference between animals in the negative control and passage-matched vehicle control cohorts (P = 0.38, log-rank; Fig. 4B). Additionally, the treatment-to-control ratio (∆T/∆C%: the ratio of change in volumes from baseline in the treatment/control arms), which corrects for baseline differences in model-passage growth rate, was significantly improved in the OncoTarget (mean ∆T/∆C% = 14.3%; 95% CI, −1.2 to 29.9, P = 0.004, two-tailed Mann–Whitney U test) and OncoTreat (mean ∆T/∆C% = 17.3%; 95% CI, −0.8 to 35.3, P = 0.014) treated cohorts, compared with the negative control cohort (mean ∆T/∆C% = 46.1%; 95% CI, 25.5–61.1), with overall statistical significance (OncoTarget vs. OncoTreat vs. negative control, P = 0.002, two-tailed ANOVA; Fig. 4C).

Pharmacodynamic (PD) studies are a critical aspect of drug development to elucidate drug MoA and to characterize primary and acquired drug resistance. PD assessment from early, on-treatment samples helps to determine: (a) whether effective, OncoTreat-predicted drugs recapitulate in vivo the TCM inversion that occurs in the cognate cell line(s) (i.e., MoA conservation) and (b) whether the failure of OncoTreat-predicted drugs corresponds to inability to conserve TCM inversion in vivo, perhaps due to pharmacokinetic factors, or occurs despite TCM inversion, for instance, due to later cell adaptation or clonal selection. The opportunity to investigate the second point was limited in the study, because most OncoTreat-predicted drugs demonstrated strong antitumor activity.

Samples for PD assessment were procured from two mice per treatment arm, for the four PDX models treated with at least three OncoTreat-predicted drugs—GIST-81050, BC-32398, CNS-16474, and PAC-05647. Mice were randomly selected for early sacrifice, independent of tumor size, 3 hours following the third dose, and were excluded from response assessment. The four selected models had the highest patient tumor fidelity (P < 10⁻¹⁰; Fig. 3A) and were thus well-suited to evaluate MoA conservation. TCM inversion was assessed by VIPER analysis of RNA-seq signatures comparing drug-treated versus vehicle control–treated PDX tumor samples.

Overall, the vast majority of OncoTreat-predicted drugs with available PD samples—i.e., 15 of 18 (83%)—significantly recapitulated in vivo (P < 10⁻⁵, by one-tailed aREA) the TCM inversion predicted from cognate cell line perturbational assays, in vitro, with P-values ranging from 10⁻⁵ (teniposide in PAC-05647) to 10⁻⁴⁰ (daunorubicin in BC-32398; Fig. 5A–D). The three drugs that failed to conserve MoA in vivo included abiraterone in CNS-16474, which was borderline for achieving disease control, belinostat in PAC-05647, which achieved disease control by end-of-study, and daunorubicin in GIST-81050, which also achieved disease control (Supplementary Table S6). As belinostat was the only HDAC inhibitor evaluated for PD effect, one is tempted to speculate if the early time point was inadequate for an epigenetic modifying drug to fully implement its in vivo effect. Intriguingly, daunorubicin demonstrated strong TCM inversion in the BC-32398 model but not GIST-81050, where it inverted the 25↓ MRs but failed to invert the 25↑ MRs, and yet had strong antitumor activity in both models. Conversely, one drug failed to achieve disease control (homoharringtonine in the CNS-16474 model) despite recapitulating TCM inversion in vivo.

As expected, four of the negative control drugs—alpelisib, selinexor, serdemetan, and pacritinib—failed to achieve significant TCM inversion in vivo (Fig. 5E) and did not achieve disease control by end of study. TAE684 did induce TCM inversion in vivo (P = 10⁻¹⁶) in the early on-treatment sample, yet still failed to achieve disease control.

In summary, 15 of 18 OncoTreat-predicted drugs, including 13 of 16 (81%) that induced disease control, recapitulated significant TCM inversion in early on-treatment PD samples. This is consistent with our hypothesis that inference of drug-induced TCM inversion in carefully selected models, whether they be cell lines or PDXs, is robust and a feature of the hyperconnected and autoregulated nature of TCMs.

The proposed PCM framework discussed in this study is uniquely suited to identify therapeutic alternatives in real-world scenarios, even for rare cancers lacking actionable mutations and standard-of-care options. Calcifying Nested Stromal Epithelial Tumor (CNSET) is an exceptionally uncommon primary hepatic tumor that occurs in children and young adults, with only about 40 cases reported in the literature (32). Although localized disease is often effectively cured with surgery, recurrent and de novo metastatic disease demonstrates chemotherapy resistance, and there are no proven therapeutic options (33–36). We thus report the observed clinical outcome for a CNSET case, where OncoTarget was used off-trial to help guide treatment selection when all other options were exhausted.

A 14-year-old male reported a 1-month history of abdominal pain and fatigue. A computed tomography (CT) scan of the chest and abdomen revealed a large hepatic mass with multiple satellite liver tumors and pulmonary metastases. The mass was biopsied, and histopathologic evaluation was consistent with CNSET. The family was initially hesitant to initiate chemotherapy, but the patient developed progressive hepatomegaly, anorexia, weight loss, constipation, and anemia in the subsequent 3 months. Memorial Sloan Kettering IMPACT (37), a targeted next-generation sequencing panel covering 468 genes, was performed on a biopsy specimen and demonstrated a CTNNB1 hotspot mutation, TERT promoter gain-of-function mutation, and an NTRK3 point mutation not known to predict response to currently available TRK inhibitors. The estimated tumor mutational burden was only 2.6 per megabase, predicting a low likelihood of response to immune-checkpoint inhibitors.

Based on the observation that the tumor shared biological features with Wilms’ tumor (CTNNB1 and TERT mutations, WT1 and β-catenin expression by IHC), the patient was treated with 10 neoadjuvant cycles of a “Wilms’ tumor like-regimen,” including five cycles of vincristine/irinotecan, four cycles of vincristine/dactinomycin/doxorubicin, and one cycle of cyclophosphamide/topotecan (38, 39). There was a PR to chemotherapy and the patient successfully underwent debulking surgery. Postoperative chemotherapy was complicated by the development of severe colitis, and the family elected to discontinue systemic therapy. Over the next 6 months, there was evidence of significant disease progression in the liver and lungs (Fig. 6A), and the patient developed biliary obstruction and transaminitis that made him ineligible for clinical trials and precluded the use of most chemotherapy agents.

Given the lack of remaining viable therapeutic options, tumor tissue was sent for the CLIA-compliant OncoTarget test. The most significantly activated targetable protein was PDGFRB (Bonferroni P = 10⁻⁷, Fig. 6B). After discussing the results with the family, including the absence of clinical data on targeting PDGFRB in this exceedingly rare malignancy, we decided sunitinib would be the best candidate drug, given its relative selectivity for PDGFRB over other kinases (RRID:SCR_003823), accessibility as a drug approved by the FDA in 2006, and safety data in the context of impaired hepatic function and pediatric patients. The patient had a PR to the first cycle of sunitinib (6 weeks) which deepened by the end of cycle 3 (Fig. 6C). Remarkably, this patient who had rapidly progressing treatment-refractory cancer has had a durable response and remains on sunitinib, now for 2 years from his original presentation with only mild side effects such as fatigue. Although OncoTarget did predict response, a caveat is that we cannot be certain of the mechanism. As pediatric patients do not routinely undergo repeat biopsies when responding to a therapy, we could not confirm the treatment effect on PDGFRB activity or changes in the activity of other canonical targets of sunitinib, including KIT and less potently VEGFRs.

The OncoTarget and OncoTreat approaches can identify multiple candidate drugs for the treatment of most tumors, a majority of which induced disease control in PDX models. Given the inherent conservation of MR and TCM activity within cancer subtypes, as identified by protein activity-based cluster analysis (7), it is reasonable to expect that subsets of patients with predicted sensitivity to the same drugs (pharmacotypes) should emerge from these analyses. Indeed, the majority of TCGA cancer cohorts were effectively stratified into 2 to 7 pharmacotype clusters by joint OncoTreat/OncoTarget analysis.

Figure 7A provides a representative example for the BLBC subtype of breast cancer in TCGA, with the pharmacotype cluster assignment of the three BLBC patients shown. Consistently, tumors are predicted to respond to distinct, pharmacotype-specific drugs. Full pharmacotype stratification of patient cohorts representing the cancers discussed in this study is provided in Supplementary Fig. S6A–S6D.

The oncogene addiction (3) and immunotherapy (4) paradigms have illuminated PCM's potential to meaningfully improve the outcomes of some patients. These approaches are straightforward to implement in the clinic. However, methodologies to predict response to the full repertoire of available antineoplastics, based on rapid and inexpensive RNA-seq profiling, are currently underdeveloped and would represent a welcome addition. Here, we present a first-in-class application of two complementary and scalable transcriptome-based tests, OncoTarget and OncoTreat, to predict drug sensitivity on an individual patient basis. Critically, both methodologies are New York and California Department of Health approved and CLIA compliant, and predictions can be generated within 2 to 3 weeks after biopsy, including for OncoTreat if perturbational profiles in a cognate cell line are available.

We encountered a number of practical challenges during this large-scale effort that are worth noting. First, as discussed, both OncoTarget and OncoTreat depend on the VIPER-based measurement of protein activities. We have published extensively on the accuracy and reproducibility of VIPER (8, 16), which compares favorably with gene and protein abundance estimation. However, as, unlike binary determination of genetic alterations, VIPER produces a fully quantitative assessment subject to inherent biological and technical noise, we decided to use a highly conservative statistical significance threshold for predicting drug sensitivity (Bonferroni P ≤ 10⁻⁵). As a result, our analyses may likely miss drugs that could have been beneficial (false negatives). In addition, protein activity measurement by VIPER, as reported (16), may be incorrect for 10% to 20% of proteins, leading to both false-positive and false-negative drug predictions.

Second, there are several recognized limitations to the use of PDX models. Tumors undergo clonal evolution under various selection pressures, including available nutrients, organ-specific environment, immune editing, pharmacologic treatment, and growth kinetics (40). The selection pressures in an immunocompromised PDX may differ from those of the original human host, thus resulting in rapid drift. To minimize this challenge, we restricted therapeutic studies to early passages and confirmed model fidelity on the basis of preserved TCM activity. Due to cost and constraints in propagating and expanding models at a single early passage, we could not test all statistically significant drugs predicted by OncoTarget and OncoTreat and had to implement practical selection criteria. We also could not test a more comprehensive list of negative control drugs, such as drugs the patient had previously or was currently receiving, or a broader panel of randomly selected drugs, blinded to the OncoTarget and OncoTreat results, across all models. Furthermore, the role of the immune system and thus the effect of immunotherapy is not evaluable in these immunodeficient models.

Additionally, the development of a xenograft from implanted patient tumor tissue may be influenced by sundry factors—including histology, mitotic rate, and fraction of cancer stem cells (41)—whereas the timeline between implantation and completion of drug testing can vary greatly (e.g., 3 to 18 months in our experience), with some tumor types rarely taking in the xenograft. Consistent with the experience of other groups, less than half of attempted PDXs could be established (Supplementary Table S1), and there was obvious tumor-type bias (42–48). Although we would not anticipate these factors to bias comparisons of response between OncoTarget/OncoTreat-predicted and negative control antineoplastic drugs, it nonetheless limits our ability to confirm the utility of our methodologies in certain tumor types. After a planned interim statistical analysis demonstrated highly significant findings, and given budgetary constraints, we tested drugs in only a fraction (the first seven) of established PDXs. Nonetheless, we believe that our systematic study of 35 algorithmically predicted drug arms and 13 negative control drug arms represents the largest such effort to date aimed at validating computational drug-sensitivity predictions in vivo.

We further note that there is no consensus on the optimal PDX endpoint that best predicts response in patients. Although clinical response assessments using RECIST are based on unidimensional measurements, bidirectional caliper measurements to estimate tumor volume are standard practice in mouse studies, but are prone to significant interobserver variability (49, 50). RECIST-defined criteria for progressive disease (PD; i.e., >20% increase in the sum of 5 target lesions) translate to a 73% increase in spherical volume (51). When considering the various factors that can contribute to variability in tumor measurements, our definition of PD (100% tumor volume increase relative to baseline) remains within reasonable margins of error (± 1–2 mm) and aligns with volumetric definitions of progression (51–53). Hence, we believe there are advantages to the volumetric end points we present here, and, reassuringly, OncoTarget and OncoTreat demonstrated significant benefit using several complementary analytic metrics.

Third, OncoTarget and OncoTreat predictions are unique and vary in number between tumors. In some cases, there are very few predictions by one or both methods and predictions may not be necessarily conserved in the PDX model, whereas in other cases there are abundant predictions. As such, it is not always possible to test an identical number of predictions of each type across models. This motivated the design of our study to evaluate the efficacy of drugs predicted by the two methods in parallel, and thus we cannot draw any robust conclusions about differences without additional studies. Further, although selecting a conservative statistical threshold generally reduces the number of predicted drugs to a reasonable number (e.g., 5 to 25), the significant variability between tumors is an important consideration for future clinical deployment.

Finally, a requisite component of generating OncoTreat predictions is the completion of one or more perturbation screens in a cognate model (e.g., cell line or organoid). Although we do not anticipate postperturbation cell viability assessments in such models to meaningfully predict clinical response, we do intuitively note that a critical aspect of assessing the ability of a drug to induce patient-specific TCM inversion is that TCM inversion cannot be assessed in the cell line model if the corresponding MRs are not consistently differentially active. Thus, identifying high-fidelity models with conserved TCM activity was a critical element of the study. In addition to their dependence on model availability, it should be noted that the generation of comprehensive, high-throughput perturbation screens is also cost and time prohibitive. As a result, perturbation screens had not yet been completed for all the tumors considered in this study. This is also a real-world consideration, as evidenced for the CNSET case study, where we could only rely on OncoTarget predictions. Fortunately, since completion of this study, we have generated a Pan-cancer Activity by Enrichment Analysis database (PanACEA), comprising perturbational profiles of about 350 drugs in 25 cancer cell lines, each representing a cognate model matched to an MR-based subtype within 15 human malignancies (Supplementary Table S7). Access to this extensive resource, which greatly extends the generalizability of the proposed approach, will be made available via a companion publication; see also ref. 54.

We also make note of several practical advantages of our methodology. First, these algorithms are tumor- and mutation-agnostic, thus supporting their application to the vast majority of human malignancies. Based on benchmarking of 11,289 primary tumors from TCGA, as well as more than 100 prospectively collected samples from patients with metastatic and treatment-refractory cancers, these tests can efficiently prioritize multiple candidate drugs for virtually every tumor, thus allowing further drug prioritization based on drug approval status, insurance reimbursement, oncologist experience, prior treatment, and toxicity. We further note that, whereas our current study was restricted to solid tumors, the approach is equally applicable to hematologic malignancies. Indeed, we have previously reported on several MR-based insights in lymphoma and leukemia (15, 55), including OncoTreat-based prediction of drug synergy (27), and these malignancies will constitute a focus of future studies.

Critically, as the number of validated OncoTreat and OncoTarget predictions increases, including in a clinical setting (56), the need to validate additional predictions in flawed preclinical models will decrease, conceivably leading to direct clinical utilization. This would significantly increase clinical utility, especially in the advanced metastatic setting. Indeed, we have already initiated and enrolled subjects to therapeutic trials based on our methodologies, including ricolinostat + nab-paclitaxel in metastatic breast cancer (NCT02632071; ref. 56), entinostat in neuroendocrine tumors (NCT03211988; ref. 13), and the HIPPOCRATES umbrella trial in pancreatic cancer (NCT04476537), without the need for PDX testing.

Second, the ability to stratify patients into a small number of well-defined pharmacotypes—i.e., tumors with shared predicted drug sensitivity—in virtually all analyzed cancer cohorts, supports prioritization and evaluation of drug predictions through standard clinical trial mechanisms. Pharmacotype analysis provides two critical insights: (i) It helps identify drugs consistently predicted across a substantial fraction of patients in a specific tumor type (e.g., BLBC), thus supporting mechanism-based hypotheses for preclinical and clinical trials, and (ii) it helps disregard singleton drugs—i.e., drugs predicted for a single patient in a cohort—as potential false positives. Basket and umbrella trials may incorporate OncoTarget and OncoTreat as companion diagnostic biomarkers. Lending support to this approach, an OncoTarget-based trial testing the HDAC6 inhibitor ricolinostat in metastatic breast cancer recently concluded with virtually complete validation of sensitivity predictions [area under the curve (AUC) = 1; ref. 56], thus paving the road to additional trials with a similar design. OncoTarget, in particular, could expand the pool of eligible patients for basket trials using targeted agents, whereas OncoTreat can be used as the basis of umbrella trials where patients with one or more cancer types can be assigned to preselected drug treatment arms based on pharmacotype classification. To illustrate this concept, we show the pharmacotypes identified by our analysis in BLBC and how this could be leveraged to design an umbrella trial (Fig. 7A and B).

Third, the resources required to perform OncoTarget and OncoTreat analysis (i.e., RNA-seq and modest computational resources) are relatively affordable and scalable. There is a critical unmet need to diversify enrollment to clinical trials, which further affects rapid adoption of new technologies in different healthcare settings (56). The ability to perform the analysis from archival formalin-fixed, paraffin-embedded tissue will allow access to clinical trial participation and clinical testing across a broad range of health care facilities. As our methodology identifies candidate drugs for the vast majority of patient tumors, we hope that it will remove at least one roadblock and source of disappointment to patients enrolling in basket/umbrella trials, due to some not matching any drug arm.

Finally, our methodology can capture longitudinal changes, such as we have shown for metastatic progression in breast cancer and neuroendocrine tumors (13, 57), histologic transformation in follicular lymphoma (15), therapy resistance in T-cell acute lymphoblastic leukemia and breast cancer (55, 58), and reprogramming of castrate-resistant prostate adenocarcinoma to a neuroendocrine state (59). OncoTarget and OncoTreat could thus be potentially applied longitudinally to adapt therapy to the dynamic nature of tumor evolution.

There are also potential limitations. Tumor heterogeneity and stroma infiltration are only partially addressed by the current approach—specifically by limiting the analysis to high tumor cellularity samples, a criterion usually met by the mostly metastatic samples collected in our study. This issue will be mitigated in the future by the increasing availability of single-cell gene-expression profiling technologies. Indeed, due to the robustness of VIPER to low sequencing depth, we have recently shown its applicability to measure protein activity in single cells, comparing favorably with antibody-based measurements and recapitulating bulk measurements (17, 60), thus supporting drug predictions for specific tumor subpopulations (bioRxiv 2022.02.28.482410), including immuno­suppressive, nontransformed populations, such as regulatory T cells (bioRxiv 2022.02.22.481404).

We conclude by stating the obvious. It would be naive to expect that all, or even a majority of, predicted drugs will be clinically effective in humans. Yet, our study and some initial clinical confirmations (56) suggest that these tests may represent a novel PCM option for mechanism-based prioritization of existing, clinically relevant drugs that is not currently available to cancer patients, and thus potentially complement successful approaches such as targeted therapy and immunotherapy.

To support context-specific regulatory protein activity inference by VIPER, we have generated comprehensive molecular interaction networks (interactomes) using the Algorithm for the Reconstruction of Accurate Cellular Networks (ARACNe; refs. 61, 62), although other suitable algorithms may be used. The networks were reverse-engineered by ARACNe from ≥100 RNA-seq profiles of human cancer tissue from (a) TCGA and (b) for meningioma and neuroendocrine tumors, from Columbia University collected datasets (Supplementary Table S8). TCGA RNA-seq level 3 data were downloaded from NCI Genomics Data Commons (63). Raw counts were normalized and variance stabilized by fitting the dispersion to a negative-binomial distribution as implemented in the DESeq2 R-package (RRID: SCR_000154; refs. 64).

ARACNe was run with 100 bootstrap iterations using an input set of candidate regulators, including (a) 1,877 transcription factors annotated in the Gene Ontology (GO; ref. 65) “Molecular Function database” as DNA-binding transcription factor activity (GO:0003700), DNA binding (GO:0003677), transcription regulator activity (GO:0030528), or regulation of transcription, DNA-templated (GO:0003677 and GO:0045449); (b) 677 transcriptional cofactors manually curated from genes annotated as transcription coregulator activity (GO:0003712), plays a role in regulating transcription (GO:0030528), or regulation of transcription (GO:0045449); and (c) 3,895 genes encoding for signal transduction proteins, dually annotated in the GO “Biological Process database” as Signal Transduction (GO:0007165) and in the GO “Cellular Component database” as either intracellular (GO:0005622) or plasma membrane (GO:0005886). The data processing inequality parameter of ARACNe was set to 0 and the mutual information (MI) threshold was set to P = 10⁻⁸. The mode of regulation was computed based on the correlation between regulator and target gene expression as previously described (16). The version of raw counts and generated networks used in our work are provided (see Key Resources Table in the supplementary material).

We have previously extensively validated the VIPER analysis algorithm as a highly robust and specific tool for the accurate inference of regulatory protein activity in a tissue context-dependent manner (8, 15, 16). VIPER leverages accurate tissue-specific gene regulatory networks, such as those produced by ARACNe (61, 62), to measure differential protein activity from bulk or single-cell gene-expression signatures (see Key Resources Table). Specifically, akin to a multiplexed gene-reporter assay, VIPER measures a protein's differential transcriptional activity through a probabilistic enrichment framework that assesses the enrichment of its activated and repressed transcriptional targets (regulon) in genes over- and under-expressed in a sample of interest compared with a set of control samples (reference model).

We considered multiple options for reference models to apply VIPER to tumor samples for the purposes of predicting drug response. Although organ-matched normal tissue expression profiles, such as those available from the genotype-tissue expression (GTEx) resource (66), and cohort-matched tumor and tumor-adjacent profiles were considered, such references would potentially overemphasize proliferative and cell-cycle signals and underemphasize lineage-specific tumor vulnerabilities, respectively. Ultimately, we prospectively opted to use the pan-TCGA tumor dataset as our reference model to run VIPER and the downstream drug prediction algorithms.

Thus, for each cancer sample, we generate a differential gene-expression signature (DGES)—computed as the gene-wise relative expression to the distribution of the expression of that gene across 11,289 TCGA samples—and expressed as its quantile relative to the reference model. Genes whose expression is not captured (quantified) in TCGA are excluded from the DGES. Next, VIPER computes enrichment scores for the targets of each regulatory protein in the DGES, using the aREA test (16), thus determining if the protein's transcriptional footprint is over-, under-, or normally represented in the tumor's DGES. Each enrichment score is computed by a two-tailed approach, rank-transforming the DGES for the positively regulated target genes and then inverting the signature to compute enrichment of the repressed target genes. The enrichment score is a weighted average based on targets’ ranks in the signature and the protein-target gene MI assignments in the network, with stronger interactions contributing more to the score. Subsequently, statistical significance for the enrichment score is computed by an analytic approach. Under the null hypothesis, the enrichment score is normally distributed with a mean of zero and variance proportional to 1/(n − 1), where n is the number of targets, scaled to the standard variable. This approximates shuffling the genes in the signature at random (67), outputting a NES. Although aREA does not directly account for biological correlation between the expression of various genes (i.e., covariance) in the given biological context, a very conservative NES/P-value threshold is empirically used for downstream analyses. Importantly, due to VIPER's use of robust reporter sets to compute the activity of each protein, moderate biases in RNA profiles based on RNA extraction method, sequencing instrument, and reference transcriptome version used to map reads, all of which can affect genes captured and abundance estimates, generally have minimal effect on protein activity assessment. Likewise, due to this use of robust reporter sets, for rare cancer types not represented in TCGA, which might have lineage-specific aberrant expression of a subset of genes, we overall do not anticipate a significant impact on protein activity assessment. When cancer type–specific networks are not available, we use an integrated network approach as implemented in metaVIPER (60, 68).

The most differentially active proteins identified by VIPER comprise candidate MR proteins. Critically, we have shown that MR proteins are ultra-conserved within each of the 112 tumor subtypes and play a key role in integrating the effect of genomic alterations to implement a stable tumor cell state, as shown for 20 tumor types in TCGA (7). VIPER reproducibility is extremely high, such that Spearman correlation of activity profiles generated from RNA-seq at 30M to as low as 50K read depth is ρ ≥ 0.8 (16) even though correlation of the underlying gene-expression profiles is low ρ ≤ 0.3. Although VIPER is uniquely suited for assessing regulatory proteins that directly control gene expression, including transcription factors, cofactors, and chromatin remodeling enzymes, we have shown that the algorithm is equally effective in monitoring the activity of signaling proteins (16) and cell-surface markers (60).

Through the use of (a) DrugBank (69), (b) the SelleckChem database (RRID:SCR_003823), (c) published literature, and (d) publicly available information on pharmaceutical company drug development pipelines, we have curated a refined list of 180 actionable proteins representing validated targets of high-affinity pharmacologic inhibitors, either FDA approved or in clinical trials (Supplementary Table S4). This manually curated target-drug database is dominated by signaling proteins and established oncoproteins, as expected from the bias in druggability assessment and past focus in drug development efforts. Pharmacologic agents with narrow therapeutic indices—such as those targeting neurotransmitters, ion channels, and vasoactive drugs—were purposefully removed from the database as less likely to be successfully repurposed in oncology.

OncoTarget simply analyzes the VIPER outputted protein activity measurements for these 180 actionable proteins and provides a multiple-testing corrected significance value for the corresponding NES. We used a conservative threshold (Bonferroni P < 10⁻⁵) to identify candidate proteins eliciting essentiality when targeted by a pharmacologic inhibitor, for in vivo validation. Based on this threshold, the average TCGA tumor yields 15 unique, significant OncoTarget predictions, ranging from an average of n = 4.5 in adrenocortical carcinoma (ACC) to n = 28.7 in renal clear cell carcinoma (KIRC). Tumors representing conserved subtypes tend to have conserved OncoTarget predictions, leading to the identification of subsets of patients whose cancer is predicted to respond to the same set of drugs (pharmacotypes; see Fig. 7; Supplementary Fig. S6A–S6D for the tumor types discussed herein).

We have previously described OncoTreat (13) as a methodology to systematically elucidate compounds capable of significantly inverting the activity of the 25↑+25↓ MRs comprising the TCM that regulates metastatic progression in enteropancreatic and rectal neuroendocrine tumors. In that study, OncoTreat identified entinostat, a class I histone deacetylase inhibitor, among 105 profiled drugs, as a highly effective TCM-inverter drug (13), leading to a clinical trial that is currently accruing patients (NCT03211988).

For the current study, we more broadly adapted OncoTreat to identify TCM-inverter compounds, by comparing individual tumor samples to the entire TCGA repository, as a reference model, as described above. Specifically, for each tumor type considered in the study, the analysis proceeds through the following steps:

1. Identifying cell lines (typically 1 or 2) jointly representing high-fidelity (i.e., cognate) in vitro models for a majority of samples in the corresponding TCGA cohort (e.g., pancreatic ductal adenocarcinoma, PAAD), based on TCM-recapitulation (Bonferroni P < 10⁻⁵), as assessed by enrichment analysis; see next section on OncoMatch. Upon acquisition of cell lines (see Key Resources Table), RNA-seq profiles of untreated cell lines were generated at CUIMC to confirm reproducibility of expression and MR/TCM profiles. All cells were inspected for Mycoplasma contamination using the ABM PCR Mycoplasma Detection Kit.

2. Generating RNA-seq profiles of each cognate cell line, from Step 1, at 24 hours following perturbation with a library of clinically relevant oncology drugs. Drugs were titrated at their maximum sublethal concentration (i.e., 48 hours EC₂₀), as determined by 10-point dose-response profiles. Profiled drugs included FDA-approved and late-stage experimental oncology drugs (in phase II and III clinical trials; see Supplementary Table S5). For completed screens, we have attempted to be all inclusive, excluding therapeutic antibodies due to a lack of availability and immunotherapy agents due to lack of appropriateness to screen in vitro. Additional miscellaneous compounds with EC₅₀ ≤ 2 μmol/L in the selected cell line were also included. Most compounds were purchased from SelleckChem or Tocris. DMSO was selected as a universal in vitro solvent (vehicle). Multiplexed, low-depth (1M to 2M reads) RNA-seq profiles were generated using 96-well plates via the PLATE-seq technology, using fully automated microfluidics for increased throughput and reproducibility (29). Eight DMSO-treated controls were included in each plate, to avoid plate-dependent batch effects and to mitigate the inherent variability of DMSO treatment.

3. Generating a subproteome-wide context-specific drug MoA for each drug, as represented by the differential activity of each protein in drug-treated versus vehicle control (DMSO)–treated cells. Differential protein activity was assessed by VIPER analysis using a tissue-matched gene regulatory network produced by ARACNe; see above.

4. Identifying sample-specific candidate MRs and the TCMs they comprise, by VIPER analysis of the sample's DGES, compared with the set of TCGA samples (reference model).

5. Finally, prioritizing pharmacologic agents based on the statistical significance of the enrichment of the tumor sample's TCM-activity signature (i.e., 25↑+25↓ MRs) in proteins inactivated and activated in drug versus DMSO-treated cells, respectively, with negative NES indicating TCM inversion (Bonferroni P < 10⁻⁵, one-tailed aREA). The number of candidate MR proteins (n = 50) used to assess TCM inversion, which for this step was restricted to only transcription factors and cofactors, was selected because we have shown that, on average, across all of TCGA, the vast majority of functionally relevant genomic events can be found upstream of the top 50 VIPER-inferred candidate MR proteins (7).

Model fidelity was assessed based on the statistical significance of the TCM-activity conservation between a human-derived sample and a model-derived sample. For computing protein activity in cell line models, we first generated an analogous DGES comparing each cell line against a large repository of cancer cell lines (reference model), which includes both the CCLE (23) and the gCSI (24). Next, we computed the enrichment of the TCM-activity signature (25 most active and 25 most inactive patient tumor-specific MRs) in differentially active and inactive proteins in the model (OncoMatch). The aREA (16) test was again used, but any suitable enrichment analysis algorithm could be substituted. The analysis was used to (a) select optimal cell lines for the generation of perturbational profiles that effectively track the activity of tested drugs on TCM proteins and (b) assess the fidelity of PDXs prior to validation of drugs predicted from the human sample. All cell lines used in these analyses were resequenced on-site to capture any potential drift effects. A conservative threshold (Bonferroni P < 10⁻¹⁰) was used to identify high-fidelity models.

All animals were maintained under barrier conditions and all experiments were performed in accordance with and approval of the Memorial Sloan Kettering Cancer Center (MSKCC) Institutional Animal Care and Use Committee (IACUC; protocol #16-08-011) and CUIMC IACUC (AAAF5850). Patient tumor tissue was collected under the CUIMC Institutional Review Board (IRB)–approved protocol AAAA7562, with written informed consent provided prospectively by the subject and conducted in accordance with recognized ethical guidelines under the U.S. Common Rule. Generation of PDX models was performed under the MSKCC IRB-approved protocol #17-387 and CUIMC IRB protocol AAAJ5811. PDX models were established as previously described (52). In summary, fresh tumor tissue was fragmented and implanted subcutaneously into nonobese/severe combined immunodeficiency IL2Rγ null, hypoxanthine phosphoribosyltransferase (HPRT)-null (NSGH) mice (Jackson Labs, IMSR catalog no. JAX:012480, RRID: IMSR_JAX:012480; see Key Resources Table) and tumor engraftment monitored by visual and manual inspection.

Engrafted tumors were measured twice weekly with calipers, and drug treatment was initiated when tumor volume (TV) reached ∼100 mm³ (TV = width² × ½ length; see Key Resources Table for sourcing of drugs). Early-passage animals (passages 1–5) were used for all therapeutic studies. Aligned with clinical response criteria, treatment response was categorized as CR (>95% volumetric reduction from baseline or no palpable tumor), PR (>50% reduction), SD (<50% reduction and no more than 100% increase), or PD (>100% increase from baseline). DCR is defined as the sum of CR, PR, and SD responses divided by the total responses assessed. ORR is defined as the sum of CR and PR responses divided by the total responses assessed. Tumor responses were assessed at the end of the therapeutic study period. Comparisons of disease control and objective response rates across treatment groups (OncoTarget alone, OncoTreat alone, Both) were performed using a Fisher exact test. The Mann–Whitney–Wilcoxon method was used to compare differences in the distribution of relative TV between OncoTarget or OncoTreat cohorts and vehicle control. Vardi's test was used to evaluate differences in AUC between treatment groups across models (70, 71). Curves for disease control using the Kaplan–Meier estimator were compared and analyzed using the log-rank test. ΔT/ΔC%, computed as the change in TV from baseline for drug-treated mice divided by the change in TV from baseline for vehicle control–treated mice, was determined for OncoTarget, OncoTreat, and negative control cohorts, and differences were evaluated using two-way ANOVA. All relevant statistical analyses were performed using R software (v3.5.0) or GraphPad Prism [v8.4.1 (RRID:SCR_002798)]. Statistical significance was defined as a P < 0.05.

Samples for pharmacodynamic assessment were procured from two mice per treatment arm, for the four PDX models treated with at least three OncoTreat-predicted drugs—GIST-81050, BC-32398, CNS-16474, and PAC-05647. Mice were randomly selected for early sacrifice, independent of tumor size, 3 hours following the third dose, and were excluded from response assessment. We performed RNA-seq and subsequent VIPER on paired drug-treated versus vehicle control–treated PDX tumor samples. TCM inversion was assessed based on the statistical significance of the enrichment of the TCM-activity signature (i.e., 25↑+25↓ MRs of the patient tumor) in proteins inactivated and activated in drug-treated versus vehicle control–treated PDX tumors, respectively, again using aREA, although alternative enrichment tests may be used. Negative NES indicates TCM inversion (Bonferroni P < 10⁻⁵).

Original tumor bulk RNA-seq data generated as part of this work have been deposited to the Gene Expression Omnibus (GEO: GSE212854; Key Resources Table). Relevant high-throughput drug perturbation data used for de novo assessment of drug MoA is available on Figshare (FS; https://figshare.com/s/d77d19eb326b8a8656b8). RNA-seq of a cohort of 102 meningioma tumors sequenced at CUIMC and used to generate a meningioma network have been deposited to GEO (GSE212377). The version of TCGA processed data (RNA-seq counts) used to generate networks and perform other analyses presented herein is available on FS (https://figshare.com/s/ad114ea4b274a523bb4a). The ARACNe-generated networks used to run VIPER are available on FS (https://figshare.com/s/5d1ffd9f8b2e86e37ed6). The RNA-seq counts and OncoTarget results for the CNSET case report are available on FS (https://figshare.com/s/6327d1857201ce657d47). Additional requests should be directed to corresponding author A. Califano.

D. Diolaiti reports grants from the NCI and CureSearch during the conduct of the study, as well as personal fees from Seagen outside the submitted work. C. Karan reports personal fees from DarwinHealth and other support from DarwinHealth outside the submitted work. R.D. Carvajal reports personal fees from Alkermes, Bristol Myers Squibb, Castle Biosciences, Delcath, Eisai, Hengrui, Ideaya, Immunocore, InxMed, Iovance, Merck, Novartis, Oncosec, Pierre Fabre, PureTech Health, Regeneron, Sanofi/Genzyme, Sorrento Therapeutics, and TriSalus outside the submitted work. B.S. Henick reports personal fees from AstraZeneca, Ideaya, Jazz Pharmaceuticals, Sorrento Therapeutics, Genentech/Roche, and OncLive, grants from NexImmune, Genentech/Roche, Johnson & Johnson, the Bristol Myers Squibb Foundation, the V Foundation, and Stand Up To Cancer outside the submitted work, as well as a patent for use of TCF1 as a blood biomarker of immunotherapy response pending. J.Y. Hou reports personal fees from IMVAX, GOG, and Foundation Medicine outside the submitted work. F.M. Iwamoto reports personal fees from AbbVie, Alexion, Gennao Bio, Guidepoint Global, Kiyatec, Massive Bio, Medtronic, Merck, MimiVax, Novocure, PPD, Regeneron, Tocagen, Xcures, Oncoceutics, Advarra, and Praesidia Biotherapeuthics, grants and personal fees from the Ben and Catherine Ivy Foundation, grants from Bristol Myers Squibb, and other support from Celldex, Forma Therapeutics, Merck, Northwest Biotherapeutics, Novocure, Sapience Therapeutics, Tocagen, and Astex outside the submitted work, as well as support from the William Rhodes and Louise Tilzer-Rhodes Center for Glioblastoma at NewYork-Presbyterian Hospital, Keep Punching, and the Glioblastoma Foundation. J.G. Jurcic reports grants from the NCI Cancer Target Discovery and Development Program (U01 CA217858), the Cancer Systems Biology Consortium (U54 CA209997), NIH Shared Instrumentation Grants (S10 OD012351 and S10 OD021764), and an NCI Cancer Center Support Grant (P30 CA013696) during the conduct of the study, as well as grants from Forma Therapeutics, Ionis Pharmaceuticals, Sumitomo Pharma, Gilead Sciences, Blueprint Medicines, Arog Pharmaceuticals, and Celularity, personal fees from AbbVie, Jazz Pharmaceuticals, Novartis, and Syros Pharmaceuticals, and grants and personal fees from Bristol Myers Squibb outside the submitted work. N. Lamanna reports grants and personal fees from AbbVie, AstraZeneca, BeiGene, Genentech, and Eli Lilly/Loxo, grants from Octapharma, Oncternal, MingSight, and TG Therapeutics, and personal fees from Janssen/Pharmacyclics and Nurix outside the submitted work. A.B. Lassman reports grants from the NIH/NCI and foundations during the conduct of the study, as well as grants, personal fees, and nonfinancial support from various sources. Specifically, support includes grants from the William Rhodes and Louise Tilzer-Rhodes Center for Glioblastoma at New York-Presbyterian Hospital, Hearst Foundations, the Michael Weiner Glioblastoma Research Into Treatment Fund, and the NIH/NCI/National Center for Advancing Translational Sciences P30CA013696, UL1TR001873, and UG1CA189960 grants. Other support includes nonfinancial support from Pfizer, personal fees from Bioclinica as an expert blinded independent reviewer of clinical and imaging data for a Bristol Myers Squibb–sponsored trial, personal fees from Affinia Therapeutics and the Gerson Lehrman Group (GLG), grants, personal fees, and nonfinancial support from Karyopharm, personal fees from Leal Therapeutics, grants, personal fees, and nonfinancial support from Novocure, grants, personal fees, and nonfinancial support from Orbus, personal fees from RadMD as an expert blinded independent reviewer of clinical and imaging data for a Roche-sponsored trial, personal fees from Sapience and Vivacitas Oncology, personal fees and nonfinancial support from Bayer, grants and nonfinancial support from Global Coalition for Adaptive Research, Novartis, and QED, nonfinancial support from Oligo Nation, Helsinn, the RTOG Foundation, the Society for Neuro-Oncology, the FDA, Aeterna Zentaris, DelMar, NextSource, Corden, Kazia, Biohaven, and Vigeo, grants and nonfinancial support from VBI, Agios/Servier, MJH Holdings (for Physicians Education Resource, PER), Genentech/Roche, AbbVie, Abbott, and Incyte, grants, personal fees, and nonfinancial support from Oncoceutics/Chimerix, and personal fees from Fondazion AIRC (Italian Foundation for Cancer Research), Guidepoint Global, Clinical Education Alliance, and Aptitude Health outside the submitted work. G.A. Manji reports grants from Plexxikon Pharm, Merck Pharm, Regeneron Pharm, BioLine Rx, Arcus Biosciences, and Pfizer, grants, personal fees, and other support from Roche/Genentech, and personal fees from Ipsen during the conduct of the study. A.I. Neugut reports personal fees from Otsuka Pharmaceuticals, grants from Otsuka Pharmaceuticals, and personal fees from GSK, Value Analytics, United Biosource Corp, and Hospira outside the submitted work, and is a member of the medical advisory board for EHE Intl. C.A. Shu reports personal fees from AstraZeneca, Genentech, Janssen, Mirati Therapeutics, Takeda, and Arcus Biosciences outside the submitted work. J.D. Wright reports grants from Merck and personal fees from UpToDate outside the submitted work. K. Kalinsky reports that his spouse has stock in EQRX, Grail, Array BioPharma, and Pfizer (prior employee), as well as advisory/consulting activities for Genentech/Roche, Immunomedics, Seattle Genetics, Oncosec, 4D Pharma, Daiichi Sankyo, Puma Biotechnology, Mersna, Menarini Silicon Biosystems, Myovant Sciences, and Takeda. P.A. Sims reports personal fees from Guardant Health outside the submitted work. M.J. Alvarez reports personal fees from DarwinHealth outside the submitted work, as well as a patent for PCT15/248,975 issued and licensed to DarwinHealth and a patent for PCT15/249,069 issued and licensed to DarwinHealth. A.L. Kung is on the scientific advisory board for DarwinHealth but has not received compensation for scientific advisory board activities over the last 36 months. A. Califano reports personal fees and other support from DarwinHealth outside the submitted work, as well as a patent for the VIPER algorithm (US 20170076035 A1) issued, licensed, and with royalties paid from DarwinHealth. No disclosures were reported by the other authors.

P.S. Mundi: Data curation, formal analysis, investigation, visualization, methodology, writing–original draft, writing–review and editing. F.S. Dela Cruz: Data curation, formal analysis, investigation, visualization, writing–original draft, writing–review and editing. A. Grunn: Data curation, investigation. D. Diolaiti: Data curation, investigation, visualization. A. Mauguen: Formal analysis, visualization. A.R. Rainey: Investigation. K. Guillan: Investigation. A. Siddiquee: Investigation. D. You: Investigation. R. Realubit: Data curation, investigation. C. Karan: Data curation, supervision, investigation. M.V. Ortiz: Writing-original draft. E.F. Douglass: Data curation, investigation. M. Accordino: Investigation. S. Mistretta: Investigation. F. Brogan: Investigation. J.N. Bruce: Investigation. C.I. Caescu: Project administration. R.D. Carvajal: Investigation. K.D. Crew: Investigation. G. Decastro: Investigation. M. Heaney: Investigation. B.S. Henick: Investigation. D.L. Hershman: Investigation. J.Y. Hou: Investigation. F.M. Iwamoto: Investigation. J.G. Jurcic: Investigation. R.P. Kiran: Investigation. M.D. Kluger: Investigation. T. Kreisl: Investigation. N. Lamanna: Investigation. A.B. Lassman: Investigation. E.A. Lim: Investigation. G.A. Manji: Investigation. G.M. McKhann: Investigation. J.M. McKiernan: Investigation. A.I. Neugut: Investigation. K.P. Olive: Investigation. T. Rosenblat: Investigation. G.K. Schwartz: Conceptualization, investigation. C.A. Shu: Investigation. M.B. Sisti: Investigation. A. Tergas: Investigation. R.M. Vattakalam: Investigation. M. Welch: Investigation. S. Wenske: Investigation. J.D. Wright: Investigation. P. Canoll: Investigation, methodology, supervision. H. Hibshoosh: Investigation. K. Kalinsky: Conceptualization, investigation. M. Aburi: Project administration. ­P.A. Sims: Resources, supervision. M.J. Alvarez: Software, validation, methodology, ­writing–original draft, writing–review and editing. A.L. Kung: Conceptualization, supervision, funding acquisition, writing–original draft, writing–review and editing. A. Califano: Conceptualization, supervision, funding acquisition, writing–original draft, writing–review and editing.

This work was supported by the NCI Cancer Target Discovery and Development Program (U01 CA217858), the Cancer Systems Biology Consortium (U54 CA209997), and NIH Shared Instrumentation Grants (S10 OD012351 and S10 OD021764), all to A. Califano. Also, this research was funded in part through the NCI Cancer Center Support Grant (P30 CA013696) to CUIMC. Also, this research was funded in part through CureSearch for Children's Cancer Acceleration Initiative Award to A. Kung and the NCI Cancer Center Support Grants to MSKCC (P30 CA008748).

The publication costs of this article were defrayed in part by the payment of publication fees. Therefore, and solely to indicate this fact, this article is hereby marked “advertisement” in accordance with 18 USC section 1734.