A multilineage screen identifies actionable synthetic lethal interactions in human cancers

Abstract

Cancers are driven by alterations in diverse genes, creating dependencies that can be therapeutically targeted. However, many genetic dependencies have proven inconsistent across tumors. Here we describe SCHEMATIC, a strategy to identify a core network of highly penetrant, actionable genetic interactions. First, fundamental cellular processes are perturbed by systematic combinatorial knockouts across tumor lineages, identifying 1,805 synthetic lethal interactions (95% unreported). Interactions are then analyzed by hierarchical pooling, revealing that half segregate reliably by tissue type or biomarker status (51%) and a substantial minority are penetrant across lineages (34%). Interactions converge on 49 multigene systems, including MAPK signaling and BAF transcriptional regulatory complexes, which become essential on disruption of polymerases. Some 266 interactions translate to robust biomarkers of drug sensitivity, including frequent genetic alterations in the KDM5C/6A histone demethylases, which sensitize to inhibition of TIPARP (PARP7). SCHEMATIC offers a context-aware, data-driven approach to match genetic alterations to targeted therapies.

Extended Data Fig. 1 |. Combinatorial gRNA library design.

(a) The dual CRISPR library targets all pairs of 67 by 176 genes across 7 cell lines. Each gene is targeted by 3 guide RNAs resulting in 9 guide pairs for each gene pair assayed in 2 replicates. (b) Design of the custom 130 base pair oligonucleotide pool used to construct the combinatorial CRISPR library. sgRNA1 and sgRNA2 can target the same gene or two different genes. hU6, human U6 promoter; sgRNA, single-guide RNA; BsmBI, BsmBI restriction enzyme recognition site. (c) Two-step cloning strategy to package the oligonucleotides in panel B into a functional combinatorial CRISPR library. mU6, murine U6 promoter.

Extended Data Fig. 2 |. Reproducibility and validation of fitness measurements.

(a–c) Scatter plots showing reproducibility across replicates of fitness measurements at the level of (a) individual guide pairs; (b) gene pairs, median over all relevant guide pairs; and (c) genes, integrating over all relevant pairwise fitnesses involving each gene. Contour lines are drawn at 5, 10, 25, 50, 75, 99 percentiles. All measurements beyond the 99 percentile are represented explicitly as black points. Note progressive increases in reproducibility (Pearson correlation $r$) with increasing integration of data. (d) Bar plot of the Pearson correlation between replicate guide pair fitness measurements (teal) or replicate gene pair fitness measurements (blue) in each of the 7 cell lines. (e) Bar plot of Pearson correlation between replicate single-gene fitness measurements from the human U6 position (hU6, blue) or the murine U6 (mU6, red) position. (f) Bar plot of the Pearson correlation between the single-gene fitness measurements from this study (hU6 position, blue; mU6 position, red) and the single-gene fitness measurements from the DepMap project. (g) Recovery of common-essential genes annotated by DepMap (area under the receiver operating characteristic curve, auROC) when scoring essential genes de novo based on sgRNAs expressed by the human U6 (hU6, blue) or murine U6 (mU6, red) promoters. (h) Scatterplot of single-gene knockout fitness measurements scored in this study versus those measured by DepMap, including data from all seven cell-line contexts. (i) Contingency table of essential genes classified by DepMap versus this study. OR, odds ratio. P-value by Fisher’s exact test. (j) Distributions of single-gene fitness measured in this study, split by DepMap status (DepMap common-essential genes in orange, non-essential genes in blue). Distributions are shown for 6/7 cell lines, as the seventh line (MCF10A) is not characterized by DepMap.

Extended Data Fig. 3 |. Reproducibility and validation of genetic interaction measurements.

(a) Decrease in fitness due to single-gene disruptions to BRCA1 or PARP1, alongside double-gene disruption to BRCA1 and PARP1. Fitnesses tracked over the course of 21 days. Error bars denote the standard error of the mean (n = 9 sgRNA pairs). (b) Volcano plot showing false discovery rate versus genetic interaction score for all gene pairs characterized in CAL27 cell line. The confidence interval contains 95% of all scores where at one sgRNA targets the AAVS1 safe-harbor locus (adeno-associated virus integration site 1). Point color shows the absolute fitness score of each gene pair. (c) Distribution of coefficients of variation (CV) of top 1000 synthetic-essential gene pairs in individual cell lines (blue) versus the corresponding distribution where the guide-pair to gene-pair mapping is randomized (orange). Dotted line shows the threshold CV that best separates the two distributions, used in the Methods to triage interactions into single cell line versus multi-lineage classifications. (d) Bar plot of the Pearson correlation between replicate genetic interaction measurements in each of the 7 cell lines. Pearson correlations are also shown after pooling measurements within each of the 3 tissues or across all tissues (multi-lineage). All measurements shown in teal, significant interactions (FDR < 30%) shown in blue. FDR, False Discovery Rate. (e) Number of significant genetic interactions in each of the 7 cell lines, in each of the 3 tissue pools, or when pooling all contexts as multi-lineage. FDR < 1% in teal, FDR < 10% in blue.

Extended Data Fig. 4 |. Heatmap of essential genes identified in this study.

Blue color indicates a human gene (columns) scoring as essential in a tumor cell line (bottom-most rows), tissue or subtype context (middle rows), or multi-lineage (top row).

Extended Data Fig. 5 |. Mapping essential systems with combinatorial CRISPR data.

(a) Same map of multi-gene systems as Fig. 3a, visualized as a vertical hierarchy. Nodes represent systems, arrows represent containment of one system within a larger one. System size (number of proteins) shown by node size. Individual proteins not shown. (b) Tests for identifying essential systems by independent gene lethality (IE), synthetic lethality within systems (SE_internal), or synthetic lethality with an external gene (SE_external). Circle nodes represent systems; diamond nodes represent genes; arrows linking one circle to another indicate hierarchical containment of the first system (child) by the second (parent). Color represents viable (gray) versus lethal (red) status of the corresponding gene or pairwise gene knockout. (c) Fraction of systems (y-axis) scoring as essential in each of four databases of subcellular systems (x-axis) revealed by the different essentiality rubrics (bar colors). Error bars, 95% confidence intervals of the sampling proportion.

Extended Data Fig. 6 |. Heatmap of all systems identified by across-system essentiality.

Each colored box indicates a system (rows) that scored as essential conditioned on knockout of an independent gene outside the system (columns). Colors denote the relevant context (tissue type or multi-lineage). Abbreviated version shown in Fig. 3b.

Extended Data Fig. 7 |. Example synthetic-essential interactions prioritized by GDSC or DepMap.

(a–d) Dose response curves of AKT1 inhibitors GSK690693 (a), afuresertib (b), uprosertib (c), and ipatasertib (d). Solid black curve shows median dose response of all cell lines not disrupted in CHEK1 (GSK690693 n = 892, afuresertib n = 936, uprosertib n = 943, ipatasertib n = 936). For each drug, three example dose response curves of CHEK1-mutant cell lines are shown in red. (e) Genetic interaction scores measured in experiments disrupting MAP2K1, MRE11, or both genes in combination. Relevant interaction measurements for all cell lines, replicates, and time points are shown (points, n = 47). Interaction measurements for single genes are measured in the combinatorial gRNA format (see fig. S1) by targeting the AAVS safe-harbor locus in addition to the gene. Box plots depict the median and the bounds of the box depict the first and third quartiles. The whiskers of the box plot mark the minima and maxima, but no further than 1.5 times the interquartile range. Significance of difference in means tested by Mann-Whitney’s U test (*, P = 8.3 × 10⁻⁹ and P = 4.1 × 10⁻¹³ for comparing MAP2K1-AAVS and AAVS-MRE11 to MRE11-MAP2K1 respectively). (f) Violin plots showing the distribution of fitness changes in tumor cell lines due to knockout of a homologous recombination gene (x-axis). $n$ represents the number of MAP2K1-disrupted or MAP2K1 wild-type cell lines profiled with CRISPR knockouts in DepMap. Blue and yellow distributions group cell lines by MAP2K1 mutation status. Box plots and asterisks follow convention of panel e.

Extended Data Fig. 8 |. Translation rate of synthetic-lethal interactions based on context.

The y-axis shows the fraction of interactions that translate to predictions of cell-line sensitivity in the DepMap screen of genome-wide CRISPR gene knockouts × cell lines^,. The x-axis considers this fraction for each of the three tissue-specific networks and the pan-essential network. A gene-gene interaction is tested in DepMap if at least 25 cell lines have an alteration in one of the two genes. Suggestive interactions refer to those at a threshold of P < 0.05 by Student’s two-sided $t$-test without multiple hypothesis correction. Stringent interactions refer to a threshold of FDR < 30%. NS = No sample interactions met the threshold.

Extended Data Fig. 9 |. Supplemental exploration of PARP7 inhibition.

(a) Violin plots of PARP7 expression (left) and fitness effects of PARP7 CRISPR knockout (right) as measured in DepMap. Both show significant differences between KDM6A loss and wildtype status (P = 0.034 and 1.4 × 10⁻⁴ respectively) by two-sided Mann-Whitney U test. Above each plot are values of $n$ providing the number of KDM6A-disrupted or KDM6A wild-type cell lines. Box plots depict the median, and the bounds of the box depict the first and third quartiles. Distributions truncated to show middle 90% of data. (b) Differences in mRNA expression associated with KDM5C loss, shown for interferon stimulated genes (left, ISGs) or PARPs (right). ISGs are significantly upregulated as a set (P = 5.4 × 10⁻⁷ by Mann-Whitney $U$ test versus a background of all human genes). (c) As for panel a, but stratifying cell lines on the genetic alteration status of the genome-integrity cluster (loss of any genes in this cluster). Both PARP7 expression and fitness of PARP7 knockout are significantly affected by loss of genes in this cluster (P = 6.4 × 10⁻⁵ and 1.6 × 10⁻⁶ respectively). (d) Violin plots showing PARP7-inhibitor dose response (area under curve), which is significantly stratified by KDM6A (left, P = 1.5 × 10⁻³) or genome-integrity cluster (right, P = 8.2 × 10⁻⁴) alteration status. Significance determined by two-sided Mann-Whitney U test. Values of $n$ represent the number of cell lines with and without alterations in the KDM6A (left) or genome-integrity cluster (right). Boxplots follow the same convention as panel A. (e–g) Dose-response curves from our PARP7 inhibitor PRISM screen, providing specific examples from tumor lines derived from pancreas (e), urinary tract (f), or colon (g). Median dose response curves for cell lines without alterations in the genome integrity cluster appear in solid black (error bars, standard errors of the mean; Pancreas n = 16, Urinary tract n = 22, Colorectal n = 30). Dose response is represented by the log₂ fold change in cell abundances between PARP7 inhibitor treatment to DMSO control. Two example dose response curves from cell lines with alterations in this cluster are depicted in red. h, PARP7i dose responses for KDM5C unaltered (blue) or KDM5 loss (red) tumor cell lines, including all colorectal, urinary tract, or pancreas cell lines on the PRISM panel. Area under the dose response curve measurements were normalized to the median of those measured in unaltered cells. *, significant difference determined by two-sided Mann-Whitney U test. Box plots depict the median and the bounds of the box depict the first and third quartiles. The whiskers of the box plot mark the minima and maxima, but no further than 1.5 times the interquartile range. Measurements beyond the whiskers are not shown.

Fig. 1 |. Strategy for discovery of actionable synthetic lethal interactions.

a, Core molecular systems of tumor cells interrogated by systematic disruption of combinations of genes across a panel of tumor-type contexts. hU6, mU6, human- and murine-U6 promoters. b, Resulting tumor cell fitness profiles analyzed to identify robust synthetic lethal interactions based on: repeated observations across tissue and cell-type contexts (penetrant interactions, top), convergence on systems (middle) and ability to predict sensitivity to targeted treatments (bottom). Drug concn, drug concentration. c, Circle-packing diagram of the NeST map of human subcellular systems, filtered to the 64 systems covered in the present study. Subcellular systems are denoted as circles; containment of one circle in another denotes a system that is a subcomponent of a larger one. The circle color denotes the fraction of genes within the system represented in the combinatorial CRISPR library, according to the color scale defined in d. EMT, epithelial–mesenchymal transition. d, Histogram showing the combinatorial CRISPR library coverage of NeST systems included in the present study. NeST systems are binned by the fraction of their genes represented by gRNAs in the CRISPR library (coverage, $x$ axis). Bar shading increases with coverage. e, Points showing all pairwise gene combinations with MSH2, with the (MSH2 × gene) double-mutant fitness plotted versus the single-mutant fitness of each gene ($y$ versus $x$ axis). The diagonal shows the least-squares-fit regression line by which a gene is determined to have a positive (above line, for example, KAT5) or negative (synthetic essential, below line, for example, BRD4) interaction with MSH2.

Fig. 2 |. Multilineage and context-specific mapping of synthetic essentiality.

a, Typing the essentiality of genes and pairwise genetic interactions. Scoring occurs first across all contexts to identify multilineage essentialities (red), then within tissue or biomarker contexts (purple) and then within individual cell lines (light blue). b,c, Piecharts showing proportions of essential genes (b) and synthetic essential interactions (c) by scoring context (colors same as a, counts in parentheses). The insets at the upper right show the overall fraction of genes scoring as essential (b) versus the fraction of gene pairs scoring as synthetic essential (c) in any context. ‘Multilineage’ denotes interaction being significant when scoring across all lineages, whereas ‘several contexts’ denotes significance in more than one context but not all. d, Heatmap of strongest synthetic essential genetic interactions based on consistent discovery across contexts (most extreme negative multilineage scores with all interactions having FDR < 0.1). Columns show interacting gene pair and row modes of interaction scoring based on (top to bottom) multilineage, tissue specific or individual cell-line analysis as per a. The blue–black–yellow color gradient represents the full range of negative–zero–positive scores. The same color gradient appears to be used for d–f. Gene pairs in red are mentioned in the text. e, Heatmap of strongest synthetic essential genetic interactions identified in specific tissue contexts (most extreme negative interactions by tissue score among interactions failing the multilineage test; all interactions with FDR < 0.1). The display is as per d. f, Heatmap of representative synthetic essential genetic interactions that are conditional on a specific biomarker (top four rows). The display is as per d. Activating gain-of-function (GOF) mutations: KRAS, PIK3CA. Loss-of-function (LOF) tumor-suppressor mutations: POLQ. Interactions dependent on TP53 are active under the TP53 wild-type status.

Fig. 3 |. Structural map of essential multigene systems.

a, Multi-scale map of 64 tumor subcellular systems, represented as a kaleidoscopic nested circle layout as per Fig. 1c. The color indicates whether a system (circle) or gene (diamond) is essential across cancer types (red), essential in specific tissue or biomarker contexts (green) or nonessential (blue). Essential genes are defined by BAGEL (BF > 5; Methods) and essential systems are determined by a one-sided hypergeometric test (FDR < 30%; Methods). Four systems are expanded at the right to show the underlying genetic data, with accompanying barplots providing odd ratios (ORs) of enrichment for IE, SE_internal and SE_external where relevant. The highlighted systems provide examples of all three effects: IE (chromosome and HR systems); SE_internal (mitosis, HR); and SE_external (G1 checkpoint). Red edges denote synthetic essential genes identified with combinatorial CRISPR. Gray edges denote protein–protein interactions. EGFR, Epidermal growth factor receptor. An alternative hierarchical view of this map is provided in Extended Data Fig. 5a. b, For systems identified by SE_external, the heatmap shows the outside gene dependencies (columns) and the contexts in which dependency is observed (colors). System names (rows) are listed next to selected interacting genes in each system. The figure shows a subset of tissue-specific and multilineage system–gene interactions considering systems catalogued by NeST and Reactome; for the full set see Extended Data Fig. 6 and Supplementary Tables 6–9. AR, androgen receptor; CDK, cyclin-dependent kinase; HDR, homology-directed repair; NER, nucelotide-excision repair; RTK, receptor tyrosine kinase.

Fig. 4 |. Prioritization of interactions that predict sensitivity to targeted therapeutic agents.

a, Chord diagram of multilineage interactions that link genes impacted by frequent somatic mutations (blue) to genes encoding druggable targets (green). Some genes have both properties (purple). b, Gene–gene interactions identified with combinatorial CRISPR (left) examined in complementary screens measuring dependency of mutated cell lines on targeted drugs (GDSC) or CRISPR sgRNA knockdowns (DepMap). Interactions are stratified into four categories (pie chart with colored slices, right). The significance was tested using two-sided Mann–Whitney $U$-test. Suggestive: P < 0.05; stringent: P < 0.05 and FDR < 30% (multiple hypothesis correction using the Benjamini–Hochberg procedure). c, Genetic interaction scores in the CRISPR experiments disrupting AKT1, CHEK1 or both genes in combination. Relevant interaction measurements for all cell lines, replicates and timepoints are shown (points, n = 47). Interaction measurements for single genes are measured in the combinatorial sgRNA format by targeting the gene in combination with the AAVS safe-harbor locus. Box plots depict the median and the bounds of the box depict the first and third quartiles. The whiskers of the box plot mark the minima and maxima, but no further than 1.5× the interquartile range. *P = 1.8 × 10⁻³ (comparing AKT1-AAVS with CHEK1-AKT’) and P = 1.5 × 10⁻² (comparing CHEK1-AAVS with CHEK1-AKT1) using two-sided Mann–Whitney $U$-test. d, Violin plots showing normalized IC₅₀ values of CHEK1-disrupted (Mut) or CHEK1 wild-type (WT) cells under separate exposures to a panel of AKT1 inhibitors or cisplatin as negative (neg.) control. Values of $n$ represent the number of CHEK1-disrupted or CHEK1 WT cell lines profiled with AKT1 inhibitors in GDSC. IC₅₀ values normalized to the median of CHEK1 WT cells. Box plots and asterisks follow the convention in c. Significance was tested using two-sided Mann–Whitney $U$-test: GSK690693 P = 0.02, afuresertib P = 5.8 × 10⁻³, AZD5363 P = 0.03, uprosertib P = 1.9 × 10⁻³, ipatasertib P = 6.4 × 10⁻³ and cisplatin P = 0.59. NS, Not significant. e, Dependency of drug response ($y$ axis) on genetic alteration of each gene targeted in the present study ($x$ axis, in rank order). The first five panels show response to AKT1 inhibitors; the last panel (far right) shows cisplatin response as a negative control. Dependency score of a gene is defined by the difference in median dose response (area under the curve (AUC)) of cells with genetic alterations in that gene compared with the median dose response of cells without that gene altered. The dependency on CHEK1 alteration is highlighted in red.

Fig. 5 |. Characterizing drug–mutation interactions that distinguish PARP family inhibitors.

a, Heatmap of multilineage, synthetic essential interactions involving PARP7, FEN1, PARG and PARP1. The blue–black–yellow color gradient represents the full range of negative–zero–positive genetic interaction scores. The gray rectangles denote clusters of genes with similar interaction patterns (hierarchical clustering of Pearson’s correlation; Methods) which are significantly enriched for genes in distinct NeST systems (colored text, hypergeometric test P < 0.02). b, Violin plots of fitness effects of PARP7 CRISPR knockout as measured in DepMap, stratified on KDM5C alteration status (red, copy number loss; blue, wild-type). Values of $n$ represent the number of KDM5C-disrupted or KDM5C wild-type cell lines profiled with PARP7 CRISPR knockout in DepMap. Box plots depict the median and the bounds of the box the first and third quartiles. Distributions are truncated to show the middle 90% of data. *P = 3.8 × 10⁻⁷ by two-sided Mann–Whitney $U$-test. c, Violin plots of PARP7 mRNA expression levels (log₂(transcripts per million)) for cell lines in the Cancer Cell Line Encyclopedia, stratified by KDM5C alteration status. Values of $n$ represent the number of KDM5C-disrupted or KDM5C wild-type cell lines profiled where PARP7 expression is measured. The display follows the convention of b. *P = 5.2 × 10⁻⁶ by two-sided Mann–Whitney $U$-test. d, Genes in the genome integrity cluster (points) each used as a biomarker to stratify DepMap cell lines into genetically altered cases (mutations, copy number aberrations) or unaltered controls. Case versus control groups are compared for differences in fitness under PARP7 CRISPR knockout ($y$ axis) and differences in PARP7 expression ($x$ axis). e, Relative mRNA expression level of KDM5C (left) or PARP7 (right) under KDM5C siRNA knockdown (black) or non-targeting control (NTC) siRNA (gray). Expression levels normalized to NTC. The horizontal line shows mean and error bars denote 95% confidence interval (CI). Points (n = 6, 2 biological replicates and 3 technical replicates each) denote replicate measurements. *P = 1.7 × 10⁻⁵ by two-sided Student’s $t$-test. f, Schema for profiling effects of PARP7 inhibition (RBN-2397) across the PRISM-barcoded pool of tumor cell lines. Cells are treated with each of seven different doses of inhibitor as well as DMSO vehicle control. Relative abundance of cells was determined by next-generation sequencing. g, Violin plots showing PARP7 inhibitor dose–response area under the curve (AUC) stratified by KDM5C alteration status. Values of $n$ represent the number of KDM5C-disrupted or KDM5C wild-type cell lines profiled with RBN-2397. The display follows the convention of b. *P = 1.1 × 10⁻⁴ by two-sided Mann–Whitney $U$-test.