Pharmacogenomic landscape of patient-derived tumor cells informs precision oncology therapy

Abstract

Outcomes of anticancer therapy vary dramatically among patients due to diverse genetic and molecular backgrounds, highlighting extensive intertumoral heterogeneity. The fundamental tenet of precision oncology defines molecular characterization of tumors to guide optimal patient-tailored therapy. Towards this goal, we have established a compilation of pharmacological landscapes of 462 patient-derived tumor cells (PDCs) across 14 cancer types, together with genomic and transcriptomic profiling in 385 of these tumors. Compared with the traditional long-term cultured cancer cell line models, PDCs recapitulate the molecular properties and biology of the diseases more precisely. Here, we provide insights into dynamic pharmacogenomic associations, including molecular determinants that elicit therapeutic resistance to EGFR inhibitors, and the potential repurposing of ibrutinib (currently used in hematological malignancies) for EGFR-specific therapy in gliomas. Lastly, we present a potential implementation of PDC-derived drug sensitivities for the prediction of clinical response to targeted therapeutics using retrospective clinical studies.

Figure 1.. Patient tumor and derived cell resources for pharmacogenomics analysis.

(a) Overview of the procedure for pharmacogenomics analysis in patient tumor-derived cells (PDCs). A total of 462 PDCs from 14 cancer types were isolated. Genomic contexts were analyzed to identify somatic variants and/or gene expression profiles. Comprehensive genomic profiling from 14 cancer lineages are summarized in a circos plot, demonstrating detailed data structure and size for each type of available molecular data. Short-term cultured PDCs have undergone drug sensitivity screening to 60 molecular targeted agents. Clinical feasibility of PDC screening-guided precision therapy was evaluated. (b) Three-dimensional bubble plot showing the frequency of somatic nonsynonymous mutations exclusively in tissue (red; left axis), exclusively in PDC (black; right axis), and in common to the two (yellow; upper axis). 122 samples with DNA sequencing data from tumor tissue and PDC samples were considered in this analysis. (c) Comparison of mRNA expression profiles between primary tissue and PDC on 24 paired samples with matched RNASeq data. Spearman correlations of mRNA expression between tissue and PDC are shown as a heat map. Paired samples are located along the diagonal.

Figure 2.. Therapeutic landscape of PDCs and lineage-specific responses

(a) Tumor lineage-specific drug association identified using 60 compounds (n=462 biologically independent samples). Two-sided wilcox rank sum test was applied to determine the relative differences of drug sensitivity between specific tumor type and the rest. Only significant associations are marked (q-value < 0.05). Drugs are ordered based on their known targets. (b) A volcano plot representation of TDA analysis showing the magnitude (Fold change, x-axis) and significance (TDA q value, y-axis) of all tumor-drug associations (n=462 biologically independent samples). Each circle represents a single tumor-drug interaction and the size is proportional to the cohort size of that tumor. (c) Distribution of gastric, glioma PDCs and BYL719 drug AUC profile over the topological representation of PDCs (n=462 biologically independent samples). Each node represents a set of PDCs with similar AUC profiles. A PDC can appear in several nodes, and two nodes are connected by an edge if they have at least one PDC in common. The P values were calculated using the pearson correlation test between the fraction distribution of gastric or glioma cell lines and mean AUC values of BYL719 drug over the nodes and they were adjusted using BH method. (d) Violin plots measures the activity level of PI3K–AKT–mTOR pathway on gastric and GBMs using TCGA RNASeq datasets (n=100 biologically independent samples for each group). We adopt the enrichment score derived from ssGSEA analysis as assessment. The P value is calculated from two-sided wilcox rank-sum test. Horizontal lines within the violin plot represent 0.25, 0.50, and 0.75 quantiles.

Figure 3.. Pharmacogenomic interactions in PDCs.

(a) A volcano plot representation of correlation analysis showing the magnitude and significance of gene-drug associations (n=360 biologically independent samples). (b) Waterfall plot enumerating significant associations between KRAS mutation and drug sensitivity (n=360 biologically independent samples). Two horizontal dashed lines indicate statistical significance. (c) SW480 was treated with DMSO (control) or trametinib (0.1 or 1 mM), followed by incubation with two EGFR inhibitors, dacomitinib (left) and gefitinib (right). Cell viability for each dose was normalized to DMSO or trametinib (0.1 or 1 mM) treatment only cells. (d) Probability distribution of drug-target families over the topological network. Each node represents a set of drugs with similar AUC profiles. A drug can appear in several nodes, and two nodes are connected by an edge if they have at least one drug in common. Colors of the nodes correspond to mean RGB values of drug families. Ibrutinib belongs to three nodes on the network encompassed with an oval. (e) Drug sensitivities to ibrutinib in 67 PDCs. The red color in the heat map represents sensitivity, while the blue color indicates resistance. EGFR alterations including genomic amplification, vIII, and expression are shown. (f) Kaplan-Meier survival plots for P2.T (EGFR amp/vIII) orthotopic mice model. Once intracranial model was established, Vehicle (0.5% methylcellulose) or iburitnib (50mg/kg/day) was administrated orally (PO) for 5 consecutive days and 2 days of resting period per each cycle (n=8 per group). P-values: a,b,e, two-sided wilcox rank sum test. P-values: f, two-sided Log-rank test.

Figure 4.. Genomic and transcriptomic correlates of panobinostat sensitivity

(a) Schematic overview for dNetFS (Diffusion kernel based Network method for Feature Selection of drug sensitivity). In brief, dNetFS integrates genomic/pharmaceutical data, protein-protein interaction network, and prior knowledge of drug-targets interaction to prioritize genetic and gene expression features of PDCs that predict drug response. (b) Predictive features of panobinostat response identified by the dNetFS are plotted for their frequency and effect size. Associations are colored in red for expression features and blue for others. Node size is proportional to the single drug–feature linear correlation. (c) Scatter plot showing linear correlation between panobinostat AUC and HDAC4 expression (left panel) and SIK2 expression (right panel) (n=69 biologically independent samples). The correlation coefficient and the P values were obtained using pearson correlation test. (d) Drug response assessment of panobinostat (1.2nM or 4.8nM) with shRNA-mediated knockdown of HDAC4 or NT (non-target) (left panel) and SIK2 or NT (right panel). Cell viability for each dose was normalized to sole transduced cells only. Data are mean ± s.d. of n=3 technical replicates. Experiments were repeated three times with similar results. *P ≤ 0.05, **P ≤ 0.01, ***P ≤ 0.001, two-tailed t-test.

Figure 5.. Predictive biomarkers for response to EGFR inhibitors in EGFR altered GBM PDCs.

(a) Mutational landscape of EGFR alterations in GBM cohort. (b) For the 10 EGFR inhibitors, top drug–feature associations identified by dNetFS are plotted for their frequency and effect size (n=49 biologically independent samples in (a). Node size is proportional to the single drug–feature linear correlation. (c) Gene expression profiles of EGFR and NRG1, and AUC drug response profiles of erlotinib and dacomitinib, over the topological representation (n=44 biologically independent samples). (d) Drug response assessment of EGFR inhibitors with shRNA-mediated knockdown of NRG1 or NT (non-target). Cell viability for each dose was normalized to shNRG1 or NT transduced cells only. (n=10 independent experiments with 3 technical replicates) (e) Comparisons of cancer pathway activities between two EGFR-altered GBM subgroups that were most sensitive and most resistance. We adopt the enrichment score derived from ssGSEA analysis as assessment. (f) Drug response assessment of EGFR inhibitors with PI3K-AKT-mTOR (PAM) inhibitors or DMSO. Cell viability for each dose was normalized to PAM or DMSO treated cells only. Mean AUC value for 4 PAM inhibitors (BYL719, BKM120, BEZ235, and AZD2014) was plotted (n=8 independent experiments with 3 technical replicates). P-values: c, pearson correlation test. P-values: d,f, two-sided wilcox rank sum test. Horizontal lines within the violin plot represent 0.25, 0.50, and 0.75 quantiles.

Figure 6.. Clinical feasibility of PDC drug screening-guided precision oncology.

(a) Bar graph represents normalized AUCs (Z-scores from pan-cancer dataset) of indicated drugs (n=31 biologically independent samples). Clinical responses were determined according to RECIST. Multiple-target drugs are classified based on corresponding representative genomic targets. (b) Representative box plot of (a). CR (complete response); PR (partial response); SD (stable disease); PD (progressive disease) per RECIST criteria. Box plot spans from the first to third quartiles and the whiskers represent the 1.5 interquartile range. (c) Receiver operating characteristic (ROC) curve was plotted by the sensitivity (%) and 100-specificity (%) values for predicting clinical response rate using z-scores in (a). (d-f) T1-weighted contrast enhanced (T1CE), T2-flare magnetic resonance images (MRI) or computed tomography (CT) images for indicated patients before and after drug treatment were demonstrated. Circles on the bar represent days obtaining presented images, where 0D refers the day starting drug treatment. Negative (−) or positive (+) D represents days before or after treatment, respectively. Red arrows indicate measurable or progressed tumors, and orange arrows represent partial response. Vertical scattered plot for AUCs of indicated drugs in pan-cancer AUC reference dataset for the indicated drug. AUCs of PDCs isolated from the illustrated patients were highlighted. P-value: b, two-sided wilcox rank sum test. P-values: c, two-sided binomial exact test.

Figure 7.. Schematic illustration of major lineage-specific and genomic associated drug interactions.

Molecular targeting agents are clustered based on drug family classification and connected to various genomic alterations including mutation, mRNA expression, copy number variation (CNV), and fusion. Edges represent sensitive or resistant gene-drug interactions. Lineage-specific drug associations are highlighted in each drug node based on specific lineage type.