Computational pipeline predicting cell death suppressors as targets for cancer therapy

Summary Identification of promising targets for cancer therapy is a global effort in precision medicine. Here, we describe a computational pipeline integrating transcriptomic and vulnerability responses to cell-death inducing drugs, to predict cell-death suppressors as candidate targets for cancer therapy. The prediction is based on two modules; the transcriptomic similarity module to identify genes whose targeting results in similar transcriptomic responses of the death-inducing drugs, and the correlation module to identify candidate genes whose expression correlates to the vulnerability of cancer cells to the same death-inducers. The combined predictors of these two modules were integrated into a single metric. As a proof-of-concept, we selected ferroptosis inducers as death-inducing drugs in triple negative breast cancer. The pipeline reliably predicted candidate genes as ferroptosis suppressors, as validated by computational methods and cellular assays. The described pipeline might be used to identify repressors of various cell-death pathways as potential therapeutic targets for different cancer types.


Figure S1. Supporting information for extraction of correlation scores (related to Figure 2) (A)
The indicated breast cancer lines were seeded in 96 well plates (maintaining a confluency of ~50%).24 hours later, the indicated FIN (erastin, RSL3, FIN56) were added in different concentrations.72 hours later MTT was performed to assess the level of cell viability.Plots show means ± sd of 3-4 independent repeats.(B) Gene ontology and pathway enrichment for the 41 genes with the highest correlation between gene expression and FIN AUCs.was performed, and the top 75 genes closest to the GCH1 node were selected for further analysis as described in Fig. 5D-E.In brief, all the genes were ranked by the correlations of their expression levels to the dependency score of GPX4 or to the AUCs of 4 FINs taken from the CTRP dataset.The enrichment of the 75 GCH1-neighboring genes in each combination was measured by GSEA against all the genes ranked by those correlations, and the normalized enrichment scores (NES) are shown in this plot.The red dot indicates the actual set of predictors used for the UMAP shown in Figure 5B (note, the NES values in Fig. 5D-E are a bit different compared to the NES of the red dot in this analysis, due to the number of signatures examined in the GSEA in each case).(C) Zoom-in of the GPX4 node in the UMAP from Fig. 5A.(D) The genes in the GPX4 neighborhood were ranked by their proximity to GPX4.PubMed search of each gene with 6 ferroptosis related terms was performed.A red labeled box indicates a citation for the gene with the indicated term in PubMed.(E) The UMAP projection shown in Figure 5A was recreated 500 times, each time with a different random seed for the UMAP creation (the seeds themselves were also randomly selected).For each UMAP, the distance between the indicated genes (which were selected for experimental validation) and the GCH1 node (for most genes) or the GPX4 node (for CTBP1 and PDAP1) was measured and ranked among other genes.(F) Correlation between the distance on the UMAP projection compared to the distance of the predictor matrix.Each point represents a gene.The xaxis shows the Euclidean distance of each gene from GCH1 in the UMAP shown in Fig. 5A.The yaxis shows the Euclidean distance between each gene and GCH1, calculated using the normalized values of the 13 predictors themselves.The blue line represents the linear correlation for the entire data; the red line shows the linear correlation for the 200 genes closest to GCH1.The red points indicate the genes from the GCH1 local neighborhood that were chosen for experimental validation.

Figure
Figure S2.Supporting information for the ferroptosis response signature (related to Figure 3) (A) UpSet plot showing the number of significant (p-value < 0.05) differentially expressed genes (DEGs, top bars) for any intersection of treatments (bottom plot).In the bottom plot, the blue and red points indicate basal and mesenchymal cell lines, respectively, and the gray and white shaded lines indicate treatment with erastin and RSL3, respectively.The annotated genes are ferroptosis to apoptosis biomarkers found in our previous publication.(B) Summary of the number of DEGs per number of treatments in each intersection.(C) Gene set enrichment analysis was performed for the 10 treatments.The pathways with the highest average normalized enrichment scores (NES)

Figure
Figure S3.Calculating the ranked transcriptomic similarity scores per gene (related to Figure 4) (A) Similarity scores for each gene and each treatment were calculated across the cell lines available from the CMAP dataset.This representative figure shows the scores for the upregulated set of 20 genes, extracted from the signatures of GPX4 and GCH1, as an example.In each panel, the red dot indicates the maximum value.The similarity score for each gene per treatment were taken to be this maximum value.(B) The maximum values from A were ranked.Representative plot showing the ranks of the scores are taken from the "BT549_E" treatment.(C) The final transcriptomic similarity scores for each gene were calculated by subtracting the downregulated gene set score (shown in B, right panel) from the upregulated gene set scores (shown in B, left panel).The red dots indicate the enrichment of 40 known ferroptosis suppressors, measured by

Figure S4 .
Figure S4.Integrating the correlation and transcriptomic similarity predictors for predicting ferroptosis targets (related to Figure 5) (A) Gene ontology (biological processes [top], molecular functions [bottom]) was performed for the GCH1 neighborhood genes shown in Fig. 5B.Only significant GO terms are shown.(B) 1023 combinations of predictors (see methods) were observed.For each combination, UMAP projection

Figure S5 .
Figure S5.Supporting information for the integration of predictors and validation (related to Figure 6) (A) Essentiality scores for the genes selected for experimental validation, across all cell lines of all lineages included in the Achilles dataset.(B) Essentiality scores for the genes selected for experimental validation, across all cell lines of all lineages included in the Achilles dataset.Boxplot depicts the essentiality scores per lineage.(C) Knock-down of the indicated genes was performed in MDA-MB-468 and HCC70 cells using shRNA.qRT-PCR analysis was performed to validate the knock-down efficiency of the shRNA.Numbers indicate fold change in the expression of the gene being knocked-down.(D) Gene ontology (molecular functions) was performed for the 8 genes used for validation.Only significantly enriched terms are shown.