Histone variant dictates fate biasing of neural crest cells to melanocyte lineage

Cell fate determination is driven by a cascade of transcription factors (TFs) that are governed by external cues (Kawakami and Fisher, 2011; Martik and Bronner, 2017). The neural crest sets a paradigm for stem cell multipotency and lineage specification, and serves as a model for understanding fate decisions. Recent data emerging from single-cell transcriptomics of neural crest cells (NCCs) demonstrates the emergence of fate bias or specification that further unfolds to culminate in fate commitment (Soldatov et al., 2019) or determination. NCCs and NCC-derived progenitor populations give rise to melanocytes during embryonic development and adult skin homeostasis and regeneration, respectively (White and Zon, 2008). The ability to trace cells visually, combined with the extensive literature on gene regulatory networks, render melanocytes as the preferred choice to trace the mechanisms behind fate decisions.

Almost all aspects of melanocyte biology, from lineage commitment to differentiation, are controlled by the central TF MITF (Goding and Meyskens, 2006; Kawakami and Fisher, 2017). WNT and BMP signalling (Jin et al., 2001; Takeda et al., 2000), and downstream TFs SOX10 (Marathe et al., 2017), PAX3 (Watanabe et al., 1998), FOXD3 (Kos et al., 2001) and LEF1 (Dunn et al., 2000), form the core of gene regulatory networks that control mitf induction and mediate melanocyte specification and determination. Although the upstream instructive signals and downstream orchestrators are well appreciated, cell-intrinsic modulators that mediate specification events remain obscure. Epigenetic mechanisms that alter chromatin accessibility could provide an intervening regulatory layer and are likely to play a key role in determining the fate of NCC lineages.

Epigenetic regulators, including the histone deacetylase HDAC1 (Ignatius et al., 2008, 2013), which is involved in locus-specific histone deacetylation; the DNA methyl transferase DNMT3b, which causes base methylation (Rai et al., 2010); and the histone variant H3.3, which replaces the core histone H3, are all involved in the derivation of specific cell types from NCCs (Cox et al., 2012). Another histone variant, H2A.Z, encoded by two genes, H2a.z.1 and H2a.z.2, has been implicated in NCC differentiation (Pünzeler et al., 2017; Sivasubbu et al., 2006). Recent studies demonstrate the crucial role of H2A.Z in human metastatic melanoma survival and proliferation (Vardabasso et al., 2015), and the involvement of the H2A.Z chaperone PWWP2A in NCC differentiation (Pünzeler et al., 2017). Another recent study demonstrates a selective role for H2A.Z.2 in craniofacial morphogenesis with implications for the congenital disorder Floating Harbor Syndrome (Greenberg et al., 2019). In addition, melanomas recapitulate neural crest fate (Kaufman et al., 2016). Therefore, it is possible that histone variants dictate the fate of melanocytes from neural crest-derived progenitors.

In this study, using zebrafish as well as mouse embryonic stem cell-derived melanocyte models, we demonstrate a selective role for H2A.Z.2 in the fate-biasing of NCCs to melanocytes. Furthermore, we show that this variant occupies the promoter of mitf and its upstream regulator sox9, and governs the inducibility of mitf under appropriate cues. CRISPR-targeted mutation in h2a.z.2 reinforces its central role during pigmentation and the related process of melanocyte regeneration in adult zebrafish. Thereby, we establish the selective role of a histone variant upstream of the TF network that enables cell-intrinsic responsiveness to external instructive signals. Our study adds a new dimension to the recent developments in tracing the lineage map of NCC derivatives and provides a molecular basis for the fate bias observed during specification.

Across vertebrates, the histone variant H2A.Z is encoded by two distinct paralogs, H2A.Z.1 (Z1) and H2A.Z.2 (Z2), that are around 97% identical at the protein level and are therefore difficult to distinguish experimentally (Dryhurst et al., 2009; Matsuda et al., 2010). In zebrafish, these two genes are named as h2az2b and h2az2a, respectively. In this study, for the sake of uniformity and simplicity, these are referred to as h2az.1 or Z1, and h2az.2 or Z2 across mice as well as zebrafish. The disruption of Z1 is not compensated for by intact Z2 and the knockout mouse embryos fail to survive (Faast et al., 2001). Therefore, to decipher their individual contribution to the differentiation of NCCs, especially to melanocytes, we adopted a morpholino (MO)-based transient knockdown strategy in developing zebrafish embryos. Knockdown of Z1 resulted in gross embryonic deformities and poor development being observed, as has also been reported previously (Greenberg et al., 2019; Madakashira et al., 2017; Murphy et al., 2018; Sivasubbu et al., 2006). Z2 morphants exhibited milder effects with regards to their gross morphology but they presented a drastic reduction in pigmentation (Fig. 1A-D).

The pigmentation phenotype of Z2 MO could be recapitulated using a splice block MO targeted to Z2, and also upon co-injection of Z2 MO with p53 MO, reiterating the specific nature of the phenotype (Fig. S1). Furthermore, analysis of recent single-cell sequencing studies (Saunders et al., 2019) confirmed the expression of Z1 and Z2 in lineages arising out of the neural crest that result in melanocytes (Fig. S2). Upon co-injection of Z2 MO along with mouse H2a.z.1 or H2a.z.2 mRNAs, a rescue of around 50% was observed only for H2a.z.2 mRNA (Fig. 1E). Thereby, it was demonstrated that Z2 has a distinct pigmentation role and is not compensated for by an increased dosage of Z1. To decipher the determinants that contribute to the Z2 selectivity, we performed single amino acid changes by site-directed mutagenesis of the three involved residues. Whereas H2AZ.2_A14T and H2AZ.2_A127V were able to rescue head melanophore numbers to a level comparable with the wild-type H2AZ.2, H2AZ.2_T38S was only capable of a partial rescue (Fig. S3).

The characteristic reduction in pigmentation of Z2 morphants could be either due to decreased pigment production or due to the absence of mature pigmented cells. At 2 days postfertilization (dpf), we observed a severe reduction in the number of ftyrp1:GFP⁺ cells in Z2 morphants, notably in the dorsal as well as the ventral pool of melanocytes above the yolk sac, clearly indicating that the number of melanocytes is affected (Fig. 1B). Additionally, at 7 dpf, the lateral line melanocytes were reduced in Z2 morphants, suggesting that the progenitor pool may be affected (Fig. 1C). We therefore used various transgenic reporter lines from the neural crest to the melanocyte lineage to estimate respective cell numbers at different stages of development (Fig. 1F). We observed that Z2 morphants displayed a 45-50% reduction in differentiating/differentiated melanophores and a 30% reduction in melanoblast numbers (Fig. 1F). However, the sox10⁺ progenitor population remained unaffected in Z2 morphants. Therefore, h2a.z.2 is likely to play a key role in the lineage commitment of melanocytes from neural crest-derived progenitors.

Multipotent NCCs progressively undergo the process of fate specification and determination to give rise to melanocytes, craniofacial cartilage, adipocytes, smooth muscle cells, Schwann cells and neurons (Donoghue et al., 2008; Le Douarin et al., 2004). Based on whole mount in situ hybridization (WISH), we observed that the early neural crest gene expression of foxd3, sox10, crestin and tfap2a were not decreased in Z2 morphants (Fig. 1G,H). This indicates that depletion of h2a.z.2 could selectively affect a subpopulation of NCC-derived progenitors from which melanocytes, and possibly a few other cell types, would be derived. Hence, we resorted to identifying gene expression changes in this seemingly homogeneous population that could explain decreased melanocyte numbers.

We created Z2 morphants in a Tg(sox10:GFP) line and sorted the GFP⁺ cells at 16-18 hpf. At this time point, neural crest-derived sox10 expression is retained and also lineage-specific differentiation markers begin to emerge and the cell fate identity is established (Curran et al., 2010; Wagner et al., 2018). Microarray and pathway enrichment analysis of the GFP⁺ cells revealed that the set of downregulated genes mapped to developmental pigmentation-related processes (Fig. 2A-C). Targeted analysis of the data indicated reduced expression of TFs with a known role in the melanocyte lineage, such as mitfa, tfap2a, tfap2e and sox9b, as well as melanocyte migration and survival factors, such as kita and ednrb1a. Similarly, we also observed decreased expression of glial lineage TFs such as nfatc4, pax7 and lmo4. However, the neuronal lineage TFs, such as foxo1, foxo3 and sox8, were upregulated in the Z2 morphants (Fig. 2D). As these are gene expression changes in the isolated sox10⁺ cells, we interpret the decreased expression to indicate alterations in the subpopulations of differentiating cells derived from the NCCs.

To verify the reduction in specific cell types, we adopted WISH for lineage-specific differentiation markers and analysis of cell type-specific transgenic lines that label NCC derivatives. dct and c-kit, which label melanocytes, showed decreased staining in the Z2 morphants (Fig. 2E). A similar pattern of decreased staining was observed for differentiated Schwann cells using the markers krox20 and mbp, which labelled the lateral line-associated glial cells (Fig. 2F). The number of undifferentiated Schwann cells, visualized by Tg(foxD3:GFP), was substantially reduced (Fig. 2F). Tg(nbt:DsRed)-marked enteric neurons remained unaffected (Fig. S4) and similarly WISH for isl2a showed no change (Fig. 2G). However, we observed the presence of ectopic neurons in the trunk, above the spinal cord neuron bundle, which suggested a local increase in the neuronal population (Fig. 2G and Fig. S4). Labelling by Tg(sox10:GFP), as well as Alcian Blue staining of craniofacial appendages and fin cartilage, showed no major changes, suggesting normal cartilage specification. However, we noticed a consistent decrease in jaw area as a result of defective ceratobranchial appendage spacing and a mild increase in the angle of the ceratobranchial cartilage (Fig. S4). Altogether, our data suggest that h2a.z.2 plays a selective role in the derivation of melanocytes and glial populations from the NCC. Interestingly, previous studies have shown that these two cell types emerge from a common progenitor (Raible and Eisen, 1994).

Based on temperature-sensitive mutants, Micropthalmia-associated TF mitfa, is seemingly dispensable for specification (Johnson et al., 2011). However, mitfa is the earliest known marker of lineage commitment that controls the expression of downstream identity genes involved in melanocyte fate determination (Goding, 2000; Lister et al., 1999). Indeed, we observed a decrease in mitfa expression in sox10⁺ cells, and this was further corroborated by the transgenic line mitfa:GFP and WISH staining for mitfa (Fig. 3A). Messenger RNA levels of mitfa were observed to be around 50% in the Z2 MO (Fig. 3B). Therefore, we set out to decipher the interplay of h2a.z.2 and mitfa in mediating melanocyte lineage commitment. A marked rescue in the melanocyte counts of mitfa mRNA-injected Z2 morphants was noted (Fig. 3C,D). To confirm the gene-selective elevation of mitfa by the two H2A.Z variants, we expressed the mouse orthologues by way of RNA injections. Western blot analysis indicated that the two variants exhibited comparable levels of expression but only Z2 enhanced the mitfa expression detected by WISH staining at 24 hpf (Fig. 3E,F). These observations confirmed the genetic interaction between h2a.z.2 and mitfa. Furthermore, h2a.z.2 likely functions upstream of mitfa and hence could modulate melanocyte fate.

Our experiments suggested that it is likely that h2a.z.2 plays a role in the specification/determination of an upstream progenitor population from which melanocytes arise. Additionally, the decrease in glia in the h2a.z.2 morphant indicated that the bipotent progenitor of melanocytes and glial populations could be affected. Furthermore, an increase in neurons is suggestive of a compensatory increase in neuronal progenitors. Our observations therefore indicate a hierarchical model of fate determination in the NCCs, which is in tune with recent studies (Le Douarin and Dupin, 2003; Raible and Eisen, 1994; Soldatov et al., 2019). We propose a model, wherein h2a.z.2 operates at the level of the tripotent glial/neuronal/melanocyte progenitor in specifying a melanocyte/glial fate (Fig. 3G).

Although the zebrafish-based experiments helped us to narrow down the function of this variant in melanocyte derivation, cellular and molecular insights are required to establish its role in specification. To understand this molecular mechanistic link, we established a method to derive melanocytes from R1/E mouse embryonic stem cells (mESCs), wherein the mESCs were converted to melanocytes through an intermediate embryoid body (EB) state (Fig. 4A). Herein, we observed ∼45% of the cells to be SOX10⁺ and around 25% to be MITF⁺, confirming the derivation of melanocytes in this model (Fig. S5). Targeted gene expression analysis of NCC, melanocyte and stemness-related genes were performed at different days of melanocyte formation (Fig. 4B). As anticipated, the expression of Mitf was induced at day 5 and further increased until day 16. Similarly, we observed that the expression of the pigmentation gene Dct was induced during the course of the experiment. We noticed that the expression of stemness factors Oct4 and Nanog were drastically reduced during the course of induction, suggesting that the cells are undergoing differentiation. Strikingly, the expression patterns of H2a.z.2 and H2a.z.1 were opposite to each other. Whereas the expression of H2a.z.2 was positively correlated with melanocyte specification, the expression of H2a.z.1 was decreased. This pattern of H2a.z.1 expression was similar to that of the stemness factors Oct4 and Nanog.

To illustrate the functional role of H2a.z.2 in mESCs, we silenced H2a.z.2 and studied melanocyte derivation. In the cells stably expressing short hairpin RNA that targets H2a.z.2 (shZ2), a reduction of ∼50% in H2a.z.2 mRNA levels could be observed, whereas H2a.z.1 levels remained unaltered (Fig. 4C). The total reduction in H2AZ protein was ∼50% upon Z2 silencing (Fig. 4D). Interestingly, previous studies report that knockdown of H2a.z.1 in embryonic stem cells leads to the formation of defective EBs with irregular structural integrity (Hu et al., 2013). We did not observe significant difference in the size and morphology of the H2a.z.2-depleted EBs (Fig. 4E), further strengthening the case for a differential role for the two proteins, at least during early developmental stages.

Upon the induction of melanocytes from R1/E cells stably expressing shZ2, the expression of the neural crest gene Sox10 remained unaltered when compared with nontargeting shRNA (shNT) (Fig. 4F). This observation is strikingly similar to the Z2 MO scenario in zebrafish embryos where a significant decrease in the mRNA levels of the master TF mitf and its target gene Tyr, in shZ2 cells, was noted when compared with shNT cells (Fig. 4G,H). Although the number of Sox10⁺ cells did not differ, we observed a 30% decrease in the Mitf⁺ cells upon Z2 silencing (Fig. 4I,J), thus, highlighting a conserved role for H2a.z.2 in melanocyte specification across zebrafish and mouse. Having established the cell-based model of melanocyte specification and the effect of H2a.z.2 in the derivation of melanocytes in this model, we set out to identify effectors and the mechanism by which H2a.z.2 could mediate melanocyte specification.

Both H2A.Z.1 and H2A.Z.2 occupy the nucleosomes in the promoter as well as enhancer regions immediate to the transcription start site (TSS) of the downstream genes. This occupancy can either facilitate (Bargaje et al., 2012; Gevry et al., 2009; Zovkic et al., 2014) or repress (Dai et al., 2017; Hardy et al., 2009) transcription in a context- and a model system-dependent manner. Partly, this discrepancy could be attributed to the near-identical proteins that are indistinguishable in chromatin immunoprecipitation (ChIP) studies and to the complex nature of chromatin changes brought about by H2A.Z. We therefore employed a silencing strategy followed by ChIP to address the molecular mechanism selective to Z2.

H2A.Z ChIP followed by qPCR for the promoter proximal nucleosomes of Mitf, and other crucial upstream effectors Sox9, Sox10 and Ckit, was performed across various stages of melanocyte specification in the mESC model. We then compared the nucleosome occupancy of H2A.Z to their mRNA level changes at these stages of differentiation. All of these four promoters showed increased occupancy during the course of melanocyte formation from mESCs. Interestingly, we observed that mRNA expression of Mitf, positively correlated with the occupancy of H2A.Z (H2A.Z.1+H2A.Z.2) in its promoter, suggesting that H2A.Z facilitates the induction of Mitf expression, possibly by modulating the accessibility of the Mitf promoter (Fig. 4K).

Although our progressive induction model provided a window to assess kinetic changes in promoter occupancy, heterogeneity in the cell populations arising during the course of melanocyte induction, limited our interpretations. Therefore, we employed a complementary approach wherein we utilized murine B16 melanoma cells that would represent specified melanocytes. ChIP with an H2A.Z antibody was carried out in B16 cells to check the occupancy in promoter regions of genes related to Mitf, and upstream regulators as well as genes downstream of it in melanocytes. Mitf and Sox9 promoters showed considerable enrichment of 5-7% over input (Fig. 5A), suggesting that these genes might be under H2A.Z control. In contrast, several pigmentation genes such as Dct, Tyr, Tyrp1 and Sox10, were not enriched for H2A.Z in B16 cells. This is in tune with previous studies on H2A.Z in human melanoma cells (Vardabasso et al., 2015) and our meta-analysis further corroborated this observation (Fig. S6).

Knockdown of H2a.z.2 in B16 cells with siRNA resulted in a consistent 70-80% decrease in mRNA levels of H2a.z.2, whereas H2a.z.1 levels remained unchanged (Fig. S7). H2A.Z ChIP-qPCR was carried out for Sox10, Sox9 and Mitf promoter regions in H2a.z.2 knockdown B16 cells (Fig. 5B,D,F). We observed decreased occupancy of H2A.Z upon Z2 knockdown. A correlated decrease in mRNA expression of Mitf and Sox9 (Fig. 5E,G) suggested that this was a direct consequence of decreased Z2 occupancy. In contrast, the Sox10 promoter, which had minimal occupancy in B16 cells, showed no change upon H2a.z.2 silencing. Correspondingly, the Sox10 mRNA levels remained unaltered (Fig. 5C). To probe the role of Z2 in the dynamic regulation of Mitf expression, we performed H2AZ ChIP of the Mitf promoter. Upon α-MSH stimulation, which increased mitf RNA levels approximately twofold, an increase in H2AZ occupancy of ∼3.5-fold was observed (Fig. 5H,I). Furthermore, Mitf promoter activity was downregulated upon silencing H2A.Z.2 under basal and α-MSH inducible conditions (Fig. 5J). The effect of H2A.Z2 occupancy on the Mitf promoter was addressed in a reporter assay with the overexpression of H2a.z.2 in B16 cells. Although there was no change in the Mitf promoter activity under basal conditions in Z2-overexpressing cells, α-MSH stimulation resulted in a sevenfold upregulation, which was much higher than the empty vector control transfected cells (Fig. 5K,L). A probable confounder could be the inefficient chromatinization of the reporter construct. To address this, we generated a stable line in B16 cells. Overexpression of Z1 and Z2 followed by stimulation with α-MSH confirmed a Z2-selective inducibility (Fig. S7). Taking these results into account, we conclude that Z2 modulates the promoter region of MITF and facilitates its induction, the absence of which leads to decreased Mitf expression resulting in a reduced number of specified melanocytes.

Occupancy of H2A.Z.2 renders the chromatin state of the promoter to be responsive to upstream activating signals. The role of H2A.Z.2 delineated in the current study, shows a striking similarity with its role in the brain, whereby the incorporation of H2A.Z mediates dynamic gene expression changes associated with memory formation (Zovkic et al., 2014). It is also interesting to note that upon Z2 silencing in the zebrafish system, sox9a and sox9b are reduced in the NCC populations (Fig. 2D). Given the upstream role of sox9 as a positive regulator of mitf (Passeron et al., 2007), it is likely it may play an additional role in melanocyte specification via mitf. Therefore, we have identified a role of the histone variant H2A.Z.2 in the specification of melanocytes from the multipotent NCCs via the activation of a key transcriptional programme downstream of Mitf (Fig. 5M).

Targeted knockout of H2a.z.1 has been attempted in several organisms and found to be embryonically lethal. However, the Z2 variant has not been systematically studied. CRISPR/Cas9-based targeted deletion of h2a.z.2 in zebrafish resulted in severe lethality and gross deformities (Fig. S8). Genotyping of the h2a.z.2 locus from these embryos revealed several indels resulting in large-scale disruptions of this gene. Although h2a.z.2 seems to be essential for early development, we isolated h2a.z.2 mutant lines by targeting the coding region downstream of the core histone domain, which maps to the C-terminal extension region beyond the αC domain (Fig. 6A). This region is unstructured and is involved in dynamic docking with histone H3/H4 during nucleosome formation (Bönisch et al., 2012). Interestingly, in the mutant line (H2a.z.2^Gln125fsX144) we observed fewer melanocytes at the embryonic stage 48 hpf (Fig. 6B). Decreased pigmentation could be traced to the juvenile stages as well (Fig. 6C), confirming the partial loss of function and dominant nature of the mutation. Therefore, the CRISPR-based targeted mutagenesis of H2A.Z.2 confirmed the role of this essential protein in specifying melanocyte numbers.

In the adult zebrafish, stripes appeared paler and on closer examination a reduction in melanophore numbers was observed (Fig. 6D). In the adult animal, the role of H2A.Z.2 in specification could be further verified by studying melanocyte regeneration, a process wherein differentiated melanocytes arise from a prespecified progenitor pool and can be experimentally studied using fin-clip experiments (Hultman et al., 2009). The adult zebrafish fin was therefore cut and allowed to regenerate for a week, during which time the quiescent progenitors were stimulated, and underwent the process of differentiation to generate pigmented cells. We noticed that the h2a.z.2 mutant animals demonstrated comparable fin growth; however, they displayed decreased melanocyte regeneration kinetics (Fig. 6E,F). The level of achieved repigmentation is in tune with the severity of the ventral stripe pigmentation observed in these mutants (Fig. 6G). The NCC-derived fate-specified progenitors established in the skin and responded to melanocyte loss by providing requisite melanocytes. The decreased recall of melanocytes indicated that the progenitors were inadequately populated, and endorsed the role of h2a.z.2 in melanocyte fate specification. Therefore, using targeted genome engineering, we unequivocally demonstrated the role of H2A.Z.2 in establishing a melanocyte progenitor population, which was evident not only at the embryonic stage but also during adult melanocyte regeneration.

The decoding of cell fate decisions has emerged as an intense area of research, primarily because of the immense application potential it has in futuristic regenerative interventions. Fate establishment is perceived to be a continuum, wherein progressive fate bias culminates in determination and is controlled by appropriate external signals. Currently, we understand deterministic events better because of the characteristic gene expression signatures orchestrated by key TFs. However, the underlying molecular mechanisms behind the fate-biasing of specification events remains unexplored. The nondeterministic nature of fate bias induced during specification necessitates a modulatory role for the specifier for which epigenetic factors are perfectly suited. Histone modifications and the deposition of variants could function in conjunction with TFs and offer a greater plasticity to modulate the cell fate. A recent study identifies a selective role for H2A.Z.2 in craniofacial cartilage formation in Xenopus (Greenberg et al., 2019). Herein, the H2A.Z.2 variant is involved in this process but H2A.Z.1 is not. Our observations confirm the conservation of this function in zebrafish.

In this study, we establish the occupancy of H2A.Z.2 on the mitf promoter that facilitates its expression by upstream activating signals. We propose a model for melanocyte lineage specification, wherein H2A.Z.2 acts as a key specifying signal. Three of the fate-biasing transcriptional programmes in the neural crest lineage are involved in suppressing the melanocyte fate. These include Sox2 (Adameyko et al., 2012) and Neurog2 (Soldatov et al., 2019) from the neural lineage, and Foxd3 (Ignatius et al., 2008; Kos et al., 2001) from the glial lineage, all of which actively suppress Mitf expression. In contrast, H2a.z.2 is a positive regulator that promotes melanocyte fate bias, as its overexpression facilitates Mitf induction. A hierarchical fate-biasing model necessitates such opposing programmes to ensure defined populations of emerging differentiated cells. The observed concomitant decrease in the glial population suggests that H2A.Z.2 could operate as a specifier at the tripotent melanocyte/glial/neural progenitor level, driving the allocation towards a melanocyte-glial bipotent fate. The regulation of Sox9 by H2a.z.2 could have a dual role of promoting Mitf expression in melanocytes and also facilitating the glial-specific programmes (Cheung and Briscoe, 2003).

Regenerative melanocytes arise from NCC-derived stem cells that are positioned along the lateral line in close juxtaposition with dorsal root ganglia in zebrafish embryos (Dooley et al., 2013). In the adult zebrafish, precursor cells directly differentiate into mature melanocytes and also divide to yield additional lineage-restricted cells (Iyengar et al., 2015). In mammalian skin, this population resembles nerve-associated melanocyte precursors (Adameyko et al., 2009) and the hair follicle bulge region progenitors that give rise to melanocytes at every hair cycle (Nishimura et al., 2002). These progenitors are thought to restore epidermal pigmentation after wound healing and during melanocyte regeneration in conditions such as vitiligo (Chou et al., 2013).

The establishment of melanocyte stem cells is thought to occur independently of Mitf (Johnson et al., 2011). However, all the downstream events involved in melanocyte fate determination are orchestrated by this master regulator. Hence the ‘melanocyte specification factor’ is likely to function upstream of Mitf. The transient epigenetic memory conferred by the incorporation of H2A.Z.2 at the Mitf promoter would not only facilitate its expression upon external cues such as WNT ligands, it would also enable the retention of this memory in transiently amplifying cells. The Mitf-inducing Wnt-signalling pathway triggers melanocyte regeneration by possibly activating this H2A.Z.2 occupied state of the Mitf promoter. Therefore, the histone variant-based specification elucidated in our study enables the retention of plasticity, a hallmark of the specified state. The identification of such upstream specification mechanisms would enable futuristic regenerative therapeutic approaches from induced pluripotent cells, by employing a selective activation of lineage determinants.

Fish experiments were performed in strict accordance with the institutional animal ethics approval (IAEC) of the CSIR-Institute of Genomics and Integrative Biology (IGIB), India (Proposal No 45a). All efforts were made to minimize animal suffering.

Zebrafish of the line ASWT were bred, raised and maintained at 28.5°C according to standard protocols (Westerfield, 2000) and were housed at the CSIR-IGIB, Mathura Road New Delhi, India. Embryos were staged on the basis of time [hours post fertilization (hpf); days post fertilization (dpf)] and morphological features, according to Kimmel et al., 1995. Embryos older than 24 hpf were treated with 0.003% 1-phenyl-2-thiourea (PTU) to inhibit pigment formation, aiding fluorescent imaging and RNA in situ hybridization analysis, as and when experiments demanded the depigmentation of the animals. The plasmid ftyrp1:GFP was a kind gift from Dr Xiangyun Wei (University of Pennsylvania, USA) (Zou et al., 2006) and the construct was injected at a concentration of 10-20 pg along with 50-75 pg of Tol2 transposase mRNA into one-cell-stage zebrafish embryos to create the transgenic line at CSIR-IGIB. Founder lines were established and propagated, which showed an expression pattern similar to pt101 as reported by Zou et al. (2006). The details of the zebrafish lines used in this study are provided in Table S1.

Single-cell transcriptomic data for 5 dpf zebrafish sox10+ cells was downloaded from GEO series GSE131136 (Saunders et al., 2019). The analysis was replicated as per their methods using Monocle v.2.99.1 (Qiu et al., 2017a,b; Trapnell et al., 2014) in RStudio v.3.6 (www.rstudio.com/) and the expression of h2afva (H2a.z.2) and h2afvb (H2a.z.1) in different cell types was plotted.

MOs for blocking translational initiation or splicing of RNA were designed by Gene Tools (Table S1). The minimum effective concentration was determined for each MO (Kaufman et al., 2016) using dosage titration experiments. Concentrations of various MOs used in this study as follows: 3-3.5 ng of h2a.z.2 translation block (Z2 MO), 4-4.5 ng of h2a.z.2 splice block (Z2 SB MO), 2 ng of h2a.z.1 splice block (Z1 MO) (Sivasubbu et al., 2006) and 0.8 ng of double-targeting (h2a.z.1 and h2a.z.2) translational block (Z1+Z2 MO). A standard control MO from Gene Tools, LLC was used as a control and dosed accordingly. The p53 MO was co-injected at a 1:1 ratio with the corresponding gene targeting MOs.

For rescue experiments, mouse H2a.z.1 (Z1-RNA) and H2a.z.2 (Z2-RNA) genes were amplified using primers listed in Table S1. Similarly, the zebrafish mitfa gene was amplified using PCR (Table S1) and cloned into the Zero Blunt TOPO vector (Thermo Scientific, K287540) according to the manufacturer's protocols. Capped and PolyA tail RNA was synthesized by in vitro transcription using the T7-Ultra mRNA Synthesis Kit (Thermo Scientific, AM1345) followed by PolyA tailing according to the manufacturer's protocol. A volume of 100 pg (for H2a.z.1 and H2a.z.2) or 50 pg (for mitfa) was co-injected into one-cell-stage embryos along with the Z2 morpholino. Site-directed mutagenesis for mH2a.z.2 A14T, T38S and A127V was carried out using the SDM II kit (Agilent) using primers listed in Table S1 according to the manufacturer's protocols. Subsequently, the mutations were verified by Sanger sequencing and cloned into the Zero Blunt TOPO vector for capped and PolyA tail RNA synthesis.

RNA in situ hybridization was performed as described by Thisse and Thisse (2008). Riboprobes of the following genes were used for this study: foxd3, sox10, crestin, krox20, sox9a (Stewart et al., 2006), mbp (Brösamle and Halpern, 2002), neuroD (van der Velden et al., 2013), isl2a (Appel et al., 1995), sox9b (Yan et al., 2005), tfap2a, mitfa, c-kit, dct (Van Otterloo et al., 2010). Zebrafish h2a.z.2 exon5-3′UTR (307 bp) was amplified using primers listed in Table S1. The product was subsequently cloned into pCR4-TOPO (Thermo Scientific, K287540) and was used as template for in vitro transcription of antisense probe using Megascript T7 (Thermo Scientific, AM1334).

Bright-field imaging for live imaging and gene expression pattern (WISH) studies of zebrafish embryos was performed using a Zeiss stereo microscope (Stemi 2000-C). A Zeiss Axioscope A1 Microscope (with AxiocamHRc) was used for fluorescence imaging. Zeiss proprietary software was used to capture images, which were processed and analysed in Adobe Photoshop CS3.

Cell counts for sox10:EGFP, mitfa:GFP, ftyrp1:GFP⁺ cells and ASWT embryos (for melanophore counts) were performed using an AMNIS ImageStream MKII flow cytometer. Briefly, embryos were dechorinated using pronase (5 mg/ml) (Sigma-Aldrich, P8811) for 10-15 min and collected in microfuge tubes. The embryos were deyolked in ice cold Ringer's solution using a micropipette tip and spun at 100 g for 2 min in a tabletop centrifuge at 4°C (Eppendorf, 5418R). The supernatant was discarded and the embryo bodies were trypsinized using TrypLE Express (Thermo Scientific, 12604039) for 15 or 30 min for 24 hpf or 48 hpf, embryos respectively, at room temperature. The cell suspension was passed through a 70 μm cell strainer and washed twice with ice cold phosphate-buffered saline. The cell suspension was analysed in the imaging flow cytometry system and at least 50,000 images were captured and processed using IDEAS (AMNIS) software.

Control MO and Z2 MO were injected into one-cell-stage Tg(sox10:GFP) embryos and allowed to develop until 15 hpf. The embryos were processed in a similar manner as that for imaging flow cytometry analysis. Finally, the trypsinized cells were resuspended in ice cold PBS+10% fetal bovine serum and subjected to FACS (BD Biosciences FACSAria™ III). The GFP⁺ cells were sorted and collected in RP1 buffer and RNA was isolated using Nucleospin RNA XS kit (Macherey Nagel, 740902) according to the manufacturer's protocol. The microarray experiment and subsequent analysis was outsourced to Genotypic Technology (Bengaluru, India).

B16 mouse melanoma cells were cultured in DMEM-High glucose media (Sigma-Aldrich, D5648) supplemented with 10% FBS (Thermo Scientific) at 5% CO₂ (Eppendorf, New Brunswick Galaxy 170S incubator). Media was changed every day and cells were passaged upon reaching 80% confluence. Mouse R1 embryonic stem cells (R1/E) were cultured in DMEM Glutamax media supplemented with sodium pyruvate (Thermo Scientific,10569 010), MEM-NEAA (Thermo Scientific, 11140 050), 2-mercaptoethanol (Thermo Scientific, 31350 010) and 20% Pansera (Pan Biotech, P30 2602) at 5% CO₂. Media was changed every day and cells were passaged every other day.

Melanocytes were generated from R1/E mESCs using a modified protocol from Yang et al. (2011). For EB formation, cells were trypsinized and 50,000 cells/ml suspensions were prepared. Subsequently, 10 μl of this suspension was placed in the lid of a 90 mm Petri dish and inverted upon the base of the Petri dish filled with 15 ml of autoclaved MilliQ H₂O. After 4 days, the EBs were transferred to 24-well or 6-well plates coated with collagen (Thermo Scientific, A10483 01) containing melanocyte conversion media (45% DMEM – High glucose, 45% Reconstituted M254) (Thermo Scientific, M254CF), 10% FBS, 50 ng/ml SCF (Peprotech, 300 07 10), 100 nM Endothelin-3 (Sigma-Aldrich, E9137), 20 pM Cholera toxin (Sigma-Aldrich, C8052), 0.5 μM Dexamethasone (Sigma-Aldrich, D1756), 50 ng/ml WNT3a (Peprotech, 315 20 10), 4 ng/ml beta-FGF (Thermo Scientific, RFGFB50), 50 nM PMA (Sigma-Aldrich, P1585), 1 N2 Supplement (Thermo Scientific, 17502 048) and 1× Antibiotic-Antimycotic (Thermo Scientific, 15240 052). Cells were collected at specific time points for downstream analysis.

Cells were plated on coverslips, washed twice with 0.01% Tween 20 in PBS, fixed with 4% paraformaldehyde for 10 min, permeabilized with 0.1% Triton X-100 in PBS for 10 min and washed three times with 0.01% Tween 20 in PBS. Coverslips were then incubated in blocking buffer (5% normal goat serum in 0.01% Tween 20) for 1 h followed by washing and overnight incubation with mouse anti-MITF (1:100, C-5, Abcam, ab12039) or rabbit anti-SOX10 (1:100, Abcam, ab155279) antibody at 4°C. This was followed by three washes with 0.01% Tween 20 and then incubation with Alexa Fluor 594 goat anti-mouse (for MITF) or anti-rabbit (for SOX10) dye for 1 h at room temperature. Coverslips were washed three times with 0.01% Tween 20 and then mounted in ProLong Gold Antifade Mountant with DAPI (Molecular Probes, Thermo Scientific). Images were visualized using a DeltaVision microscope (GE Healthcare).

2×10⁵ B16 melanoma cells were seeded in 6-well plates before the day of transfection. The following day, 100 nM of nontargeting control siRNA, H2a.z.2 siRNA (siRNA pool of five siRNAs) (Dharmacon ON-TARGETplus; L-063612-01-0005) or H2a.z.1 siRNA (siRNA pool of five siRNAs) (Dharmacon ON-TARGETplus; L-042994-01-0005) was transfected according to the manufacturer's protocol using Dharmafect (Dharmacon, T-2001) transfection reagent. Cells were collected at specified time points for downstream analysis. Plasmid transfections were performed using Lipofectamine 2000 (Thermo Scientific, 11668019).

RNA was isolated from zebrafish embryos,B16 melanoma cells or R1/E cells using Nucleospin Triprep (Macherey Nagel, 740966) according to the manufacturer's protocols. cDNA was synthesized using the Superscript III First-Strand cDNA Synthesis Kit (Thermo Scientific, 1800051). Quantitative real-time PCR was carried out as described (Pfaffl, 2001) using a Roche Lightcycler Real-Time PCR System. For quantification, the relative standard curve method was used (as described by the manufacturer) to generate raw values representing arbitrary units of RNA transcripts. Data analysis was performed using the 2^−ΔΔCT method.

Lentiviral packaging and transduction was performed as previously reported by Motiani et al. (2018). Briefly, lentiviral packaging was performed by co-transfecting pVSVG, pdR8.2 and H2a.z.2 shRNA, or nontargeting shRNA in HEK cells. Two days post-transfections, lentiviral particles were collected from the cell supernatant. Next, 1:100 HEPES pH 7.25 1 M stock solution and 1:100 Polybrene 0.4 mg/ml were added to the lentiviral supernatant. Embryonic stem cell media [ESCM; DMEM GlutaMax media supplemented with sodium pyruvate (Thermo Scientific; 11140-50), β-mercaptoethanol (Thermo Scientific; 31350-010) and 20% Pansera ES (Pan Biotech; P30-2602)] and lentiviral cocktail were mixed in a 1:1 ratio and added to R1/E ESCs. The mixture was spun for 30 min at 37°C and at 1076 g. The cells were incubated at 37°C for 10-12 h and then the media was changed to fresh ESCM. Cells were checked for GFP fluorescence at 48 h after transduction in R1/E cells and were selected using 0.5 µg/ml puromycin for a week. Experiments were performed with FACS-sorted GFP⁺ cells propagated as a pool.

ChIP assays were performed according to the protocol provided by Upstate Biotechnology with modifications as suggested in the Fast ChIP protocol. The anti-H2A.Z antibody (Abcam, ab4174) was used for the assay. Anti-rabbit IgG (Thermo Scientific, 02-6102) was used for isotype control in all of the cell lines. Briefly, B16 melanoma cells or R1/E ESCs were fixed with 10% formalin (Sigma-Aldrich, HT501128) and incubated at 37°C for 10 min. Next, 2.5 M glycine (Sigma-Aldrich, 50046) was added to the cells and again incubated at 37°C for 10 min. Cells were washed with ice cold 1× PBS containing protease inhibitors. Cells were then scraped and centrifuged at 86 g for 5 min at 4°C. The cell pellet was lysed in SDS lysis buffer [1% SDS (Sigma-Aldrich, L3771), 10 mM EDTA (Sigma-Aldrich, E6758), 50 mM TRIS (Sigma-Aldrich, T6066) (pH 8.1)] on ice for 30 min. The cells were then sonicated in a Bioruptor (Diagenode) in ice. The chromatin lysate was then estimated for protein content using a BCA kit (Thermo Scientific, 23225). The samples were decrosslinked overnight at 65°C and deproteinized with Proteinase K (Sigma-Aldrich P4850). Sheared DNA was column-purified using a PCR Purification Kit (Qiagen, 28104) and estimated using a QUBIT dsDNA HS Kit (Thermo Scientific, Q32851). ChIP was performed using 3-5 μg of the respective antibody added to an equal amount of chromatin across samples and incubated overnight at 4°C (10% of chromatin was kept separately as input). Next day, 100 μl of Protein A agarose beads (G-Bioscience, 786 283) were added to the mixture of chromatin and antibody, and incubated for 4-6 h at 4°C. After incubation the beads were washed twice with low salt buffer [0.1% SDS, 1% Triton X-100 (Sigma-Aldrich, T8787), 2 mM EDTA, 20 mM Tris HCl (pH 8), 150 mM NaCl (Sigma-Aldrich, S3014)], high salt buffer [0.1% SDS, 1% Triton X-100, 2 mM EDTA, 20 mM Tris HCl (pH 8), 500 mM NaCl] and LiCl buffer [0.25 M LiCl (Sigma-Aldrich, D62476), 1% Igepal CA-630 (Sigma-Aldrich, I8896), 1 mM EDTA, 10 mM Tris HCl (pH 8), 1% deoxycholate (Sigma-Aldrich, D6750)]. Finally, the beads were washed with Tris EDTA buffer. The beads were then incubated with 100 µl elution buffer [1% SDS, 0.75% sodium bicarbonate (Sigma-Aldrich, S5761)] and 1 µl of 20 mg/ml proteinase K for 15 min; elution was performed again with 100 µl for another 15 min. Subsequently, the samples were kept for 12-14 h at 65°C for reverse crosslinking. After incubation, the samples were column-purified using the PCR Purification Kit; the inputs from each of the samples were also included in the purification step. SYBR green- (Kapa Biosystems) based qRT-PCR was set up using eluted DNA and graphs were plotted as percentage input. Primer sequences for target genes are provided in Table S1.

Zebrafish embryos were dechorinated, deyolked and made into a single-cell suspension. Cell lysis was performed by adding 1 ml hypotonic lysis buffer to the cell pellet and gently mixing by pipetting to dissolve the cell pellet. The cells were incubated on a rotator at 4°C for 30 min and then centrifuged at 10,000 g for 10 min at 4°C. Nuclear lysis was performed by removing the supernatant and 400 µl of 0.4 N H₂SO₄ was added to the cells. The pellet was dissolved without any clumps and incubated at 4°C in a rotator overnight at 10 rpm. The samples were then centrifuged at 16,000 g for 10 min at 4°C. The supernatant containing histones was transferred to a fresh epitube. Histone precipitation was performed by adding TCA (trichloroacetic acid) (stock 66%; final concentration 33%) drop by drop to the tube while vortexing constantly at slow speed. The sample was then incubated on ice for 30 min until the solution turned milky white. The samples were centrifuged at 16,000 g for 10 min at 4°C. The supernatant was removed and the pellet was washed with 500 µl of ice cold acetone. The samples were centrifuged again at 16,000 g for 5 min at 4°C. The acetone wash step was repeated again. The pellet was air dried and then dissolved in nuclease-free water. Histone quantification was performed by running 10 µl of samples on a 12% SDS gel and stained with Coommassie Brilliant Blue R-250 overnight. The gel was scanned after destaining and densitometry was performed using ImageJ software. According to the ImageJ analysis, the sample volumes were recalculated for equal loading to perform western blotting.

Mouse mitf-luciferase (a kind gift from Dr Krishnamurthy Natarajan, Jawaharal Nehru University, India) and sox9-luciferase (a kind gift from Dr Peter Koopman, The University of Queensland, Australia) constructs were co-transfected (at a ratio of 1:10) with Renilla luciferase (pRL-TK vector, Promega) into nontargeting control- or H2a.z.2 siRNA-treated B16 melanoma cells, and B16 cells transfected with pIRES empty vector or pIRES-mH2a.z.2 vector. After 24 h the cells were lysed with passive lysis buffer and luciferase assays were performed according to the manufacturer's protocols (Dual Luciferase Reporter Assay System; Promega). Alpha MSH (Sigma-Aldrich, M4135) treatment was provided 6 h post-transfection of the luciferase constructs. Stable B16 cells containing the Mitf promoter luciferase were generated by transfecting the construct followed by selection with G418 for a week. This chromatinized reporter line was subsequently used for downstream experiments.

sgRNA targeting the zebrafish h2a.z.2 gene was selected using the ECRISP (www.e-crisp.org/E-CRISP/) online tool using default parameters. Primers were designed for generating the complete sgRNA (related to STAR methods: oligonucleotides) using annealing PCR. In vitro transcription was performed on these PCR products using T7 MEGAshortscript Transcription Kit (Thermo Scientific) according to the manufacturer's protocols. Next, 100 pg of sgRNA was injected along with 500 pg of spCAS9 protein (kind gift from Dr Debojyoti Chakraborty, CSIR-IGIB, India). The F0 embryos were grown to adulthood and were then outcrossed with ASWT to give rise to F1 animals. In the F1 generation, genomic DNA was isolated from fin clips of putative mutants and the target region was amplified using primers provided (related to STAR methods: oligonucleotides). The PCR products were subjected to Sanger sequencing to confirm mutations.

We used the F2 generation h2a.z.2 mutants, namely h2a.z.2^Gln125fsX144, for this study. In the H2a.z.2^Gln125fsX144 mutants we observed a 10-base deletion in exon 5, which leads to skipping of the stop codon and adds an additional 14 amino acids to the H2a.z.2 protein sequence.

Fin amputation experiments were performed on heterozygous F1 mutant animals. Mature zebrafish adults (5-10 months) were anaesthetized and the distal two-thirds of the caudal fin was amputated with a scalpel. The amputated fin was imaged immediately using a stereozoom microscope (Zeiss Stemi 2000-C). The fish were then returned to fresh water at 25°C on a regular feeding schedule for the duration of the experiment. The fin was allowed to regenerate for a week, during which the quiescent melanocyte stem cells are stimulated and undergo the process of differentiation to generate pigmented cells. The tail fin of the fish was imaged again after 7 days and melanocyte regeneration was assessed.

A Student's t-test was performed to determine statistical significance of the data. Asterisks on the error bars correspond to *P≤0.05, **P≤0.01, ***P≤0.001, ****P≤0.0001 and ns (not significant; P>0.05). Graphs were plotted using GraphPad Prism.

We acknowledge Dr C. Sachidanandan (CSIR-IGIB, New Delhi, India) for providing several zebrafish lines and WISH probes; and B. Sharma for assistance with WISH experiments. We acknowledge the Imaging and FACS facility of CSIR-IGIB.