A Library of Phosphoproteomic and Chromatin Signatures for Characterizing Cellular Responses to Drug Perturbations

Abstract

Although the value of proteomics has been demonstrated, cost and scale are typically prohibitive, and gene expression profiling remains dominant for characterizing cellular responses to perturbations. However, high-throughput sentinel assays provide an opportunity for proteomics to contribute at a meaningful scale. We present a systematic library resource (90 drugs × 6 cell lines) of proteomic signatures that measure changes in the reduced-representation phosphoproteome (P100) and changes in epigenetic marks on histones (GCP). A majority of these drugs elicited reproducible signatures, but notable cell line- and assay-specific differences were observed. Using the "connectivity" framework, we compared signatures across cell types and integrated data across assays, including a transcriptional assay (L1000). Consistent connectivity among cell types revealed cellular responses that transcended lineage, and consistent connectivity among assays revealed unexpected associations between drugs. We further leveraged the resource against public data to formulate hypotheses for treatment of multiple myeloma and acute lymphocytic leukemia. This resource is publicly available at https://clue.io/proteomics.

Keywords: GCP; L1000; LINCS project; P100; drug discovery; epigenetics; mass spectrometry; mechanism of action; proteomics; signaling.

Figure 1. Analyte-Centric Analysis of the Resource

(A) Experiment schematic. Cells were treated with one of 90 small-molecule perturbagens with a minimum of three biological replicates. After 3 (P100) or 24 (GCP) hr of exposure to treatment, cells were lysed and profiled in each of the two assays. After quality control filtering, over 3,400 individual profiles constituted the resource. (B) JNK inhibitors downregulate phosphorylation of S100 on JUND. Observed abundance of JUND, S100 phosphorylation across all P100 samples (level 4 data) is shown as a single horizontal row of a heatmap. JNK inhibitor profiles are marked with black ticks along the top of the heatmap. Computing enrichment of the JNK inhibitors across the 1,683 samples with data present yields a GSEA FDR result of 0.00. (C) MEK inhibitors downregulate phosphorylation of S3426 on AHNAK. Computing enrichment of the MEK inhibitor profiles across 1,390 P100 samples with data present yields an FDR of 0.00. (D) Decitabine increases ubiquitination of lysine 18 on histone H3. Observed abundance of H3K18ub1K23ac0 across all GCP samples (level 4 data) is shown as a single horizontal row of a heatmap. Computing enrichment of decitabine profiles across 1,159 GCP samples with data present yields an FDR of 0.00. (E) Number of drug-cell type conditions with at least one “dominant” analyte. (F) Number of analytes dominant in at least one condition. (G) Percentage of signal derived from dominant analytes (average of all profiles). (H) Number of analyte-drug combinations that transcend cell type.

Figure 2. Connectivity Framework

(A) Each sample is represented as a profile of analyte measurements. Spearman correlations are computed between all profiles within a cell line. Finally, we compute connectivity scores by comparing the observed correlations with a background of correlations. Computing connectivity collapses replicates. Connectivity maps may be represented as matrices or networks. Different color scales are used for profile and similarity matrices to make them more distinguishable. (B) Vorinostat versus belinostat: example of a positive connectivity score close to 1. The background distribution (blue) consists of the correlations between the replicates of belinostat and all other samples. The test distribution (green) consists of the correlations between the replicates of belinostat and the replicates of vorinostat. (C) Vorinostat versus DMSO: example of a connectivity score close to 0. The background distribution (blue) consists of the correlations between the replicates of DMSO and all other samples. The test distribution (green) consists of the correlations between the replicates of DMSO and the replicates of vorinostat. (D) Vorinostat versus EX527: example of a negative connectivity score close to −1. The background distribution (blue) consists of the correlations between the replicates of EX527 and all other samples. The test distribution (green) consists of the correlations between the replicates of EX527 and the replicates of vorinostat. (B) to (D) show data for the GCP assay in A375 cells. See also Figure S1.

Figure 3. Majority of Compounds Are Reproducible in Both Assays, but Assays Show Different Sensitivities to MoA Classes

(A) Distributions of all Spearman correlations among replicates (green) and among non-replicates (blue). (B) Bar chart showing the number of compounds considered reproducible in each cell line-assay combination. A compound was considered reproducible if its replicates were significantly correlated compared with a permutation null (q value < 0.05). The algorithm was rerun 10 times to generate error bars; the center is the median, error bars represent 95% confidence interval of 1,000 bootstrapped iterations. The shaded component indicates the overlap of reproducible compounds between GCP and P100. (C) Bar charts showing the median connectivity of compounds annotated with the same mechanism of action (MoA). The center is the median, error bars represent 50% confidence interval of 1,000 bootstrapped iterations. (D) Heatmaps showing connectivity among compounds belonging to the BRD inhibitor, sirtuin modulator, and JNK inhibitor MoA classes. The BRD inhibitor class shows high connectivity in both P100 and GCP; the sirtuin modulator class shows high connectivity in P100 but low connectivity in GCP; and the JNK inhibitor class shows low connectivity in both P100 and GCP. Each square is the median of six cell-specific connectivity scores. The labels of the matrices are symmetric; that is, the columns (left to right) have the same annotations as the rows (top to bottom). Color scale applies to both (C) and (D). See also Figures S2 and S3.

Figure 4. Connectivity Profile Analysis

(A) t-SNE projection of P100 and GCP connectivity profiles for connections within each cell type (distance metric = Pearson correlation, perplexity = 60, learning rate = 10). (B) Schematic representation of cutting a dendrogram at a fixed percentage of its height and counting resulting clusters, for illustration only. (C) Number of connectivity clusters formed as a result of cutting dendrograms as depicted in (B). Individual data points (six per assay) are overlaid on the box plots and jittered on the y axis for clarity. The center is the median, error bars represent the 25th and 75th percentiles. (D) Connectivity flows from P100 connectivity clusters to GCP connectivity clusters. For each cell line, only compounds reproducible in both assays are included. Each cluster is annotated by the major pertinent mechanistic classes for each assay, with the number of drugs in each class shown in parentheses. Colors are arbitrary. Because the analysis was restricted to the reproducible compounds in each cell type and single member clusters were eliminated, the number of clusters at the 60% cut for P100 and GCP may differ slightly from (C). See also Figures S4 and S5.

Figure 5. Comparison with Transcriptomic Data Demonstrates Assay-Specific Sensitivities

(A) Schematic of the comparison of connectivity matrices in three assays, including L1000 transcriptomic data. (B) Percent overlap of the top 5% of connections in the P100, GCP, and L1000 assays. The light-gray bars show percent overlap for cell-specific connectivities, and the dark-gray bars show percent overlap for aggregated connectivities. The dashed line indicates the percent overlap expected by chance. The center is the median, error bars represent 95% confidence interval of 1,000 bootstrapped iterations. (C) Recall of connectivity profiles across assays. The y axis indicates the percent of compounds (n = 90) that have recall greater than 0.95 for a pairwise comparison. Recall of 0.95 means that the connectivity profile for a particular compound in one assay had higher similarity to its corresponding connectivity profile in another assay than to 95% of other connectivity profiles (see Figure S7E for a schematic of this algorithm). The center is the median, error bars represent 95% confidence interval of 1,000 bootstrapped iterations. Shading as in (B). (D) Network views of the top 0.5% of connections in each assay. All connectivity scores are positive. Compounds are represented by nodes, and MoA is encoded by the color of the node. See also Figure S7.

Figure 6. Multi-Assay Data Integration Reveals Cell-Specific Vulnerabilities

(A) Schematic of the integration of connectivity matrices to create AVG data. (B) Network view of the top 0.5%of connections for AVG data, which is an average of the connectivity scores in P100, GCP, and L1000. All connectivity scores are positive. Compounds are represented by nodes, and MoA is encoded by the color of the node. TG101348 (circled) has unexpected connectivity to the BRD inhibitors. (C) Heatmap view of the connectivity scores between vemurafenib and the two MEK inhibitors. The connectivity scores in A375 (0.91 and 0.93) are considerably higher than connectivity scores in any other cell line. (D) Heatmap view of the connectivity scores between SCH 900776 and the two JNK inhibitors. The connectivity score between SCH 900776 and SP600125 in MCF7 (0.88) is an unexpected cell-specific connection. (E) Results of a 5-day follow-up viability experiment. The y axis shows GR values in A375 (red rectangles) and MCF7 (green diamonds) for SP600125 (left), SCH 900776 (middle), and CC-401 (right).GR values quantify drug cytotoxicity and are insensitive to different cell growth rates. The x axis shows drug concentration on a log₁₀ scale. See also Supplemental Information.

Figure 7. Connectivity Query and Perturbation Set Analysis of a Diverse Set of Cancer Lineages Validates Genetics and Identifies Potential Therapeutic Avenues

(A) Connectivity query of chromatin signatures from EZH2 loss-of-function cell lines from the CCLE. Results are sorted by the median connectivity to the perturbation across the five EZH2 loss-of-function cell lines. (B) Adaptation of the GSEA algorithm to test for enrichment of MoA classes in connectivity results. The top ranked set is EZH2 inhibitors; all hits to this set are at the top of the list sorted by average connectivity. (C) Stratification of two sets of NSD2 gain-of-function classes via hierarchical clustering of connectivity query results. (D) The most highly anti-connected perturbations to the t4;14 and NSD2:p.Glu1099Lys gain-of-function classes of cell lines, when ranked by connectivity for each class. Enriched perturbation sets with FDR of <0.05 are shown for each class. The HDAC inhibitors are the most anti-connected perturbations to the t4;14 subtype while BRD, CDK, and mTOR inhibitors are all anti-connected to the NSD2:p.Glu1099Lys subtype. See also Figure S8.