Analysis of CRISPR‐Cas9 screens identifies genetic dependencies in melanoma

Abstract Targeting the MAPK signaling pathway has transformed the treatment of metastatic melanoma. CRISPR‐Cas9 genetic screens provide a genome‐wide approach to uncover novel genetic dependencies that might serve as therapeutic targets. Here, we analyzed recently reported CRISPR‐Cas9 screens comparing data from 28 melanoma cell lines and 313 cell lines of other tumor types in order to identify fitness genes related to melanoma. We found an average of 1,494 fitness genes in each melanoma cell line. We identified 33 genes, inactivation of which specifically reduced the fitness of melanoma. This set of tumor type‐specific genes includes established melanoma fitness genes as well as many genes that have not previously been associated with melanoma growth. Several genes encode proteins that can be targeted using available inhibitors. We verified that genetic inactivation of DUSP4 and PPP2R2A reduces the proliferation of melanoma cells. DUSP4 encodes an inhibitor of ERK, suggesting that further activation of MAPK signaling activity through its loss is selectively deleterious to melanoma cells. Collectively, these data present a resource of genetic dependencies in melanoma that may be explored as potential therapeutic targets.

metastatic melanoma, targeted therapy using BRAF and MEK inhibitors can lead to significant tumor regression. Almost invariably melanoma cells acquire resistance to these targeted treatments, and disease relapse occurs. Consequently, there is a need to identify additional genetic dependencies that might serve as therapeutic targets.
The application of genome-wide screens in human cancer cell lines has the potential to identify genetic dependencies that may be targeted therapeutically (Barretina et al., 2012;Thompson, Adams, & Ranzani, 2017). Genome editing with CRISPR-Cas9 technology has improved the identification of genetic dependencies due to its high precision and limited off-target effects. In CRISPR-Cas9 dropout screens a population of cells is transduced with a pooled sgRNA library and following culture selective depletion of sgRNAs is measured to identify genes associated with a growth disadvantage or lethal phenotype, designated as fitness genes (Liu & Li, 2019;Shalem et al., 2014). Core fitness genes are involved in essential processes that cells depend on for survival and proliferation. In addition, context-dependent fitness genes that are specific for cell lineage or genotype are distinguished. Recently, genome-wide CRISPR-Cas9 screens have been performed in a range of cancer cell lines, yielding cancer-specific fitness genes (Bakke et al., 2019;Behan et al., 2019;Dempster et al., 2019;Picco et al., 2019;Tzelepis et al., 2016;T. Wang et al., 2015). Importantly, the majority of cancer-specific fitness genes were found to be limited to only one or two tumor types . Specific genetic dependencies in cancer cells may constitute targetable therapeutic vulnerabilities. The objective of this study was to identify novel genetic dependencies that may serve as potential therapeutic targets in melanoma cells through analysis of CRISPR-Cas9 screen data.
Among the 33 fitness genes that we define in melanoma, there are multiple genes that have not previously been associated with melanoma growth, including inhibitors of MAPK signaling activity.

| CRISPR-Cas9 screen data
The generation of CRISPR-Cas9 screen data at Broad Institute, Cambridge Massachusetts, available online at https://depmap.org/ ceres/, was described previously (Meyers et al., 2017). Briefly, 341 tumor cell lines including 28 melanoma lines were engineered to express Cas9 and subsequently screened using the human Avana4 library composed of 70,086 sgRNAs, targeting 17,670 protein-coding genes (4sgRNAs per gene) and 995 non-targeting control sgRNAs (Meyers et al., 2017). Cancer cell lines were transduced at a multiplicity of infection (MOI) of 0.3 to ensure that each cell expresses only one sgRNA. Genomic DNA was purified from transduced cells cultured under puromycin selection at day 1 and day 21 for nextgeneration sequencing. The cell lines included in the screen expressed Cas9.

| Pan-cancer analysis to determine melanoma fitness genes
The CRISPRcleanR package was applied to process the CRISPR-Cas9-derived essentiality profiles and to correct for copy-number amplifications, associated with exacerbated vulnerability scores (Meyers et al., 2017;Iorio et al., 2018). CRISPR-Cas9 screen data were analyzed using an R implementation of the BAGEL (Bayesian analysis of gene essentiality) algorithm, generating a scaled Bayesian factor (BF) score per gene Hart & Moffat, 2016).
A 5% false discovery rate (FDR) cutoff was applied. The mutation annotation for each melanoma cell line was derived from the cell line encyclopedia (CCLE) (Ghandi et al., 2019). Gene-level BFs were computed by calculating the average of the BFs across sgRNAs targeting a gene. This algorithm uses reference sets of predefined essential and non-essential genes. Each gene was assigned a scaled BF computed by subtracting the BF at the 5% FDR threshold (obtained from classifying reference essential/non-essential genes using BF rankings) from the original BF. Those genes with a statistically significant depletion at 5% FDR had a scaled BF above zero. Fitness genes were determined by comparing the average dropout of sgRNAs targeting the same gene, to that of reference essential and non-essential genes Hart & Moffat, 2016). Scaled BFs were binarized to 0 (scaled BF<0) and 1 (>0). A Fisher's exact test was performed on a two-way contingency table of fitness and non-fitness genes with binarized scaled BF scores from melanoma and the other tumor cell lines with the resulting p-values corrected for multiple testing using the Benjamini-Hochberg procedure (adjusted p-value<.01). Gene expression data for the 28 melanoma cell lines were available for analysis from the CCLE data portal. The data sets are available at https://data.broad insti tute.org/ccle/CCLE_RNAseq_081117.rpkm.

Significance
Genome-wide CRISPR-Cas9 genetic screens in tumor cell lines allow systematic identification of genetic dependencies that may be cancer type-specific. Here, we analyzed available data from such screens aimed at identifying genes, inactivation of which impairs proliferation of melanoma cells. A subset of the melanoma fitness genes encode proteins that can be inhibited with available small molecule inhibitors. Our data show that for growth melanoma cells not only depend on activators but also on inhibitors of MAPK signaling activity such as DUSP4. The melanomaspecific genetic dependencies identified might be explored as novel specific therapeutic targets.
The sgRNAs used for validations differed from those used in the Broad Avana4 library to provide orthogonal validation. Lentivirus stocks were produced following transfections into HEK293T cells using polyethylenimine (PEI) (Carlotti et al., 2004). Viral titers were determined by antigen capture ELISA measuring HIV p24 (ZeptoMetrix Corp.). Following lentiviral transduction with the Cas9 expression vector to reach a MOI of 3 (100% infection efficiency) in A375, IGR1, and WM983B cells, the Cas9-editing efficiency was tested . Capillary sequencing analysis was performed using the human U6 promoter forward primer (GACTATCATATGCTTACCGT) to align the sgRNA sequence to the backbone vector LV04 (Sigma-Aldrich) and ensure that the sgRNA sequence was correct. For lentiviral transductions with sgRNA expression vectors, 10 5 cells were seeded in 6-well plates or 2.5 × 10 5 cells in 6-cm dishes. Cells were transduced next day in 2 ml or 5 ml culture medium supplemented with 8 μg/ml polybrene to reach a MOI of 3. Cells were incubated overnight at 37°C and cultured for 4 days in fresh medium supplemented with blasticidin-S (5 μg/ml) and puromycin (2 μg/ml) for selection and further analysis.
Cells were fixed for 5-10 min in 4% paraformaldehyde (PFA) and stained for 30 min with crystal violet (0.05%) when control wells transduced with negative control sgRNAs targeting SSX3 reached 80% confluency. Scanned images of the wells were obtained before solubilization of retained dye with 100% methanol to measure absorbance at OD540.

| Cell confluency assay
IGR1-Cas9 and A375-Cas9 cells were seeded in 96-well plates to monitor cell confluency using the IncuCyte live-cell analysis system (Essen BioScience). For this assay cells were seeded at a density of 50 cells/well (A375) and 200 cells/well (IGR1) in 6 replicates and scanned images of each well were taken every 12 hr. Data were analyzed using IncuCyte software and confluency was calculated over 6 days (A375) and 8 days (IGR1) and normalized to the confluency on day 1 to correct for seeding variation. Colony formation and cell confluency assays were performed in biological duplicates, and combined data were analyzed using GraphPad Prism version 8.
A 2-way ANOVA and Bonferroni's multiple comparisons test were performed to detect statistically significant differences in cell confluency (p-value<.01).

| CRISPR-Cas9 screen analysis for fitness genes in melanoma
To identify genes that are specifically required for the fitness of melanoma cells, data available from CRISPR-Cas9 screens performed at the Broad Institute (Cambridge, Massachusetts) were extensively processed and analyzed. In these genome-wide CRISPR-Cas9 dropout screens, the depletion of sgRNAs from a genome-wide library was measured in populations of Cas9expressing tumor cells following 21 days of culture after transduction (Meyers et al., 2017). Fitness genes were determined by comparing the average dropout of sgRNAs targeting the same gene to the profile of reference essential and non-essential genes using a supervised approach called BAGEL Hart & Moffat, 2016). Of all targeted genes, 4,423 (25%) negatively affected the fitness of one or more melanoma cell lines. We found an average of 1,494 fitness genes in each of the 28 melanoma cell lines. There were 1,396 genes that reduced the fitness in half or more of the melanoma cell lines, a number similar to what was reported for other cancer types . Our analysis yielded 193 genes that reduced fitness in all melanoma cell lines, 175 of which (91%) had been found to be core fitness genes in the haploid cell line HAP1 and are not specific to melanoma (Table S2) (Blomen et al., 2015).
To identify tumor type-specific genetic dependencies for melanoma, we performed comparative analysis of the scaled BFs rep-  Table S4) (Hemesath et al., 1994;Shakhova et al., 2012). A second cluster of fitness genes related to melanoma consisted of components of the MAPK signaling pathway, such as BRAF and MAPK1 (Table S5). The finding of these established fitness genes in this gene set confirms the sensitivity of the CRISPR-Cas9 screen analysis. A third group of genes, to which MDM2 belongs, is involved in regulating p53 activity (Table S4). For several genes, a role in melanoma progression has been reported; FERMT2 for instance has been found to impact melanoma metastasis (Karras et al., 2019). Moreover, activation of AHR was reported to promote resistance to BRAF inhibitors in melanoma (Corre et al., 2018). In addition, a set of fitness genes related to melanoma with undefined roles in melanocyte biology or melanoma pathogenesis was identified, such as MTMR6 and CRTC3. Pharmacological compounds are available that may halt melanoma growth for the products of 12 of the identified fitness genes related to melanoma, including AHR and MDM2 (Table S4).
F I G U R E 1 Identification of genetic dependencies in cutaneous melanoma cells. (a) Heatmap representation of genes significantly associated with cellular fitness in melanoma cell lines indicated in a box. The scale bar represents scaled Bayesian factors (BFs) within each screen, calculated by subtracting the BF at the 5% FDR threshold, obtained when classifying prior known essential/non-essential genes based on their BF's rank . Red color indicates genes that are likely to be fitness genes and, therefore, have a positive scaled BF. Blue color indicates genes less likely to be important for cell fitness with a negative scaled BF. Tumor types are shown and clustered. Genes were ranked according to the noted Fisher exact test-adjusted p-values. (b) String protein interaction network for 33 significant fitness genes in melanoma cell lines. Color coding: Red-proteins involved in the MAPK signaling pathway; blue-melanocytelineage-specific proteins; green-p53-regulatory pathway components; purple-unknown interaction network, miscellaneous function; dark-core protein components; light-regulatory components. For proteins indicated in underlined-bold text pharmacological inhibitors are currently available  MITF  MAPK1  DUSP4  ZEB2  CRTC3  SOX9  TCEB3  PEA15  PPM1G  IRF4  MAP2K2  TFAP2A  FERMT2  EGLN1  PAX3  INTS12  AHR  PPP2R2A  TXNL1  GTF2H5  MTMR6  MDM2  FAM25C  MUC12  NFATC2  FBXW11  RFWD2  CHMP4B  ARHGAP11A BPTF COASY

| Negative regulators of MAPK signaling are fitness genes in melanoma
Remarkably, the identified MAPK signaling genes encoded not only activating but also inhibitory components such as DUSP4 and PEA15. DUSP4 dephosphorylates ERK1/2, in addition to p38 and JNK (  Proteins encoded by the target genes BRAF, DUSP4 and PPP2R2A were significantly depleted by CRISPR-Cas9-mediated inactivation using two independent sgRNAs, as was confirmed by immunoblot analysis 4 days after transduction ( Figure S2a). The sgRNAs against BRAF also led to a decrease in DUSP4 protein levels, in accordance with regulation of DUSP4 by MAPK signaling as part of a negative F I G U R E 2 Significant fitness genes in 28 melanoma cell lines. Heatmap representation of fitness genes specifically important in cutaneous melanoma cell lines and indication of established driver mutations (BRAF, NRAS, NF1) for each cell line. Supervised clustering was performed on columns to distinguish clusters of BRAF, NRAS, and NF1-mutant and wild-type cell lines. The genes are ranked according to fitness effect that is denoted by the scaled BF score within the set of 28 melanoma cell lines. Genes that are high on the heatmap had a more positive scaled BF score (red) indicating that are more likely to be important for fitness, whereas genes that are low on the heatmap (blue) had a more negative scaled BF score indicating that are less likely to be important for fitness  (Hutchinson et al., 2015). Upon CRISPR-mediated inactivation of DUSP4 and PPP2R2A, a significant effect on WM983B cell viability was observed as measured using the crystal violet assay ( Figure S2b). Loss of proliferation caused by DUSP4-mediated depletion was similar to BRAF-mediated depletion in WM983B cells ( Figure S2c). The effects of inactivation of PPP2R2A using both sgR-NAs were slightly less pronounced in this cell line. Combined, these data suggest that genetic depletion of DUSP4 and PPP2R2A has a significant effect on the proliferation of WM983B cells.
Next, we examined the effects of DUSP4 and PPP2R2A depletion in two additional melanoma cell lines, both harboring BRAF mutations. Whereas the IGR1 melanoma cell line was also included in the CRISPR-Cas9 screens, the A375 melanoma cell line was not.
Immunoblot analysis confirmed depletion of PPP2R2A and DUSP4 upon sgRNA transduction (Figure 3a). We could verify that PPP2R2A inactivation affects the viability of IGR1 cells 8 days after seeding according to the colony formation assay (Figure 3b

| D ISCUSS I ON
Here, we present the results of analysis of genome-wide CRISPR- have an as yet undetermined role in melanoma biology. For 12 proteins encoded by genes identified as essential in the present study, pharmacological compounds are available, implying that certain existing drugs might be efficacious in the treatment of melanoma.
Strengths of this study are the sensitivity and robustness of the CRISPR-Cas9 screening methodology by targeting all protein-coding genes using 4 sgRNAs per gene and data analysis through normalization of copy number-associated effects, as well as comparative analysis of genetic dependencies in 28 melanoma cell lines with those of 313 other tumor cell lines. Evaluation of CRISPR-Cas9 screens in human melanocytes would have allowed further delineation of genes that are essential in melanoma from lineage-specific fitness genes, but such data are not yet available. As our identified hits include the melanocytic lineage transcription factors MITF and SOX10, it is probable that some of the other genes affect melanocyte fitness as well. Eleven of the identified 33 genes have been previously identified as fitness genes in haploid, CML-derived HAP1 cells (BRAF,MAPK1,ZEB2,TCEB3,FERMT2,PPP2R2A,BPTF CHMP4B,INTS12,FBXW11,and COASY) (Blomen et al., 2015). Inactivation of these genes is significantly more detrimental to melanoma cells than to other tumor cell types, but they are likely to be involved in essential cellular processes. The composition of the fitness gene set related to melanoma may be determined in large part by dependencies associated with mutant BRAF. In a previous study, gene dependency associations with BRAF mutation have been identified in 16 cancer types (Dempster et al., 2019). Mutant BRAF was present in nine of 146 of cancer cell lines, and this oncogenic mutation was found to be associated with dependency on 50 genes in these cancer types.
We demonstrate genetic dependency on multiple MAPK signaling components in melanoma. Interestingly these include not only activators but also inhibitors of MAPK signaling, such as DUSP4 and PEA15. DUSP4 regulates phosphorylation of ERK, but also of p38, JNK and other proteins (Mazumdar et al., 2016). DUSP4 depletion was recently reported to diminish the negative effects on melanoma cell viability induced by MEK inhibitors through increasing MAPK activity (Gupta et al., 2020). Our results support the notion that further activation of MAPK signaling through loss of inhibitors such as DUSP4 in melanoma cells that already harbor activating BRAF or NRAS mutations is detrimental, a hypothesis that will need to be explored with further experiments. Accordingly it has been reported that ERK1 and ERK2 overexpression results in cell death in BRAF and NRAS-mutant melanoma cells (Leung et al., 2019). Most melanoma cell lines included in the CRISPR-Cas9 screens carried the BRAF mutation, but screen data suggested that some NRAS-mutant melanoma cells might also be sensitive to inactivation of DUSP4 and PPP2R2A. Targeting these proteins may provide an alternative approach to treatment of metastatic melanoma, particularly after relapse from immunotherapy or targeted therapy. Acquisition of resistance to BRAF inhibitors commonly involves reactivation of MAPK signaling. Discontinuation of BRAF inhibitor treatment in resistant melanoma cells results in hyperactivated MAPK signaling, which may be detrimental to these cells (Kong et al., 2017;Sun et al., 2014). We hypothesize that once melanoma cells have acquired resistance to treatment with BRAF and MEK inhibitors by upregulation of MAPK activity they will become more sensitive to DUSP4 inhibition. Targeting negative regulators of MAPK signaling inhibitors such as DUSP4 and PEA15 further activates MAPK F I G U R E 3 Validation of CRISPR-mediated inactivation of DUSP4 and PPP2R2A on proliferation of IGR1 and A375 cutaneous melanoma cell lines. (a) Immunoblot analysis of BRAF, PPP2R2A, and DUSP4 protein expression 4 days after transduction with sgRNA expression vectors in IGR1 and A375 cells. Expression of vinculin was investigated as a loading control. (b) Images of crystal violet assays of a control depletion (SSX3), BRAF, and two independent sgRNAs for PPP2R2A and DUSP4 in IGR1 and A375 cells. (c) Comparison of two independent techniques, live-cell imaging in 96-well plates and cell viability assay in 12-well plates for different knockout IGR1 and A375 lines. A 2-way ANOVA and Bonferroni's multiple comparisons test was performed between the two techniques (*p<.01). Error bars represent SE of the mean. Two independent experiments were performed. (d) Graphical representation of cell proliferation up to 8 days after seeding of IGR1 cells and 6 days after seeding of A375 cells. sgCtrl represents depletion of SSX3, BRAF depletion was used as a positive control, and 2 independent sgRNAs were used for PPP2R2A and DUSP4. The normalized confluency (%) of IGR1 and A375 cells was corrected based on day 1 measurements. A 2-way ANOVA and Bonferroni's multiple comparisons test was performed between the control and all other lines (*p<.01) (a)  (Koutsioumpa et al., 2018), but also to promote MAPK signaling by regulating RAF and KSR (Hein et al., 2016;Ory, Zhou, Conrads, Veenstra, & Morrison, 2003). The PP2A complex has broader cellular functions, including regulation of oxidative stress signaling and DNA repair response (Martina & Puertollano, 2018;Omerovic et al., 2010;Yan et al., 2012). We have not investigated whether the fitness effects of DUSP4 and PPP2R2A are strictly dependent on their function as inhibitors of oncogenic kinome signaling.
The generation of effective treatment strategies in melanoma remains a challenge. Here, we present an analysis of CRISPR-Cas9 screen data aimed at melanoma cell lines, identifying 33 genes that specifically affect the fitness in this tumor type. In vitro experiments in human melanoma cell lines confirmed that inactivation of DUSP4 and PPP2R2A results in decreased cell proliferation. Collectively, these data present a resource of genetic dependencies in melanoma that may be explored as potential therapeutic targets.

ACK N OWLED G EM ENTS
We would like to thank the Broad Institute DepMap team for gen-

CO N FLI C T O F I NTE R E S T
The authors declare that there is no conflict of interest to disclose.

S U PP O RTI N G I N FO R M ATI O N
Additional supporting information may be found online in the Supporting Information section.