SATB2 induction of a neural crest mesenchyme-like program drives invasion and drug resistance in melanoma

Recent genomic and scRNA-seq analyses of melanoma identified common transcriptional states correlating with invasion or drug resistance, but failed to find recurrent drivers of metastasis. To test whether transcriptional adaptation can drive melanoma progression, we made use of a zebrafish mitfa:BRAFV600E;tp53-/- model, in which malignant progression is characterized by minimal genetic evolution. We undertook an overexpression-screen of 80 epigenetic/transcriptional regulators and found neural crest-mesenchyme developmental regulator SATB2 to accelerate aggressive melanoma development. Its overexpression induces invadopodia formation and invasion in zebrafish tumors and human melanoma cell lines. SATB2 binds and activates neural crest-regulators, including pdgfab and snai2. The transcriptional program induced by SATB2 overlaps with known MITFlowAXLhlgh and AQP1+NGFR1high drug resistant states and functionally drives enhanced tumor propagation and resistance to Vemurafenib in vivo. Here we show that melanoma transcriptional rewiring by SATB2 to a neural crest mesenchyme-like program can drive invasion and drug resistance in endogenous tumors.

Both sequencing of patient precursor lesions (Shain et al., 2015) and modeling in transgenic animals (Perez-Guijarro, Day, Merlino, & Zaidi, 2017;van Rooijen, Fazio, & Zon, 2017) showed that acquisition of oncogenic BRAF or NRAS mutations alone is insufficient to drive melanoma initiation. Reconstruction of the clonal history of advanced human melanoma suggests that few driver mutations, including loss of tumor suppressors, are acquired in early transformed melanocytes, and that high mutational burden and genomic instability are later events in malignant progression (Birkeland et al., 2018;Ding et al., 2014). A pan-cancer whole-genome analysis of metastatic solid tumors, including 248 melanoma patients, showed a surprisingly high degree of similarity in mutational landscape and driver genes between metastatic tumors and their respective primaries, as well as across multiple metastatic lesions within the same patient (Priestley et al., 2019). The lack of recurrently mutated genes in metastasis, and the limited genetic diversity across metastatic sites suggests that genetic selection of mutated drivers is unlikely to be responsible for metastatic progression. Despite the high degree of mutational burden and genetic intra-and inter-patient heterogeneity observed in melanoma, several bulk RNA-seq and scRNA-seq analyses of metastatic melanoma patient samples have identified common recurrent transcriptional states correlating with invasion (mesenchymal signature) (Verfaillie et al., 2015), drug resistance (MITF low /AXL high ) (Tirosh et al., 2016), a neural crest-like state correlating with minimal residual disease persistence during targeted therapy with MAPK inhibitors (MITF low /NGFR1 high /AQP1 high ) (Rambow et al., 2018), and even response to immunotherapy (NGFR1 high ) (Boshuizen et al., 2020). The existing literature on these often less proliferative cell states has been limited to in vitro perturbations, or relied on transplant models, and as such their role in oncogenesis outside of iatrogenic drug treatment remains less clear.
Using a zebrafish transgenic melanoma model, we have previously shown that in a cancerized field of melanocytes carrying driving oncogenic mutations, such as BRAF V600E and loss of tp53, only a few cells or subclones will go on to develop a malignant tumor, and that this event is marked by re-expression of developmental NC-markers sox10 and crestin McConnell et al., 2018). Similar to the observation in human metastatic patient samples, malignant progression from transformed melanocytes in experimental animal models is not explained by the very limited genetic evolution observed in these tumors (Yen et al., 2013).
Taken together, these observations in human patients and animal disease models raise the question whether epigenetic adaptation, rather than genetic selection, can contribute to tumor progression by altering the transcriptional state of the tumor.
Here, we address this hypothesis by high throughput genetic perturbation of 80 epigenetic/transcriptional regulators in endogenous primary tumors, by leveraging a transgenic zebrafish melanoma model uniquely suited to precisely and rapidly perturbate tumor development in vivo (Ceol et al., 2011), at a scale and speed much greater than possible with genetically engineered murine melanoma models. We identify Special AT-rich Binding protein 2 (SATB2) as a novel accelerator of melanoma onset driving an aggressive phenotype in primary tumors suggestive of metastatic spreading. SATB2 acts as a transcriptional regulator by recruiting members of the Acetyl Transferase and Histone demethylase complexes to target genes and altering the local chromatin organization and activation state (L. Q. Zhou et al., 2012). Its expression has been shown to correlate with patient outcome in several tumor types (Naik & Galande, 2019;Yu, Ma, Shankar, & Srivastava, 2017). In different tissues and tumor types SATB2 has been shown to play a role in oncogenic transformation and proliferation, or Epithelial-to-Mesenchymal Transition (EMT), migration, and self-renewal (Gan et al., 2017;Naik & Galande, 2019;Nayak et al., 2019;Wang et al., 2019;Wu et al., 2016;Xu et al., 2017;Yu, Ma, Ochoa, Shankar, & Srivastava, 2017)). SATB2 is a transcription factor and chromatin remodeler with a well-conserved structure and expression pattern across chicken, mouse and zebrafish during the development and migration of the cranial neural crest (CNC) and neuronal development. It is required for the development of the exo-mesenchymal lineages of the CNC (Sheehan-Rooney, Palinkasova, Eberhart, & Dixon, 2010;Sheehan-Rooney, Swartz, Lovely, Dixon, & Eberhart, 2013), and neuronal axon formation (McKenna et al., 2015;Shinmyo & Kawasaki, 2017). In facts, SATB2 inactivating mutations in humans have been associated with cleft palate, intellectual disability, facial dysmorphism, and development of odontomas, defining a neurocristopathy referred to as SATB2-associated syndrome (Kikuiri et al., 2018;Zarate et al., 2019;Zarate & Fish, 2017). Through a combination of zebrafish in vivo allotransplants and validation in human melanoma cell lines, we show that SATB2 drives enhanced invasion via invadopodia formation and an EMT-like phenotype. Mechanistically, chromatin and transcriptional characterization of primary zebrafish SATB2 tumors vs. EGFP controls via ChIP-seq and RNA-seq shows SATB2 to bind and induce transcriptional activation of neural crest regulators, including snai2 and pdgfab. The transcriptional program induced by SATB2 overexpression is conserved between zebrafish and human melanoma, and overlaps with the aforementioned MITF low /AXL high (Tirosh et al., 2016) and neural crest-like MITF low /NGFR1 high/ AQP1 high drug resistant states (Rambow et al., 2018). Finally, we show SATB2 transcriptional rewiring to functionally drive enhanced tumor propagation and resistance to MAPK inhibition by Vemurafenib in zebrafish tumor allografts in vivo.

In vivo overexpression screen of epigenetic factors identifies SATB2 as melanoma accelerator
To interrogate whether epigenetic reprogramming can accelerate melanoma development, we utilized a genetic discovery driven approach and undertook an in vivo overexpression screen.
We tested 80 chromatin factors: 15 pools of 5 factors, and 6 additional single factors including previously published positive controls SETDB1, SUV39H1 (Ceol et al., 2011) (see Supplementary   Table 1). EGFP was used as a negative control, and we tested CCND1 as an additional positive control (known driver often amplified in melanoma (Cancer Genome Atlas, 2015)). As a screening platform, we leveraged a zebrafish melanoma model driven by tissue-specific expression of human oncogenic BRAF V600E in a tp53 and mitfa-deficient background (Ceol et al., 2011). Tg(mitfa:BRAF V600E ); tp53-/-; mitfa-/zebrafish lack melanocytes and do not develop melanoma. Mosaic integration of the transposon-based expression vector MiniCoopR (MCR), rescues melanocyte development by restoring mitfa, while simultaneously driving tissuespecific expression of candidate genes (Fig. 1A). Thus in mosaic F0 transgenics all melanocytes express the candidate gene tested (described in Ceol et al. (Ceol et al., 2011)). We identified 6 candidate pools, of which Pool F ( Figure 1B) had the strongest acceleration. Single factor validation of Pool F ( Figure 1C) showed that SATB2-overexpression (Supplemental Figure 1A Figure 2B-C). The acceleration phenotype is not due to increased cellular proliferation, since SATB2overexpression in primary zebrafish tumors (Supplemental Figure 3A) or in a panel of human melanoma cell lines via a TETon tetracycline inducible lentiviral vector (here referred to as iSATB2) did not result in increased proliferation (Supplemental Figure 3B,C). Our data demonstrate that SATB2-overexpression accelerates melanoma malignant progression without affecting proliferation, but SATB2 is neither necessary nor sufficient for melanoma initiation.

SATB2 overexpression leads to invadopodia formation and increases migration and invasion in vitro and in vivo
Given the lack of proliferation changes, we next asked whether SATB2's aggressive tumor phenotype and internal tumors might be explained by an increase in tumor migration and/or invasion. Using primary zebrafish melanoma in vitro cell cultures we performed scratch migration assays, that confirmed a heightened migration potential of 3 independent MCR:SATB2 cell lines compared to 3 independent MCR:EGFP cell lines (Supplementary Figure   4A). Cytoskeletal staining of zebrafish melanoma in vitro cultures (Zmel1 MCR:EGFP and 45-3 MCR:SATB2) with phalloidin revealed the presence of strong F-actin positive foci in MCR:SATB2 cells (Figure 2A), reminiscent of invadopodia. Invadopodia are cell protrusions involved in metastatic spreading by facilitating anchorage of cells to, and local degradation of the ECM (Murphy & Courtneidge, 2011), and are known to be regulated by the PDGF-SRC pathway and EMT. During normal development, invadopodia-like physiologically equivalent structures called podosomes are utilized by neural crest cells to migrate (Murphy & Courtneidge, 2011;Murphy et al., 2011). Co-localization of F-actin with invadopodium structural component Cortactin in MCR:SATB2 melanoma cell lines confirmed these foci to be invadopodia ( To further validate whether the internal organ involvement in MCR:SATB2 primary tumors was due to an increased migratory potential in vivo, we used an orthotopic allotransplantation model in sub-lethally irradiated transparent casper zebrafish (Heilmann et al., 2015;P. Li, White, & Zon, 2011;White et al., 2008). First, we allotransplanted zebrafish melanoma-derived cell lines Zmel1 (MCR:EGFP) vs. 45-3 (MCR:SATB2) into irradiated casper recipients. This showed SATB2 to cause an invasive histological phenotype and reduction in the recipient's overall survival compared to EGFP (45-3; n=31, median survival=25 days vs. Zmel1; n=31, median onset not reached at the experimental end point) (Supplemental Figure 5A-C). We then tested whether this difference was also present in primary tumors using an established in vivo migration assay (Heilmann et al., 2015) where we transplanted 300,000 primary pigmented melanoma cells into casper recipients, and monitored the formation of distant metastasis  Figure 6B). Human melanoma patient genomic datasets publicly available on cBio portal show SATB2 to be infrequently but recurrently amplified in ~4-8% of patients (Supplementary Figure   6C) in 3 independent datasets of metastatic melanoma (Hugo et al., 2016;Van Allen et al., 2015), and the high mRNA expression level correlate with poor survival in two independent metastatic melanoma patient datasets available on TIDE portal (GSE22153 and GSE8401, Supplementary Figure 6D).

SATB2 binds and transcriptionally regulates neural crest-and EMT-associated loci
To gain insight into the mechanism underlying the SATB2-overexpression phenotype in melanoma, we performed ChIP-seq and RNA-seq on primary zebrafish tumors. We conducted ChIP-seq on MCR:SATB2 primary zebrafish melanomas to identify SATB2-bound target genes (antibody validation in Supplementary Figure 7A). HOMER motif analysis of ChIP-seq significant peaks for de novo motif underlying SATB2 (Special AT-rich Sequence-Binding protein 2) binding showed the top motif to be AT-rich (Supplemental Figure 7B) as expected (Hassan et al., 2010;Naik & Galande, 2019;Savarese et al., 2009;L. Q. Zhou et al., 2012). GREAT analysis of GOterms associated with ChIP-seq peaks, showed SATB2 binding to be enriched at loci associated with NC development and migration, and Epithelial-to-Mesenchymal Transition (EMT) ( Figure   3A). Indeed, ~38.8% (127/327) of known NC-associated genes (Tan et al., 2016), including sox10, snai2, pax7, chd7, and semaforin (sema3fa) are bound by SATB2 (Supplementary Table   2). We next performed ChIP-seq for H3K27Ac and H3K9me3 histone marks, in MCR:SATB2 and MCR:EGFP control tumors to investigate the effect of SATB2 overexpression on the chromatin state of SATB2-bound targets. Globally, genome wide GO term analysis of loci with increased H3K27Ac deposition in MCR:SATB2 vs. MCR:EGFP also showed an increased activation of neural crest development (Supplementary Figure 7C). Furthermore, HOMER motif analysis of known motifs shows SATB2 ChIP-seq to be enriched for TFAP2A and RXR motifs, which are both transcription factors involved in neural crest specification (Supplementary Figure 7D). Locally, H3K27Ac deposition in MCR:SATB2 vs. MCR:EGFP around transcriptional start sites (TSS) of SATB2-bound genes suggested an increased chromatin activation state of SATB2-bound NC targets ( Figure 3B). SATB2's binding pattern is consistent with its role as a necessary specifier of cranial migratory neural crest differentiation in exo-mesenchymal tissues in the pharyngeal arches (which will develop into the jaw and teeth), and the consequent cleft palate defects observed in patients with mutated SATB2 in the human SATB2 syndrome (Zarate & Fish, 2017). Consistent with these prior literature findings we find satb2 to be highly expressed in a previously published dataset of migrating crestin:EGFP + -NC cells sorted from zebrafish embryos at the 15 somite stage , when EMT-mediated delamination initiates migration of the neural crest, including cranial (Supplementary Figure 7E) . Furthermore, knock-down of satb2 during embryonic development via microinjection of a validated splicing morpholino (Ahn et al., 2010) in the transgenic neural crest reporter zebrafish line Tg(sox10:mCherry) showed a reduction of sox10 reporter expression and craniofacial abnormalities, but did not affect melanocyte development (Supplementary Figure 7F).
To define the transcriptional effect of SATB2-overexpression we performed RNA-seq on 3 MCR:SATB2 and 3 MCR:EGFP tumors, and correlated SATB2-bound loci with their transcriptional changes in MCR tumors to identify SATB2 target genes that might underlie SATB2's phenotype ( Figure 3C). Given the inter-tumor transcriptional variability of primary zebrafish tumors (evident by qPCR on SATB2 itself in MCR:SATB2 tumors, Figure 3D), we conducted extended validation on a subset of these bound and transcriptionally altered genes, plus some additional manually curated SATB2-bound genes by qPCR on a large set of primary tumors ( Figure 3D). Genes that are SATB2-bound [within 3kb of the transcriptional start or end site, and the gene body] or SATB2-associated (predicted nearest gene association), which were altered by RNA-seq revealed that NC, cytoskeleton and extracellular matrix-associated gene programs are highly stimulated by SATB2-overexpression ( Figure 3D, Supplementary Table 2). chd7 and snai2 in particular have been described to be key regulators of the migratory NC (Bajpai et al., 2010;Okuno et al., 2017;Schulz et al., 2014), while pdgfab has been implicated in podosome and invadopodia formation (Ekpe-Adewuyi, Lopez-Campistrous, Tang, Brindley, & McMullen, 2016;Murphy & Courtneidge, 2011;Murphy et al., 2011;Paz, Pathak, & Yang, 2014). SATB2-binding to targets such as pdgfab (Supplemental Figure 7G Overall, SATB2 binds and activates chromatin at a subset of neural crest-related loci including snai2, resulting in an EMT-like phenotype. Overexpression of the human ortholog of downstream target snai2 partially recapitulates SATB2's phenotype in the MCR model.

SATB2-induced program is conserved across zebrafish and human melanoma and overlaps with known drug resistance transcriptional states
To ascertain whether SATB2-induced transcriptional changes were conserved across species, we inducibly overexpressed SATB2 in primary human melanocytes and human melanoma cell lines, and performed qPCR for selected orthologues of genes validated in zebrafish. qPCR analysis of SATB2 target genes after 48 hours of culture in the presence of doxycycline on iSATB2 melanoma cell lines (SKMEL2, A375 and SKMEL28 iSATB2), and iSATB2 untransformed primary human melanocytes (HEMA-LP) showed a similar increase in SNAI transcription factors ( Figure 4A). Cell line specific differences in overall SNAI1 or SNAI2 (both snai2 orthologs) induction levels may reflect context dependence (e.g. oncogene or cell line baseline status) (Nieto, Huang, Jackson, & Thiery, 2016). To analyze global transcriptional changes induced by SATB2 in human melanoma, we selected SKMEL2 iSATB2, since this cell line had the lowest endogenous SATB2 level and showed the strongest invadopodia induction ( Figure  To assess how the SATB2 program relates to previously described transcriptional states in melanoma we conducted GSEA analysis in SKMEL2 iSATB2 for known signatures: (1) a mesenchymal signature (Verfaillie et al., 2015), (2) MITF low /AXL high drug resistance state (Tirosh et al., 2016), and (3) neural crest state and other transcriptional signatures in different tumor cell subpopulations (Rambow et al., 2018 ) . This analysis showed a significant overlap of the SATB2-induced transcriptional program with the less differentiated MITF low /AXL high state (NES=1.54, p=0.000, q=0.17), and the Rambow neural crest-like/minimal residual disease MITF low /NGFR1 high /AQP1 high state driven by RXRG (NES=1.64, p=0.000, q=0.09) ( Figure 4D), which have both been previously described in metastatic human melanoma to correlate with a self-renewing MAPK inhibitor drug resistant state (Rambow et al., 2018;Tirosh et al., 2016). The key marker genes of these states (i.e. AQP1, RXRG, NGFR, AXL, Figure 4C drives a transcriptional induction of invadopodia, EMT and neural crest regulators in zebrafish and human melanoma alike, consistent with the transcriptional changes and phenotypes described above and its known role in development. Furthermore, the SATB2-induced transcriptional program shows strong overlaps with known drug resistance transcriptional states in melanoma.

SATB2 increases tumor propagating potential and resistance to MAPK inhibition in vivo
Given the overlap between the SATB2-induced program and the known minimal residual disease and MAPK inhibition resistant states, we asked whether SATB2 functionally confers enhanced self-renewal or resistance to MAPK inhibition. To address the former, we performed a limiting dilution tumor propagation assay injecting 300,000, 3,000, 300 or 30 primary Our chromatin and transcriptional analysis via ChIP-seq, and RNA-seq followed by qPCR validation on MCR:SATB2 vs. MCR:EGFP control tumors showed SATB2 to directly bind and transcriptionally induce several genes related to neural crest development (snai2 and chd7) (Bajpai et al., 2010;Betancur, Bronner-Fraser, & Sauka-Spengler, 2010;Mayor & Theveneau, 2013;Schulz et al., 2014), EMT (vim and fn1) (Betancur et al., 2010), and invadopodia formation (pdgfab and rhobtb2) (Maheswaran & Haber, 2015;Mayor & Theveneau, 2013;Murphy & Courtneidge, 2011;Nieto et al., 2016;Paz et al., 2014). Functionally, an EMT-like phenotype and ECM remodeling in MCR:SATB2 vs. MCR:EGFP tumors was also apparent at the protein level  Figure 3H). Taken together, these data suggest that downstream SATB2 target snai2 is partially responsible for the SATB2 phenotype, but that other genes in the program contribute to the full phenotype spectrum.
To further assess the conservation and relevance of our findings to human melanoma, we conducted qPCR analysis on a panel of iSATB2 human melanoma cell lines and primary cultured melanocytes (HEMa-LP). qPCR of the human orthologs of snai2, pdgfab, and sema3f overall showed a similar induction as zebrafish tumors, with some variation across cell lines ( Figure 3D and 4A). For an unbiased global transcriptional comparison, we conducted RNA-seq on human melanoma cell line SKMEL2 iSATB2 and identified transcriptional differences induced by SATB2 overexpression (with/without Dox induction for 48 hours). Significant differentially regulated genes from SKMEL2 iSATB2 were compared with zebrafish primary tumor RNA-seq (MCR:SATB2 vs. MCR:EGFP) using both Gene Set Enrichment Analysis (GSEA) and Ingenuity Pathway Analysis (IPA). Despite 1) the ortholog prediction constraints between human and zebrafish datasets that may have resulted in some expected loss of information, 2) the variation due to species genetic and functional conservation differences, and 3) the in vivo vs. in vitro conditions, we see strong conservation across the two models in the top altered pathways ( Figure 4B). Overall, our study shows a robust transcriptional and phenotypic conservation of SATB2 overexpression effects in melanoma across human and zebrafish ( Interestingly, SATB2 has been reported as a primary hit in an unbiased overexpression screen for MAPK inhibitor resistance drivers in human melanoma cell lines in vitro (Johannessen et al., 2013). We sought out to functionally validate SATB2 as a resistance driver by conducting in vivo limiting dilution transplants and drug treatments with Vemurafenib using established assays in zebrafish allografts (Dang et al., 2016;Heilmann et al., 2015) and showed MCR:SATB2 tumors to have increased tumor propagating potential ( Figure 5A-C), and primary resistance to MAPK inhibitors in vivo ( Figure 5D-F).
Analysis of publicly available genomic datasets showed SATB2 to not be recurrently altered in the TCGA dataset (Consortium, 2020), but its amplifications is observed in ~4-8% of patients in 3 independent datasets of metastatic melanoma (Supplementary Figure 6C) Van Allen et al., 2014), with a similar fraction of SATB2 high-expressing patients having a higher risk of metastasis-related poor outcome in two additional datasets (Supplementary Figure 6D).
In summary, our work identifies SATB2 as a driver of invasion and resistance to Vemurafenib treatment in melanoma. Yet, further studies will be needed to assess the prevalence of SATB2 expression changes in targeted therapy and immunotherapy clinical settings. Furthermore, we show in an endogenous transgenic tumor model that SATB2 transcriptional rewiring of melanoma towards a state similar to the less proliferative neural crest-like state described by Rambow and colleagues (Rambow et al., 2018)

CRISPR/Cas9 inactivation of satb2
To specifically inactivate satb2 in melanocytes, we engineered the MiniCoopR vector to express

Morpholino knockdown of satb2
To address whether SATB2 is necessary for neural crest and melanocyte development, 4 pg of a splicing morpholino against satb2 (GCAGTGTTGAACTCACCATGAGCCT, Ahn et al., 2010) was injected into one-cell stage TG(sox10:mCherry) embryos. At 3 dpf, injected and uninjected control embryos were scored for melanocyte and cranio-facial abnormalities using light and fluorescence microscopy. The experiment was repeated 3 times, with an average clutch size of 40 embryos per experiment.

Zebrafish primary tumor and cell line transplants
Pigmented primary melanomas were excised from euthanized MCR/MCR:EGFP and MCR:SATB2 zebrafish, and dissociated in 50% Ham's-12/DMEM medium containing 0.075 mg/mL liberase (Promega) for 30 minutes at RT, with periodical manual disaggregation using a razor blade. The dissociation medium was inactivated by the addition of 3 x 5 ml 50% Ham's-12/DMEM with 15% heat-inactivated Fetal Bovine Serum. To obtain a single cell suspension, cells were passaged through a 40 µm filter into a 50 mL falcon tube, and pelleted by centrifugation at 500 x g for 5 minutes. Cell numbers were determined, and a cell suspension of 100,000 cells/µl in PBS was made, and kept on ice. casper recipients received a split-dose irradiation of 15 Gy over two consecutive days prior to transplantation (Heilmann et al., 2015). During transplantation, irradiated casper recipients were anesthetized in 0.4% MS222, mounted in a moist sponge. The response rate for experimental cohorts was depicted via waterfall plots, and t-test statistics were applied for significance as previously described (Dang et al., 2016).

Whole-mount tumor immunohistochemistry
Zebrafish where euthanized on ice, melanomas were surgically removed, and fixed in either 4%

Western Blot
Adherent cells were scraped, and zebrafish tumors were mechanically homogenized in RIPA lysis buffer containing 1:100 protease inhibitors (P8340, Sigma-Aldrich) and 20 µM N-ethylmaleimide (Sigma-Aldrich). Lysates were incubated for 20 minutes on ice, and spun down for 10 minutes at 14,000 RPM at 4°C. Samples were denatured by adding Laemmli sample buffer (BioRad) with 5% β-mercaptoethanol (Sigma-Aldrich), and boiled at 95°C for 5 minutes prior to loading.

RNA extraction, quantitative RT-PCR analysis and RNA-seq sample preparation
Zebrafish tumors were isolated and mechanically homogenized in 350 μl RTL buffer (Qiagen) containing β-mercaptoethanol, on ice for 20 seconds. Adherent cells (zebrafish and human cell culture) were washed twice with ice cold PBS on ice, and cells were scraped in RTL buffer containing β-mercaptoethanol (Sigma-Aldrich). Tumor and cell lysates were next transferred onto a QiaShredder column (Qiagen). RNA isolation was performed using the RNA micro plus kit (Qiagen), according to the manufacturers instruction. RNA quality was determined using a Nanodrop. For RNA-seq on primary zebrafish melanomas, additional quality control of the total RNA was performed on and Fragment Analyzer. Total RNA was depleted of ribosomal RNA with the RiboZero gold kit (Epicentre), and enriched mRNA was applied to library preparation according to manufacturer's protocol (NEBNext Ultra). After repeated quality control for an average DNA input size of 300 base pairs (bp), samples were sequenced on a HiSeq Illumina sequencer with 2 × 100-bp paired-end reads. For qPCR, cDNA was synthesized with the SuperScript III Kit (Life Technologies) using a 1:1 mixture of random Hexamers and OligoDT.
Quantitative PCR was performed on a BioRad iQ5 real-time PCR machine using the Ssofast EvaGreen Supermix (BioRad). The ∆Ct or ∆∆Ct methods were used for relative quantification.
The average of 2 independent experiments is shown for SKMEL2 and a single experiment is shown for primary human melanocytes. qPCR primers were designed using GETprime (Gubelmann et al., 2011) or qPrimerDepot (https://primerdepot.nci.nih.gov). For primer sequences see Supplementary Table 5.

Zebrafish primary melanoma cell culture
Zebrafish primary melanoma cell lines were generated as described (Heilmann et al., 2015).

Matrix degradation assay
Oregon green 488-conjungated gelatin covered 12 mm coverslips, or glass bottom 6-well plates were coated as described . Human melanoma cell lines were induced with doxycycline for 48 hours prior to plating. Human and zebrafish cells were plated at a density of 30,000 cells/well in a 12 well plate, and grown on the fluorescent gelatin for 23-25 hours at 37°C (human cells) or 28°C (zebrafish cells). Cells were fixed in 4% PFA for 10 minutes at room temperature, washed 3 times with PBS, permeabilized with 0.4% Triton-PBS for 4 minutes, followed by 30 minute blocking (10% lab serum, 1% DMSO, 0.1% Tween in PBS). Coverslips were next stained with Alexa 650-or 568-conjungated Phalloidin (1:50) in blocking buffer, overnight at 4°C. Coverslips were washed extensively and nuclei were stained with DAPI, 10 minutes in PBS/0.1% Tween at RT. Coverslips were inversely mounted in Slow Fade Gold (Invitrogen) on a slide, and imaged on a Nikon C2si Laser Scanning Confocal, using ElementsX software.
Experiments were seeded in triplicate, in at least 3 independent experiments. Four to six high power (40-60x) areas per coverslip were imaged to determine the fraction of cells with degraded gelatin. An unpaired two-tailed t-test was used to compare significance between groups.

Time-lapse video of fluorescent gelatin degradation after SATB2 induction
Induced human melanoma cell lines were seeded onto an Oregon green 488-conjungated gelatin coated glass bottom 6-well plate, and allowed to adhere for 3 hours at 37°C, prior to imaging. Time-lapse movies were recorded on a Nikon Eclipse Ti Spinning Disk Confocal with a 10X objective (70um, 4um step size), for 16 hrs in 20 minute intervals at 37°C, 5% CO2. Images were processed using Photoshop, ImageJ, or Imaris.

Proliferation assays
Proliferation rates were determined using the CellTiter-Glo (Promega) luminescent cell viability assay, according to the manufacturer instructions. Cells were seeded at an initial density of 5,000 or 10,000 cells/well (96-well plate), in triplicate or quadruplicate wells, and experiments were repeated at least 3 independent times. An unpaired two-tailed t-test was used to compare significance between groups.

ChIP-seq tumor sample preparation and sequencing
ChIP-seq was performed as previously described (Lee et al., 2006). Zebrafish were euthanized on ice and melanomas were excised, finely minced using a scalpel blade in 5 ml cold PBS in a petri dish, transferred to a 50 ml Falcon tube, and the petri dish was rinsed once with 5 ml PBS.

ChIP-sequencing bioinformatics
All ChIP-Seq datasets were aligned to build version danRer7 of the zebrafish genome using Bowtie2 (version 2.2.1) (Langmead and Salzberg, 2012) with the following parameters: --end-toend, -N0, -L20. We used the MACS2 version 2.1.0 (Zhang et al., 2008) peak finding algorithm to identify regions of ChIP-Seq peaks, with a q-value threshold of enrichment of 0.05 for all datasets (Langmead and Salzberg, 2012). Genome track images were generated using the UCSC browser. SATB2-bound loci were determined by overlapping all regions 3 kb upstream of the transcriptional start site (TSS), gene body (GB) and 3 kb downstream of the transcriptional end site (TES) of all transcripts in the danRer7 assembly, with all significantly bound SATB2-peaks (P<10 -7 ). SATB2 bound peaks were compared to SATB2 enhancers peaks to determine the SATB2-bound enhancers. GREAT (Hiller et al., 2013) (2) fold change > 1.5 and a P-value < 0.05. Human orthologs were predicted using DIOPT33. Data sets are deposited to the GEO Gene Expression Omnibus, accession number GSE77923. RNA-seq was performed on three biological replicates of SKMEL2 transduced with pINDUCER20-SATB2 after 48 hours of induction with 2ug/mL of doxycycline as outlined above, except for the mechanical homogenization. Pathway analyses were conducted using GSEA (v4.0) (Subramanian et al., 2005) and IPA (Qiagen). Pathway and gene list overlaps were plotted using BioVenn (Hulsen et al., 2008).

Statistics
Comparison of Kaplan-Meier survival curves was performed by a log-rank (Mantel-Cox) test.
Statistical difference in qPCR expression analysis between large groups of MCR:EGFP and MCR:SATB2 primary tumors, and difference in estimated fraction of tumor propagating cells between MCR:EGFP and MCR:SATB2 donors were determined by an unpaired two-tailed Mann-Whitney test. The rest of the statistics were performed with an unpaired two-tailed t-test. An Ftest was used to determine similar variation between compared groups. Graphs show the median with s.e.m. No statistical methods were used to predetermine sample size. Experiments to quantify proliferation by PH3 immunohistochemistry, and cells with degraded 488conjungated gelatin were scored blindly. Irradiated casper zebrafish were randomized between transplantation groups. All statistical analyses were performed with GraphPad Prism. NS, not significant, P > 0.05; *P ≤ 0.05; **P ≤ 0.01; *** P ≤ 0.001; **** P ≤ 0.0001