Metabolic modelling-based in silico drug target prediction identifies six novel repurposable drugs for melanoma

Despite high initial response rates to targeted kinase inhibitors, the majority of patients suffering from metastatic melanoma present with high relapse rates, demanding for alternative therapeutic options. We have previously developed a drug repurposing workflow to identify metabolic drug targets that, if depleted, inhibit the growth of cancer cells without harming healthy tissues. In the current study, we have applied a refined version of the workflow to specifically predict both, common essential genes across various cancer types, and melanoma-specific essential genes that could potentially be used as drug targets for melanoma treatment. The in silico single gene deletion step was adapted to simulate the knock-out of all targets of a drug on an objective function such as growth or energy balance. Based on publicly available, and in-house, large-scale transcriptomic data metabolic models for melanoma were reconstructed enabling the prediction of 28 candidate drugs and estimating their respective efficacy. Twelve highly efficacious drugs with low half-maximal inhibitory concentration values for the treatment of other cancers, which are not yet approved for melanoma treatment, were used for in vitro validation using melanoma cell lines. Combination of the top 4 out of 6 promising candidate drugs with BRAF or MEK inhibitors, partially showed synergistic growth inhibition compared to individual BRAF/MEK inhibition. Hence, the repurposing of drugs may enable an increase in therapeutic options e.g., for non-responders or upon acquired resistance to conventional melanoma treatments.

viability databases with more tested concentrations, allow determining more precise potency (such as IC 50 ) but on a smaller drug set. The three sets of drugs (predicted, NO-based, and anti-melanoma) were compared using the IC 50 measures in high-throughput pan-cancer cell viability databases. Four databases with IC 50 values were merged: Secondary PRISM (19Q4) (Corsello et al., 2020), GDSC1000 (Release 8.4) (Yang et al., 2013), GDSC2000 (Release 8.4) (Yang et al., 2013), and Genentech Cell Line Screening Initiative (Haverty et al., 2016). As fotemustine was the only anti-melanoma missing in the merged IC 50 database, fotemustine IC 50 in melanoma cell lines were added manually from the literature. The melanoma cell lines in the secondary database were further classified based on the BRAF and NRAS mutation status to identify sensitive candidate drugs regardless of the cell line mutation profile. Mutation status retrieved from the DepMap website ("Binary Calls for Copy Number and Mutation Data") was used to classify the melanoma cell lines in the secondary PRISM database based on BRAF and NRAS mutations.
Dose-response curves and determination of IC 50 values: Melanoma cells were seeded at a density of 0.5×10 cells/well in 96-well black µclear plates (Greiner). 3-fold dilution series of each drug (ranging from 0.05 to 100000nM, depending on the drug) were assayed in technical triplicates for 72h (Margue et al., 2019).

Proliferation assay
Briefly, 624mel and SKMel30 cells were seeded at a density of 5000 cells/well in a 96-well µclear plate (655090, Greiner). Cells were treated with either the drugs alone or combinations thereof for 72h at concentrations shown in Supplementary Table S  Fluorescence intensity was measured on Cytation 5 (Biotek). Quality Control concerning the essential genes predictions

Models capture metabolic variations between conditions
Models of the same tissue or condition are expected to cluster together as sharing a higher similarity (high quality). If the quality of the models is low, these are not able to capture metabolic variation among the conditions. Hence, to assess the quality of the context-specific reconstruction process, the sample were clustered based on the Jaccard Similarity Index.
Overall, models separated mainly along with their tissue of origin and then in function of their status (cancer or control), as previously described for the TCGA dataset (Pacheco et al., 2019), showing that the models were able to capture metabolic variations between the different tissues and also between healthy and cancer models. For the CCLE models, no clear separation in clusters was observed, as all models are derived from melanoma cell lines and do not include replicates (Supplementary Figure 1). However, for the IN-HOUSE melanoma models, the replicates in the sample models clustered together and in the consensus model, three major clusters composed of the control melanocytes, the A375 cell lines a low metastatic melanoma cell line), together with the metastatic samples, and the remaining cell lines could be observed. The third cluster could further be separated for the brain metastasis, which formed a group with WM1346 and WM1366 and another composed of MeWO, SKMel5, Malme3M. Furthermore, the A375 cell line shared metabolic similarities with both the metastatic samples and to some extent with the control melanocytes (Supplementary  The similarity between the different models is greater than 75%. For each model pair in the sample-specific and consensus selection, the Jaccard similarity score was determined to assess how similar the models are, based on the reaction presence. The data is then represented in a clustergram that allows finding clusters of similar models. For the melanoma models, three clusters can be observed for the sample-specific and consensus models. One cluster comprises the control melanocytes, the second cluster includes mainly the patient-derived Tumel models and the A375 cell line with the IZI treated, and the third cluster is formed by the remaining cell lines. 30 essential genes are common to TGCA, SKCM, INHOUSE models 30 essential genes are common to TGCA, SKCM, IN-HOUSE models but the knockout of these genes does no affect ATP maintenance in the liver and kidney models. Figure S 4: Conserved essential genes across the different melanoma, cancer cell types and healthy tissues. Venn diagram presenting the overlaps of predicted essential genes between the consensus IN-HOUSE models, the TCGA, the SKCM, TCGA CONTROL LIVER and TCGA CONTROL KIDNEY. Common versus context-specific fitness essential genes How common or how specific a gene is, was estimated by counting the number of the sample models that predicted this gene to be essential. Common melanoma or cancer genes, if druggable, might be interesting drug targets as are likely to have a low NNT. In the main text, the analysis was shown for the 40 predicted essential genes by IN-HOUSE consensus melanoma models. Supplementary Figure 7 displays the same analysis for the top essential genes in the IN-HOUSE sample models. The same analysis was performed on the TCGA dataset (Supplementary Figure 5) to assess if a gene is a common essential gene, but also on SKCM and CCLE sample models (Supplementary Figure 6). 0 5000 10000  In silico gene deletion analysis that shows (in x-axis) the number of models that are affected by the deleted gene. The y-axis represents the essential genes, sorted by growth ratio in the TCGA models. The growth ratio is presented in the colour scale from dark (complete reduction of the objective function) to white (no effect on the objective function). The deletion of a set of genes (blue band) completely shuts down the biomass reaction. The deletion of the gene set within the green band only shut down biomass in some models. In the purple band are genes that are implicated in oxidative phosphorylation. The genes in the green, blue and purple bands are mostly conserved across the TCGA, SKCM, CCLE models. The arrows show the essentially for some genes of interest.   Number of models

Number of models
Conserved essential genes across the IN-HOUSE models. The upper part of the left panel includes genes whose deletion completely reduces the biomass reaction to 0 in all of the 32 cancer models of the IN-HOUSE melanoma dataset (blue band). Below are sets of genes that reduce the biomass reaction to 0 in most cancer models but not all (green band). Differential effects can be observed for some genes (purple and turquoise bands) that affect only a few cancer samples and to different degrees. Some sets of genes that reduce the biomass production in the cancer models also reduce the ATP production in the control models (purple band). The arrows show the essentially for some genes of interest.

Quality control concerning the drug candidate predictions
The drug predictions of the sample models mostly overlap with the consensus models The number of shared drugs is displayed for every cancer consensus and sample-specific model. The largest number of drugs was found for the SKCM consensus model whereas only 31 drugs were shared with the sample-specific models. The melanoma consensus and sample-specific models predicted the same drugs. The lowest number of drugs was predicted by the CCLE model, but they were also predicted by every other model. Some predicted drug predictions were already FDA approved for cancer As an additional quality control, cancer data from drug databases were mined. 12 of the 28 candidate drugs have already been approved for antineoplastic therapy in at least one of the Figure S 10: Shared predicted drugs across the different cancer tissues of the TCGA based on the consensus models. The number of shared drugs between the different tissue-specific consensus models is shown. A cluster of cancer tissues (LAML, BLCA, UCS, SKCM, and CESC) with a high number of predicted drugs (54 to 148 drugs) that are also partially shared can be observed as well as another cluster (LGG, ACC, COAD, GBM, PRAD, LIHC, OV, THCA, READ and KICH) that share between 38 and 52 drugs can be observed, suggesting similar cancer liabilities in these clusters.
investigated databases (see Supplementary File 2). Further, only one database (SEER*RX) returned a significant enrichment for known antimetabolites (p-value: 0.0091), chemotherapeutic agents (p-value: 0.0300), and both (p-value: 0.0025), while no significant enrichment was found in other tested databases (see Methods). One possible explanation for this is the specificity of the predicted drugs that only target metabolic genes, whereas the databases include all the possible cancer drugs, which is further shown by the low number of found to target metabolic genes included in the Recon 2.04 reconstruction. Interestingly, no enrichment was found for melanoma-related anti-cancer agents in the SEER*Rx database, suggesting that we have predicted novel drugs for repurposing in melanoma, see Supplementary Table   S3. Predicted drugs and their targets show higher viability reduction and dependency, respectively, than anti-melanoma drugs and their targets, in both resistant and metastatic melanoma cell lines The viability reduction were retrieved for 25 predicted, 9 anti-melanoma and 29 NO-based drugs from the primary PRISM database. Aminomethyltransferase enzyme (https://clue.io/repurposing-app?q=Name:aminomethyltransferase) was mislabeled as a compound in the Drug Repurposing Hub: thus, it was discarded from the NO-based drugs. Three predicted drugs (cladribine, fluvastatin, and gemcitabine) ranked higher than anti-melanoma drugs in both metastatic (Supplementary Figure S11) and resistant cell lines (Supplementary Figures S12 and S13). Among the 29 NO-based drugs, diphenyleneiodonium is the only drug with a median viability reduction higher than 50% in either metastatic or resistant cell lines.
The high ranking of the non-anti-cancer drug, fluvastatin, among the predicted anti-cancers (gemcitabine and cladribine) in both metastatic and resistant cell lines, underlines fluvastatin as a potential anti-melanoma drug.
In both metastatic and resistant cell lines, predicted drug main targets have stronger dependency than anti-melanoma targets, with RRM1 and RRM2 scoring the highest, see Supplementary Figures S14 and S15. Main drug targets (anti-melanoma and predicted drugs) were selected from the drug targets using using the Drug Repurposing Hub (Corsello et al., 2017) for its high manual curation. Unlike the primary PRISM database, DepMap offers eight drug-induced resistance cell lines that are corresponding to four wildtype cell lines. ASS1 was the highest NO-related gene with a 30% dependency probability in metastatic cell lines. Rank of the predicted drugs in the Primary PRISM database in the melanoma cell lines out of 4604 drugs Figure S 11: Cladribine, fluvastatin and gemcitabine and the NO-based drug diphenyleneiodonium rank better than anti-melanoma drugs in viability reduction assays. The viability reduction relative to DMSO for our 28 predicted, anti-melanoma and NO-based drugs were gathered from the primary PRISM database for metastatic, primary, and uncategorized cell lines. The drugs were ranked by their median viability reduction in metastatic cell lines. X-axis represents the median viability reduction relative to DMSO (%) for metastatic (left), primary (middle) and uncategorized (right) cell lines. The rank of each drug out of the 4606 tested drugs in the primary PRISM database was printed beside each bar. Drug sensitivity of the anti−melanoma drugs in the Primary PRISM database Figure S 12: Over a quarter of the melanoma cell lines displayed resistance to anti-melanoma drugs in the primary PRISM database. Effect of drugs ranges from viability reduction (in blue) to increase of proliferation (up to more than 2-fold increase, in red). Cell lines with at least 50% increased proliferation (-50% viability reduction) to any anti-melanoma drug, were considered resistant cell lines (names highlighted in red) for visualisation in Figure S13. The x-axis represents the anti-melanoma drugs and the y-axis represents the melanoma cell lines. Missing drug-cell line experiments are marked in grey tiles.  The knockout of the predicted drug targets and essential genes induced a stronger viability reduction than the melanoma targets and NO-related genes. Predicted drug targets and essential genes were compared to both the anti-melanoma targets and NO-related genes based on gene dependency probability scores (probability to induce cell death or stop cell growth upon knock-out) from DepMap in metastatic melanoma cell lines. X-axis represents the median dependency probability across three types of melanoma cell lines (metastasis, primary, and uncategorized). Essential genes are highlighted in bold and the non-druggable essential genes (according to DrugBank v5.1.3) have yellow bars. Genes are ranked on the y-axis by the median dependency probability in the metastatic cell lines, with the gene rank displayed on the right of each bar.  Sensitivity analysis using DepMap dependency probability for the predicted drug targets and essential genes in the drug−resistant melanoma cell lines Figure S 15: The knockout of predicted essential genes induce a higher viability reduction in resistant cell lines. The dependency probability of the predicted drug targets and essential genes were compared to anti-melanoma targets and NO-related genes using the drug-induced resistance (violet, brown, green and blue bars) and wildtype (red bars) cell lines from the DepMap database. Genes (predicted drug targets, essential genes, anti-melanoma targets and NO-related genes) were sorted by the median dependency probability across the eight resistant cell lines.
anti-melanoma drugs in melanoma cell lines from high-throughput drug viability databases.