We present cisTopic, a probabilistic framework to simultaneously discover co-accessible enhancers and stable cell states from sparse single-cell epigenomics data (http://github.com/aertslab/cistopic). On a compendium of single-cell ATAC-seq datasets from differentiating hematopoietic cells, brain, and transcription-factor perturbations, we demonstrate that topic modelling can be exploited for a robust identification of cell types, enhancers, and relevant transcription factors. cisTopic provides insight into the mechanisms underlying regulatory heterogeneity within cell populations.

Melanoma cells can switch between a melanocytic and mesenchymal-like state. Scattered evidence indicates that additional, intermediate state(s) may exist. To search for such states and decipher their underlying gene regulatory network (GRN), we studied ten melanoma cultures by single-cell RNA-seq, and 26 additional cultures by bulk RNA-seq. Although each culture exhibited a unique transcriptome, we identified shared GRNs that underlie the extreme melanocytic and mesenchymal states, and the intermediate state. This intermediate state is corroborated by a distinct chromatin landscape and governed by the transcription factors SOX6, NFATC2, EGR3, ELF1 and ETV4. Single-cell migration assays confirmed its intermediate migratory phenotype. By time-series sampling of single cells after knockdown of SOX10, we unravelled the sequential and recurrent arrangement of GRNs during phenotype switching. Jointly, these analyses indicate that an intermediate state exists and is driven by a distinct and stable "mixed" GRN rather than being a symbiotic, heterogeneous mix of cells.

Single-cell techniques are advancing rapidly and are yielding unprecedented insight into cellular heterogeneity. Mapping the gene regulatory networks (GRNs) underlying cell states provides attractive opportunities to mechanistically understand this heterogeneity. In this review, we discuss recently emerging methods to map GRNs from single-cell transcriptomics data, tackling the challenge of increased noise levels and data sparsity compared with bulk data, alongside increasing data volumes. Next, we discuss how new techniques for single-cell epigenomics, such as single-cell ATAC-seq and single-cell DNA methylation profiling, can be used to decipher gene regulatory programmes. We finally look forward to the application of single-cell multi-omics and perturbation techniques that will likely play important roles for GRN inference in the future.

Deciphering the genomic regulatory code of enhancers is a key challenge in biology as this code underlies cellular identity. A better understanding of how enhancers work will improve the interpretation of non-coding genome variation, and empower the generation of cell type specific drivers for gene therapy. Here we explore the combination of deep learning and cross-species chromatin accessibility profiling to build explainable enhancer models. We apply this strategy to decipher the enhancer code in melanoma, a relevant case study due to the presence of distinct melanoma cell states. We trained and validated a deep learning model, called DeepMEL, using chromatin accessibility data of 26 melanoma samples across six different species. We demonstrate the accuracy of DeepMEL predictions on the CAGI5 challenge, where it significantly outperforms existing models on the melanoma enhancer of IRF4. Next, we exploit DeepMEL to analyse enhancer architectures and identify accurate transcription factor binding sites for the core regulatory complexes in the two different melanoma states, with distinct roles for each transcription factor, in terms of nucleosome displacement

Chromatin accessibility refers to the level of physical compaction of chromatin, a complex formed by DNA and associated proteins consisting mainly of histones, transcription factors (TFs), chromatin-modifying enzymes and chromatin-remodelling complexes 1-3. Although eukaryotic genomes are generally packed into nucleosomes, which comprise ~147 bp of DNA wrapped around an octamer of histones 4,5 , nucleosome occupancy is not uniform in the genome, and varies across tissues and cell types. Nucleosomes are typically depleted at genomic locations that represent cis-regulatory elementsenhancers and promoters, among others-that interact with transcriptional regulators (for example, TFs), resulting in accessible chromatin 6-10. Profiling chromatin accessibility on a genome-wide scale is an excellent tool to map putative regulatory elements in a cell type or cell state. Post-translational chemical modifications of chromatin, including DNA methylation (in vertebrates) and histone methylation and acetylation, are dynamic and change between different cell states, similar to nucleosome positioning. These post-translational modifications are often correlated with chromatin accessibility and can reflect specific functionalities of genomic regions related to the regulation of gene expression 11,12. Changes in these post-translational modifications, such as increased or decreased histone methylation and acetylation, are affected by a large set of chromatin-modifying enzymes that can be recruited to chromatin regions by TFs. These modifications alter the physico-chemical properties of the chromatin, which in turn can influence the formation of transcriptional condensates 13,14. In addition, active chromatin remodelling impacts nucleosome occupancy; for example, the SWI/SNF complexes use ATP hydrolysis to alter histone-DNA contacts, thereby repositioning or removing nucleosomes 15. Dynamic changes in the chromatin structure, chemical modifications and nucleosome positioning form a crucial interplay with the TFs that drive differentiation of cells during development 16,17. Initial changes in chromatin accessibility are caused by the binding of TFs, which outcompete histones and recruit cofactors, including ATP-dependent chromatin remodellers 18,19 , or by TFs that preferentially bind to their recognition sequence in nucleosomal DNA 20,21. The binding of these initial TFs, known as pioneer factors, can recruit other TFs to co-bind and further stabilize the nucleosome-depleted

Genomic sequence variation within enhancers and promoters can have a significant impact on the cellular state and phenotype. However, sifting through the millions of candidate variants in a personal genome or a cancer genome, to identify those that impact cis-regulatory function, remains a major challenge. Interpretation of non-coding genome variation benefits from explainable artificial intelligence to predict and interpret the impact of a mutation on gene regulation. Here we generate phased whole genomes with matched chromatin accessibility, histone modifications, and gene expression for 10 melanoma cell lines. We find that training a specialized deep learning model, called DeepMEL2, on melanoma chromatin accessibility data can capture the various regulatory programs of the melanocytic and mesenchymal-like melanoma cell states. This model outperforms motif-based variant scoring, as well as more generic deep learning models. We detect hundreds to thousands of allele-specific chromatin accessibility variants (ASCAVs) in each melanoma genome, of which 15-20% can be explained by gains or losses of transcription factor binding sites. A considerable fraction of ASCAVs are caused by changes in AP-1 binding, as confirmed by matched ChIP-seq data to identify allele-specific binding of JUN and FOSL1. Finally, by augmenting the DeepMEL2 model with ChIP-seq data for GABPA, the TERT promoter mutation as well as additional ETS motif gains can be identified with high confidence. In conclusion, we present a new integrative genomics approach and a deep learning model to identify and interpret functional enhancer mutations with allelic imbalance of chromatin accessibility and gene expression.

Background: Hypoxia is pervasive in cancer and other diseases. Cells sense and adapt to hypoxia by activating hypoxia-inducible transcription factors (HIFs), but it is still an outstanding question why cell types differ in their transcriptional response to hypoxia. Results: We report that HIFs fail to bind CpG dinucleotides that are methylated in their consensus binding sequence, both in in vitro biochemical binding assays and in vivo studies of differentially methylated isogenic cell lines. Based on in silico structural modeling, we show that 5-methylcytosine indeed causes steric hindrance in the HIF binding pocket. A model wherein cell-type-specific methylation landscapes, as laid down by the differential expression and binding of other transcription factors under normoxia, control cell-type-specific hypoxia responses is observed. We also discover ectopic HIF binding sites in repeat regions which are normally methylated. Genetic and pharmacological DNA demethylation, but also cancer-associated DNA hypomethylation, expose these binding sites, inducing HIF-dependent expression of cryptic transcripts. In line with such cryptic transcripts being more prone to cause double-stranded RNA and viral mimicry, we observe low DNA methylation and high cryptic transcript expression in tumors with high immune checkpoint expression, but not in tumors with low immune checkpoint expression, where they would compromise tumor immunotolerance. In a low-immunogenic tumor model, DNA demethylation upregulates cryptic transcript expression in a HIF-dependent manner, causing immune activation and reducing tumor growth. Conclusions: Our data elucidate the mechanism underlying cell-type-specific responses to hypoxia and suggest DNA methylation and hypoxia to underlie tumor immunotolerance.

Understanding how enhancers drive cell type specificity and efficiently identifying them is essential for the development of innovative therapeutic strategies. In melanoma, the melanocytic (MEL) and the mesenchymal-like (MES) states present themselves with different responses to therapy, making the identification of specific enhancers highly relevant. Using massively parallel reporter assays (MPRA) in a panel of patient-derived melanoma lines (MM lines), we set to identify and decipher melanoma enhancers by first focusing on regions with state specific H3K27 acetylation close to differentially expressed genes. An in-depth evaluation of those regions was then pursued by investigating the activity of overlapping ATAC-seq peaks along with a full tiling of the acetylated regions with 190 bp sequences. Activity was observed in more than 60% of the selected regions and we were able to precisely locate the active enhancers within ATAC-seq peaks. Comparison of sequence content with activity, using the deep learning model DeepMEL2, revealed that AP-1 alone is responsible for the MES enhancer activity. In contrast, SOX10 and MITF both influence MEL enhancer function with SOX10 being required to achieve high levels of activity. Overall, our MPRAs shed light on the relationship between long and short sequences in terms of their sequence content, enhancer activity, and specificity across melanoma cell states.