Ludmila Recoules

Alexandre Heurteau

Flavien Raynal

Fatima Moutahir

Fabienne Bejjani

Isabelle Jariel-Encontre

Olivier Cuvier

Anne-Claire Lavigne
email: anne-claire.lavigne@ibcg.biotoul.fr

Kerstin Bystricky
email: kerstin.bystricky@ibcg.biotoul.fr

The histone variant macroH2A1.1 regulates gene expression by direct association with their transcription start site

The histone variant macroH2A1 (mH2A1) is involved in cellular growth, differentiation and reprogramming, but the underlying molecular mechanisms are a matter of debate. Different roles of mH2A1 in gene expression may relate to functional differences of its two splicing isoforms, mH2A1.1 and mH2A1.2. Here, we map for the first time genome-wide localization of endogenous mH2A1.1 and link the distribution of mH2A1.1 to control of gene expression in human breast cancer cells.

In addition to localization shared with mH2A1.2 to facultative heterochromatin, mH2A1.1 specifically associates with regulatory elements required for gene activation, super-enhancers and promoters of highly expressed genes. Depending on the recruitment profile of mH2A1.1 to these elements, selective depletion of mH2A1.1 up-or downregulates its target genes. mH2A1.1 represses transcription when its binding is spread over the entire gene and promoter, and activates transcription when its binding is strictly confined to the transcription start site (TSS). Notably, RNA Polymerase II was frequently in pause at mH2A1.1-activated genes. Functionally, mH2A1.1-dependent regulation of a subset of paused genes impedes mammary tumor cell migration. Molecular mechanisms of mH2A1.1 function at the TSS uncovered by our study define an intriguing new mode of transcription regulation in cancer cells.

Introduction

Compaction of DNA into chromatin modulates DNA accessibility, thereby regulating key molecular processes such as transcription [START_REF] Venkatesh | Histone exchange, chromatin structure and the regulation of transcription[END_REF][START_REF] Luger | New insights into nucleosome and chromatin structure: an ordered state or a disordered affair?[END_REF]. Histone posttranslational modifications, chromatin-remodeling enzymes, DNA-binding factors and architectural proteins fine-tune genome organization and dynamics [START_REF] Venkatesh | Histone exchange, chromatin structure and the regulation of transcription[END_REF][START_REF] Luger | New insights into nucleosome and chromatin structure: an ordered state or a disordered affair?[END_REF]. In addition, histone variants replace canonical histones in a locus-specific manner, which endows chromatin with additional properties required to regulate DNA accessibility and functions [START_REF] Buschbeck | Variants of core histones and their roles in cell fate decisions, development and cancer[END_REF].

Among the histone variants, macroH2A1 (mH2A1) is a vertebrate-specific [START_REF] Pehrson | Evolutionary conservation of histone macroH2A subtypes and domains[END_REF][START_REF] Rivera-Casas | The characterization of macroH2A beyond vertebrates supports an ancestral origin and conserved role for histone variants in chromatin[END_REF] histone H2A variant composed of an N-terminal H2A-like domain (64 % identical to H2A) and a C-terminal 25 kDa 'macro' domain. These two domains are joined by an unstructured 41 amino acid long linker that positions the macro domain outside of the nucleosome [START_REF] Gamble | Multiple facets of the unique histone variant macroH2A: From genomics to cell biology[END_REF]. Expression of the highly conserved H2AFY gene produces two splicing isoforms, mH2A1.1 and mH2A1.2, whose sequences differ in a 30 aminoacid region within the macro domain [START_REF] Gamble | Multiple facets of the unique histone variant macroH2A: From genomics to cell biology[END_REF]. mH2A1 was originally found to be enriched on the transcriptionally silent X chromosome [START_REF] Costanzi | Histone macroH2A1 is concentrated in the inactive X chromosome of female mammals[END_REF]. mH2A1 is also present at autosomes, forming large domains in association with histone marks associated with heterochromatin, such as H3K27me3 and H3K9me3 [START_REF] Gamble | The histone variant macroH2A1 marks repressed autosomal chromatin, but protects a subset of its target genes from silencing[END_REF][START_REF] Douet | MacroH2A histone variants maintain nuclear organization and heterochromatin architecture[END_REF][START_REF] Sun | Transcription-associated histone pruning demarcates macroH2A chromatin domains[END_REF]. In vitro studies have demonstrated that nucleosomal mH2A1 interferes with binding of the transcription factor NFkB, and inhibits nucleosome sliding by the remodeling complex SWI/SNF and initiation of RNA polymerase II (Pol II) transcription [START_REF] Doyen | Mechanism of polymerase II transcription repression by the histone variant macroH2A[END_REF][START_REF] Angelov | The histone variant macroH2A interferes with transcription factor binding and SWI/SNF nucleosome remodeling[END_REF]. Therefore, mH2A1 is believed to play a role in transcriptional repression. However, in a few cases the presence of mH2A1 has been seen to correlate with active transcription of a subset of genes involved in a variety of processes such as cell differentiation, lipid metabolism and cell-cell signaling [START_REF] Gamble | The histone variant macroH2A1 marks repressed autosomal chromatin, but protects a subset of its target genes from silencing[END_REF][START_REF] Dell'orso | The Histone Variant MacroH2A1.2 Is Necessary for the Activation of Muscle Enhancers and Recruitment of the Transcription Factor Pbx1[END_REF][START_REF] Chen | MacroH2A1.1 and PARP-1 cooperate to regulate transcription by promoting CBP-mediated H2B acetylation[END_REF][START_REF] Podrini | MacroH2A1 isoforms are associated with epigenetic markers for activation of lipogenic genes in fat-induced steatosis[END_REF][START_REF] Wan | MacroH2A1.1 cooperates with EZH2 to promote adipogenesis by regulating Wnt signaling[END_REF][START_REF] Changolkar | Genome-wide distribution of macroH2A1 histone variants in mouse liver chromatin[END_REF]. Thus, the roles of mH2A1 in regulating gene expression are seemingly contradictory. Studying mH2A1 isoforms separately may help to better understand their roles.

The two mH2A1 splice variants exhibit tissue-and cell-specific expression patterns [START_REF] Marjanović | MacroH2A1.1 regulates mitochondrial respiration by limiting nuclear NAD+ consumption[END_REF]. In normal cells, the mH2A1.2 isoform appears ubiquitously expressed [START_REF] Sporn | Histone macroH2A isoforms predict the risk of lung cancer recurrence[END_REF][START_REF] Cantariño | MacroH2A -An epigenetic regulator of cancer[END_REF][START_REF] Sporn | Differential regulation and predictive potential of MacroH2A1 isoforms in colon cancer[END_REF]. Mainly incorporated into heterochromatin (X-inactive chromosome and . CC-BY 4.0 International license perpetuity. It is made available under a preprint (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in The copyright holder for this this version posted January 29, 2020. ; https://doi.org/10.1101/2020.01. [START_REF] Xu | The histone variant macroH2A1.1 is recruited to DSBs through a mechanism involving PARP1[END_REF].924704 doi: bioRxiv preprint autosomic heterochromatin) [START_REF] Costanzi | Histone macroH2A1 is concentrated in the inactive X chromosome of female mammals[END_REF][START_REF] Gamble | The histone variant macroH2A1 marks repressed autosomal chromatin, but protects a subset of its target genes from silencing[END_REF][START_REF] Douet | MacroH2A histone variants maintain nuclear organization and heterochromatin architecture[END_REF][START_REF] Mermoud | Histone MacroH2A1 . 2 Relocates to the Inactive X Chromosome after Initiation and Propagation of X-Inactivation[END_REF], this splicing isoform is generally associated with gene repression. However, a few recent studies discovered that mH2A1.2 was required for gene expression. Indeed, in mouse muscle cells [START_REF] Dell'orso | The Histone Variant MacroH2A1.2 Is Necessary for the Activation of Muscle Enhancers and Recruitment of the Transcription Factor Pbx1[END_REF], mH2A1.2 binds muscle-specific enhancers necessary for the activation of the myogenic regulator network. In human embryonic stem cells (hESCs) [START_REF] Pliatska | The Histone Variant MacroH2A Blocks Cellular Reprogramming by Inhibiting Mesenchymal-to-Epithelial Transition[END_REF], binding of mH2A1.2 to gene promoters seems to control gene expression both positively and negatively.

However, the binding profile of mH2A1.2 in hESC was not conserved in mesodermderived human fibroblasts [START_REF] Pliatska | The Histone Variant MacroH2A Blocks Cellular Reprogramming by Inhibiting Mesenchymal-to-Epithelial Transition[END_REF].

In contrast to mH2A1.2, mH2A1.1 is only expressed in differentiated cells with low proliferation rates [START_REF] Sporn | Histone macroH2A isoforms predict the risk of lung cancer recurrence[END_REF][START_REF] Cantariño | MacroH2A -An epigenetic regulator of cancer[END_REF][START_REF] Sporn | Differential regulation and predictive potential of MacroH2A1 isoforms in colon cancer[END_REF]. In tumors, expression of the mH2A1.1 isoform is frequently reduced, as compared to normal tissues, suggesting that this isoform is a tumor suppressor [START_REF] Cantariño | MacroH2A -An epigenetic regulator of cancer[END_REF][START_REF] Sporn | Differential regulation and predictive potential of MacroH2A1 isoforms in colon cancer[END_REF][START_REF] Lavigne | Increased macroH2A1.1 expression correlates with poor survival of triplenegative breast cancer patients[END_REF]. In highly metastatic cancers such as triple-negative breast cancers however, expression levels of mH2A1.1 are increased and correlate with poor prognosis [START_REF] Lavigne | Increased macroH2A1.1 expression correlates with poor survival of triplenegative breast cancer patients[END_REF]. Only mH2A1.1 can bind NAD + metabolites through its macro domain [START_REF] Kustatscher | Splicing regulates NAD metabolite binding to histone macroH2A[END_REF] and interact with the DNA-damage repair and chromatin remodeling factor PARP1 (Poly-(ADP-Ribose) Polymerase 1) [START_REF] Chen | MacroH2A1.1 and PARP-1 cooperate to regulate transcription by promoting CBP-mediated H2B acetylation[END_REF][START_REF] Chaudhuri | The multifaceted roles of PARP1 in DNA repair and chromatin remodelling[END_REF][START_REF] Ouararhni | The histone variant mH2A1.1 interferes with transcription by downregulating PARP-1 enzymatic activity[END_REF][START_REF] Marjanović | MacroH2A1 . 1 regulates mitochondrial respiration by limiting nuclear NAD + consumption[END_REF]. Interaction between mH2A1.1 and PARP1 seems to be key for the capacity of mH2A1.1 to regulate DNA damage responses [START_REF] Xu | The histone variant macroH2A1.1 is recruited to DSBs through a mechanism involving PARP1[END_REF][START_REF] Kim | The histone variant macroH2A1 is a splicing-modulated caretaker of genome integrity and tumor growth[END_REF], mitochondrial respiration [START_REF] Marjanović | MacroH2A1 . 1 regulates mitochondrial respiration by limiting nuclear NAD + consumption[END_REF] and gene transcription [START_REF] Chen | MacroH2A1.1 and PARP-1 cooperate to regulate transcription by promoting CBP-mediated H2B acetylation[END_REF][START_REF] Ouararhni | The histone variant mH2A1.1 interferes with transcription by downregulating PARP-1 enzymatic activity[END_REF][START_REF] Chen | MacroH2A1 and ATM Play Opposing Roles in Paracrine Senescence and the Senescence-Associated Secretory Phenotype[END_REF]. In particular, it was proposed that mH2A1.1 and PARP1 recruit CBP (CREB-binding protein) to mediate acetylation of H2BK12/K120 and to regulate, in part, mH2A1-target gene expression in IMR90 cells, without really demonstrated the specific action of mH2A1.1 with respect to mH2A1.2 [START_REF] Chen | MacroH2A1.1 and PARP-1 cooperate to regulate transcription by promoting CBP-mediated H2B acetylation[END_REF]. In other hand, recruitment of mH2A1.1 to repressed signal-inducible genes is required for their stress response [START_REF] Chen | MacroH2A1 and ATM Play Opposing Roles in Paracrine Senescence and the Senescence-Associated Secretory Phenotype[END_REF][START_REF] Ouararhni | The histone variant mH2A1.1 interferes with transcription by downregulating PARP-1 enzymatic activity[END_REF].

Despite this knowledge about functions of mH2A1.1, neither its chromatin association, independently of the one of mH2A1.2, nor its specific mechanisms of action during gene activation, was clearly analyzed.

To gain a better understanding of the roles of mH2A1.1 in regulating gene expression, we determined 1) the genomic localization of mH2A1.1 and 2) its effect on gene expression. To that end, we generated a ChIP-grade mH2A1.1-specific antibody. Chromatin immunoprecipitation followed by sequencing (ChIP-seq) of mH2A1.1, of total mH2A1 and of a variety of histone marks, shows that both mH2A1 isoforms bind to heterochromatin and to enhancers. Surprisingly, only the mH2A1.1 isoform binds TSSs (Transcription Start Site) of active genes. Interestingly, we show that mH2A1.1-activated genes are mostly Pol II paused genes. Frequently, these genes are negatively implied in cell migration in mammary tumor cells. We demonstrate that mH2A1.1 impedes cell migration whereas mH2A1.2 promotes it.

Our work describes a novel gene activation pathway in cancer cells dependent on the selective recruitment of mH2A1.1 to the TSS of paused genes. Surprisingly, 80% of H3K27me3 overlapped with H3K9me3 (Fig 1C), although often presented as non-overlapping heterochromatin histone marks [START_REF] Allshire | Ten principles of heterochromatin formation and function[END_REF]. However, within peaks common to H3K27me3 and H3K9me3, the contribution of each mark was inversely proportional (Figs 1D and S5A). In addition, a high H3K27me3 to H3K9me3 ratio was more associated with genes whereas low H3K27me3 to H3K9me3 ratio was more associated with gene-poor genomic regions (Fig 1D andS5B). Hence, in this cell line, levels of relative enrichment between H3K27me3 and H3K9me3 could be used to distinguish "facultative-like" heterochromatin (with high levels of H3K27me3 and low levels of H3K9me3), from "constitutive-like" heterochromatin (with low levels of H3K27me3 and high levels of H3K9me3 To assess whether mH2A1.1 binding at the TSS correlates with transcription levels, we generated RNA-seq data from MDA-MB231 cells. We then stratified genes into four equal categories according to their expression levels (from silent to highly expressed genes). We determined the distribution of both mH2A1 isoforms at the TSS (+/-2 kb) and over the gene body (region from the TSS to the Transcription End Site (TES)) for those four groups of genes (see Materials and Methods). We found that the amount of mH2A1.1 bound to TSSs increased with expression levels while relative amounts recruited onto gene bodies decreased (Figs 3B, 3C and S7A-S7D). On the contrary, mH2A1 bound uniformly to TSSs and over the gene body of silent genes and its presence was restricted to the proximal promoter regions (around 1kb upstream ). Altered gene expression was confirmed by RT-qPCR on a subset of genes using two different siRNAs directed against the mH2A1.1 isoform (S3 Fig) . We found that all mH2A1.1 regulated genes, activated-as well as repressed-genes, were active or even highly expressed genes in WT conditions (S8B and S8C Fig) . Strikingly, we revealed that mH2A1.1 and mH2A1 binding differed between mH2A1.1-activated and repressed genes (Figs 4B-4D and S9A). Indeed, mH2A1.1 was mainly bound upstream of the TSS and within the gene body of mH2A1.1-repressed genes. At mH2A1.1-activated genes, however, mH2A1.1 exclusively associated with the TSSs and upstream (up to 1kb) (Figs 4B-4D and S9A). mH2A1 was bound upstream and downstream of the TSS of mH2A1.1-repressed genes, whereas its level was at its lowest at the TSS. In contrast, its relative enrichment was very low at mH2A1.1-activated genes, with its highest levels upstream TSS (up to 1 kb) (Figs 4B-4D and S9A). Overall, these results suggest that mH2A1.1 promotes transcription when its binding is restricted to TSSs and represses transcription when its binding is spread over genes. mH2A1.1 regulates expression of paused genes. We noted that Pol II was bound to the TSSs of mH2A1.1-activated genes but was barely detected on the corresponding gene bodies (Figs 4B-4D and S9A). This pattern of Pol II distribution raised the possibility that mH2A1.1-activated genes may be in pause. To confirm this, we calculated the Pol II pausing index (PI) for transcribed genes using Pol II ChIPseq data as described in [START_REF] Adelman | Promoter-proximal pausing of RNA polymerase II: emerging roles in metazoans[END_REF] (see Materials and Methods). To see if the binding of mH2A1.1 at TSSs was related to the level of Pol II pausing, we first plotted the mH2A1.1 ChIP-seq signal around the TSS +/-10 kb for each gene ranked by their PI (Fig 5A). We observed that confinement of mH2A1.1 to the TSS and its absence from the gene body is a characteristic of genes with a high PI. Genes where mH2A1.1 was spread along the entire gene body were characterized by low PIs (Fig 5A). In agreement, we found that mH2A1.1 peaks, as well as Pol II ones, were significantly sharper at the TSS (TSS +/-1kb) of genes with a high PI as opposed to genes with a low PI (Figs 5B, 5C and S9B). Moreover, 40% of genes bound by mH2A1.1 only at the TSS have a PI > 2 and 31% of paused genes (defined as genes with PI>2, n=6,821) have mH2A1.1 only peaks at their TSS (Fig 5D). Furthermore, the correlation between mH2A1.1 binding and activation of paused genes appeared to be functionally relevant because the majority of mH2A1.1-activated genes (68%) have a pausing index greater than 2 and their PI are significantly higher than that of any other gene category tested (Fig 5E and5F). In contrast, the majority of mH2A1.1-repressed genes (69%) was characterized by a PI < 2. These results suggest that mH2A1.1 spreading over the gene body represses genes whose expression is independent of Pol II pausing, whereas mH2A1.1 restricted to the TSS activates genes whose expression is dependent of Pol II pausing. Altogether, our results provide the first evidence of a functional link between the binding of mH2A1.1 at TSSs and Pol II pausing at promoter proximal regions. We propose that the association of mH2A1.1 with the TSSs of a subset of paused genes enhances their transcription.

Results

mH2A1.1 inhibits cell migration in MDA-MB231 cells. We found that genes deregulated by the loss of mH2A1.1 were involved in four main processes: cell cycle, DNA repair, cytoskeleton organization and cell adhesion (S10A, S10B Fig and S4 Table ). The two first processes were expected based on earlier studies [START_REF] Sporn | Differential regulation and predictive potential of MacroH2A1 isoforms in colon cancer[END_REF][START_REF] Xu | The histone variant macroH2A1.1 is recruited to DSBs through a mechanism involving PARP1[END_REF][START_REF] Kim | The histone variant macroH2A1 is a splicing-modulated caretaker of genome integrity and tumor growth[END_REF][START_REF] Novikov | QKImediated alternative splicing of the histone variant MacroH2A1 regulates cancer cell proliferation[END_REF].

However, the relationship between mH2A1.1 and expression of genes involved in cytoskeleton organization or cell adhesion was undocumented. Using a standard wide-field microscope, we observed that cells became more elongated after transfecting two different siRNA against mH2A1.1 (Fig 6A and6B). We used a siRNA against mH2A1.2 (S3A-S3C Fig) to test whether this effect was specific of mH2A1.1. Interestingly, we observed that mH2A1.2 KD cells became rounder compared to control cells (Fig 6A -6B). Importantly, we noticed that numerous mH2A1.1 KD deregulated genes modulate cell migration (examples of anti-migratory genes: ARRDC3 [START_REF] Draheim | ARRDC3 suppresses breast cancer progression by negatively regulating integrin beta4[END_REF], SOCS4 [START_REF] Mei | Interleukin-23 Facilitates Thyroid Cancer Cell Migration and Invasion by Inhibiting SOCS4 Expression via MicroRNA-25[END_REF], HACE1 [START_REF] Castillo-Lluva | The tumour suppressor HACE1 controls cell migration by regulating Rac1 degradation[END_REF] and FBXL4 [START_REF] Stankiewicz | Identification of FBXL4 as a Metastasis Associated Gene in Prostate Cancer[END_REF] -and of pro-migratory genes: EIF6 [START_REF] Pinzaglia | EIF6 over-expression increases the motility and invasiveness of cancer cells by modulating the expression of a critical subset of membrane-bound proteins[END_REF], MT1E [START_REF] Ryu | Role of metallothionein 1E in the migration and invasion of human glioma cell lines[END_REF], JUND [START_REF] Selvaraj | Extracellular signal-regulated kinase signaling regulates the opposing roles of JUN family transcription factors at ETS/AP-1 sites and in cell migration[END_REF] and DAPK3 [START_REF] Kake | Death-associated protein kinase 3 controls the tumor progression of A549 cells through ERK MAPK/c-Myc signaling[END_REF]). Thus, we investigated the effect of mH2A1.1 and mH2A1.2 depletion on the migratory capacity of MDA-MB231 cells using a Boyden Chamber assay (see Materials and Methods). Upon depletion of mH2A1.1, the migratory capacity of MDA-MB231 cells was significantly increased compared to control cells. In contrast, depleting the mH2A1.2 isoform led to a decrease in cell migration (Fig 6C and6D). Strikingly, mH2A1.1-activated genes involved in cytoskeleton organization and cell adhesion were also amongst genes with a high Pol II pausing index (Fig 6E ; see Materials and Methods). Taken together, these results suggest that mH2A1.1 inhibits cell migration by in part enhancing the expression of paused genes involved in cytoskeleton organization and cell adhesion.

Discussion

In this study, we present the first map of the histone variant mH2A1.1 genomic distribution over the entire human genome in breast cancer cells. Integration of this map with chromatin features confirms that mH2A1.1, together with its splice variant mH2A1.2, localizes to "facultative-like" heterochromatin domains and enhancers. Importantly, for the first time in cancer cells, we identify a direct link between mH2A1.1 recruitment of transcribed genes and their transcriptional regulation. The impact of this recruitment on the transcriptional rate is functionally dichotomous, with positive or negative effects. This bivalence correlates with a differential distribution of mH2A1.1 around the TSS as well as to the level of dependency of the transcribed gene to the Pol II pausing.

The results of our mH2A1 ChIP-seq analysis are overall consistent with previously published reports identifying this histone variant as a widespread chromatin feature on autosomes, covering large domains, with a preferential recruitment to H3K27me3-decorated "facultative-like" heterochromatin (Fig 1C) [START_REF] Gamble | The histone variant macroH2A1 marks repressed autosomal chromatin, but protects a subset of its target genes from silencing[END_REF][START_REF] Douet | MacroH2A histone variants maintain nuclear organization and heterochromatin architecture[END_REF][START_REF] Chen | MacroH2A1.1 and PARP-1 cooperate to regulate transcription by promoting CBP-mediated H2B acetylation[END_REF][START_REF] Lavigne | Composite macroH2A/NRF-1 Nucleosomes Suppress Noise and Generate Robustness in Gene Expression[END_REF]. Our mH2A1.1 ChIP-seq analysis allows us to show that even if mH2A1.1 distribution is mainly reminiscent of that of mH2A1.2, each variant has also specific heterochromatin localizations. Genome mapping of mH2A1. [START_REF] Luger | New insights into nucleosome and chromatin structure: an ordered state or a disordered affair?[END_REF] is not yet possible due to the poor specificity or lack of precipitation efficacy of commercial antibodies.

We further demonstrate that mH2A1.1-bound chromatin co-localizes significantly with the H3K9me3 histone mark (Fig 1C ). A fraction of these sites are devoid of H3K27me3 and could correspond to the recently identified macroH2A localization at constitutive heterochromatin 8 . However, the vast majority of mH2A1.1bound H3K9me3-decorated chromatin contained also tri-methylated H3K27 (Fig 1C). This difference may be a feature of the MDA-MB231 cell line, a high migratory capacity cancer cell line in which there is an abnormal expansion of frequently overlapping H3K9me3 and H3K27me3 marks along the genome [START_REF] Segal | The Heterochromatin Landscape in Migrating Cells and the Importance of H3K27me3 for Associated Transcriptome Alterations[END_REF][START_REF] Yokoyama | Cancer-associated upregulation of histone H3 lysine 9 trimethylation promotes cell motility in vitro and drives tumor formation in vivo[END_REF].

Moreover, we noted that even though mH2A1.1 associates with heterochromatin domains, its partial loss alone is not sufficient to reactivate silenced genes present in these domains (S8 Fig) . Similarly, even if mH2A1 binding was shown to overlap with H3K27me3 modified chromatin in primary human cells, no enrichment of H3K27me3 at mH2A1-regulated genes was observed [START_REF] Chen | MacroH2A1.1 and PARP-1 cooperate to regulate transcription by promoting CBP-mediated H2B acetylation[END_REF]. Furthermore, mH2A1.2-occupied and repressed target genes are not reactivated upon mH2A1.2 knock-down [START_REF] Dell'orso | The Histone Variant MacroH2A1.2 Is Necessary for the Activation of Muscle Enhancers and Recruitment of the Transcription Factor Pbx1[END_REF]. These observations suggest that enrichment of both mH2A1 isoforms in heterochromatin domains may serve as a lock to conserve heterochromatin stability and architecture [START_REF] Douet | MacroH2A histone variants maintain nuclear organization and heterochromatin architecture[END_REF].

In parallel, we show that the depletion of mH2A1.1 alone is sufficient to alter the expression of many genes present in chromatin domains devoid of H3K27me3 (Figs 4A andS8C). Previous work identified a role for mH2A1.1 in the transcriptional rate of genes encompassed in H2B acetylated regions but this mechanism was described as specific of primary cells and absent from cancer cells [START_REF] Chen | MacroH2A1.1 and PARP-1 cooperate to regulate transcription by promoting CBP-mediated H2B acetylation[END_REF]. Here, we demonstrate that the effect of mH2A1.1 depletion correlates with the distribution of mH2A1.1 along transcribed genes. Indeed, mH2A1.1 binding is restricted to the TSSs of mH2A1.1-activated genes while mH2A1.1 binds both upstream and on the gene body of mH2A1.1-repressed genes (Fig 4B). Our findings reinforce the notion that mH2A1.1 has a dual role in gene regulation and this seems to depend on its distribution on mH2A1.1-regulated genes [START_REF] Chen | MacroH2A1.1 and PARP-1 cooperate to regulate transcription by promoting CBP-mediated H2B acetylation[END_REF]. Molecular mechanisms determining mH2A1.1 differential localization and the consequent effects on gene regulation should now to be determined. Perhaps, mH2A1.1 post translational modifications favor binding to the TSS of active genes. Indeed, mH2A1-S137 phosphorylation excluded mH2A1 from the inactive X chromosome [START_REF] Bernstein | A phosphorylated subpopulation of the histone variant macroH2A1 is excluded from the inactive X chromosome and enriched during mitosis[END_REF]. Specific protein partners could also be involved in the recruitment of mH2A1.1 to the TSSs of transcribed genes, for instance PARP1. Indeed, association of mH2A1. Therefore, Ab αmH2A1 is certainly less efficient in immunoprecipitating mH2A1.1 in this context. Having a high affinity antibody led us to investigate new localizations and functions of the mH2A1.1 isoform. Now, it will be interesting to have a high affinity antibody against mH2A1.2. Indeed, through native-ChIP-seq experiments using a mH2A1.2 custom-made specific antibody, Dimitris Thanos and his collaborators recently showed that mH2A1.2 binds to the TSSs of transcribed genes in human embryonic stem cells [START_REF] Pliatska | The Histone Variant MacroH2A Blocks Cellular Reprogramming by Inhibiting Mesenchymal-to-Epithelial Transition[END_REF]. Thus, adapted sequencing techniques, isoformspecific antibodies and appropriate extractions of mH2A1 chromatin sub-complexes will enable studying mH2A1 splicing isoforms in different chromatin-and celldependend contexts.

Interestingly, we remarked that mH2A1.1 TSS-binding is distributed in a uniform manner around the TSS placing the maximum binding to the nucleosome free region (NRF) (Fig 3E). Its localization around the TSSs was also observed for mH2A1.2 using Native ChIP-seq [START_REF] Pliatska | The Histone Variant MacroH2A Blocks Cellular Reprogramming by Inhibiting Mesenchymal-to-Epithelial Transition[END_REF], eliminating the hypothesis of a technical bias.

The long-and unstructured-linker domain places the macro domain of mH2A1 outside nucleosomes [START_REF] Gamble | Multiple facets of the unique histone variant macroH2A: From genomics to cell biology[END_REF]. Because antibodies against mH2A1 frequently recognize the linker or the macro domain (as Ab αmH2A1.1 and Ab αmH2A1, respectively), we hypothesize that mH2A1.1 is incorporated at the two adjacent nucleosomes of the TSSs but their macro domain could be joined at the NFR, which could explain the higher signal at the NFR compared to the two adjacent nucleosomes. Future analysis will be necessary to confirm this hypothesis.

At the cellular level, we observed that the silencing of mH2A1.1 promotes cell migration, whereas mH2A1.2 silencing increases migration in the MDA-MB231 breast cancer cell line. mH2A1 has been shown to be involved in cell migration in mouse [START_REF] Dardenne | Splicing switch of an epigenetic regulator by RNA helicases promotes tumor-cell invasiveness[END_REF][START_REF] Kapoor | The histone variant macroH2A suppresses melanoma progression through regulation of CDK8[END_REF] and human [START_REF] Dardenne | Splicing switch of an epigenetic regulator by RNA helicases promotes tumor-cell invasiveness[END_REF][START_REF] Kapoor | The histone variant macroH2A suppresses melanoma progression through regulation of CDK8[END_REF][START_REF] Borghesan | DNA Hypomethylation and Histone Variant macroH2A1 Synergistically Attenuate Chemotherapy-Induced Senescence to Promote Hepatocellular Carcinoma Progression[END_REF][START_REF] Li | QKI5-mediated alternative splicing of the histone variant macroH2A1 regulates gastric carcinogenesis[END_REF] models. The opposite roles of both splice variants on the migratory capacities that we see are in agreement with observations in previous studies using gastric cancer cells [START_REF] Li | QKI5-mediated alternative splicing of the histone variant macroH2A1 regulates gastric carcinogenesis[END_REF] and MDA-MB231 cells [START_REF] Dardenne | Splicing switch of an epigenetic regulator by RNA helicases promotes tumor-cell invasiveness[END_REF]. At the molecular level, the selective incorporation of mH2A1.2 along the SOD3 gene (playing a part in cell migration) has been shown to be directly linked to gene repression in mouse metastatic cells [START_REF] Dardenne | Splicing switch of an epigenetic regulator by RNA helicases promotes tumor-cell invasiveness[END_REF]. Here, we show that the migratory capacities of human metastatic breast cancer cells are directly regulated by selective recruitment of mH2A1.1. Among mH2A1.1-activated genes were genes responsible for inhibiting cell migration. Thus, it is tempting to speculate that the negative effect of mH2A1.1 on cell migration is due to mH2A1.1-dependent activation of genes inhibiting migration. Interestingly, even if mH2A1.2 is very similar to mH2A1.1 in terms of amino-acid sequence, it has an opposite function to mH2A1.1 in the regulation of cell migration. Possibly the NAD+ metabolite binding pocket of mH2A1.1 plays a role in directing function. Here, we decipher in part the molecular mechanism by which mH2A1.1 could inhibit cell migration. Taken together, our data reveals antagonistic cellular functions of both mH2A1 isoforms and highlight the need to distinguish these isoforms when studying the role of mH2A1.

We demonstrate that mH2A1.1 favors expression of paused genes (70% of mH2A1.1-activated genes) (Fig 5E). It will be now interesting to dissect the molecular mechanisms that allow mH2A1.1 to promote transcription of paused genes. Perhaps, mH2A1.1 binding to the TSS could stimulate the recruitment of P-TEFb (Positive Transcription Elongation Factor) and consequently allow Pol II release [START_REF] Adelman | Promoter-proximal pausing of RNA polymerase II: emerging roles in metazoans[END_REF]. PARP1, well known partner of mH2A1.1, is also involved in Pol II pausing release by mediating ADP-ribosylation of NELF (Negative elongation factor) [START_REF] Gibson | Chemical genetic discovery of PARP targets reveals a role for PARP-1 in transcription elongation[END_REF] and could be recruited by mH2A1.1. This hypothesis could also explain why we observed mH2A1.1 at the TSS of highly transcribed genes (Fig 3B and3D). The presence of mH2A1.1 at SEs could also induce Pol II release (Fig 2C and2D). SEs are involved in cellular identity through the regulation of key genes involved in cellular identity [START_REF] Whyte | Master Transcription Factors and Mediator Establish Super-Enhancers at Key Cell Identity Genes[END_REF][START_REF] Hnisz | Super-Enhancers in the Control of Cell Identity and Disease[END_REF]. Furthermore, SEs are known to play an important part in many diseases, including several cancers 51-53 in which they drive expression of oncogenes [START_REF] Lovén | Selective inhibition of tumor oncogenes by disruption of super-enhancers[END_REF][START_REF] Donati | BRD4 and Cancer: going beyond transcriptional regulation[END_REF].

To note, several clinical trials already utilize SE blockers (bromodomain and extraterminal motif (BET) inhibitor and CDK7i) [START_REF] He | Targeting Super-Enhancers as a Therapeutic Strategy for Cancer Treatment[END_REF]. Among them, BRD4, one of the BET protein family members, was targeted. BRD4 binds acetylated histones at TSSs and SEs, brings them together, and mediates transcriptional activation and elongation by RNA Pol II [START_REF] Lovén | Selective inhibition of tumor oncogenes by disruption of super-enhancers[END_REF][START_REF] Donati | BRD4 and Cancer: going beyond transcriptional regulation[END_REF]. Interestingly, we observed that 74% of SEs are bound by mH2A1.1 and BRD4 in MDA-MB231 cells (data not shown). Could the presence of mH2A1.1 at SEs and TSSs be involved in the recruitment of BRD4 to enable elongation of RNA Pol II on paused mH2A1.1-activated genes ?

Overall, we demonstrate for the first time that mH2A1.1 activates transcription of paused genes by its selective recruitment to their TSSs. It remains to demonstrate whether mH2A1.1 is directly involved in the recruitment of factors acting on Pol II pausing events and/or on the establishment of 3D-chromatin structures (e.g chromatin looping mechanisms), facilitating Pol II pausing release.

Materials and Methods

Cell culture. MDA-MB231, HEK-293T and MCF7 cell lines were purchased from ATCC, and were maintained and amplified in Dulbecco's Modified Eagle's (DMEM) for HEK-93T and MDA-MB231 cells, and in DMEM-F12 for MCF7 cells, supplemented with gentamycin (50 µg/ml) (Gibco), fetal bovine serum (10%, Gibco) and sodium pyruvate (100 mM, Sigma). Cells were maintained in a humidified incubator at 37°C with 5% CO 2 . Cells lines were regularly tested for mycoplasma infection (MycoAlert, Lonza). In Montpellier, MDA-MB231 cells were cultured in DMEM supplemented with 10% fetal calf serum and penicillin/streptomycin (100 µg/ml each) and regularly tested for mycoplasma infection. Western blotting. Cells were lysed and subjected to western blot analysis as previously described [START_REF] Yang | Western blot: Technique, theory and trouble shooting[END_REF]. Briefly, proteins extracts were separated in 10% polyacrylamide (1:125 bisacrylamide:acrylamide) SDS gels, transferred onto nitrocellulose membrane (Bio-Rad) and blocked with PBS-Tween 0.4% -Milk 5% for 1h at RT with rotation. Membranes were then incubated with primary antibodies overnight (O/N) at 4°C in PBS-Tween 0.4% -Milk 5% with rotation (or 1h30 at RT).

Primary antibodies are described in the S1 Table Membranes were next incubated with secondary antibody in PBS-Tween 0.4% -Milk 5% 1h at RT with rotation and the signal was detected using chemiluminescence. Secondary antibodies are described in the S1 Table . Signal quantifications were carried out with Image Lab software (Bio-Rad).

RNA extraction, reverse transcription and quantitative real-PCR (qRT-PCR).

Total RNA was isolated using the RNAeasy midi kit (Qiagen). Purified RNA was reversed transcribed to cDNA using Maxima H Minus first Strand cDNA synthesis kit (Promega). The sequences of the primers used are available in S6 Table . RT-PCR was performed using iTAq Universal SYBR Green (Bio-Rad) according to manufacturer's instructions. At least two independent experiments were performed for each condition. The relative expression levels of mRNA were normalized to RPLP0 mRNA expression and evaluated according to the 2 -ΔΔCt method [START_REF] Rao | An improvement of the 2ˆ(-delta delta CT) method for quantitative real-time polymerase chain reaction data analysis[END_REF].

Fluorescence microscopy. Two or three days post-transfection, cells were fixed with 4 % paraformaldehyde for 15 min for MDA-MB231 cells and 10 min for HEK-293T at RT. Cells permeabilization was carried out using 0.1 % Triton X-100 in PBS for 10 min at RT. Cells were then blocked with 5 % BSA-0.15% Tween in PBS for 1h at RT. Next, cells were incubated with primary antibody O/N at 4°C. Cells were then incubated with Alexa conjugated secondary antibody for 1h at RT. Actin was labelled using cytoPainter Phalloidin iFluor diluted 1:1000 with secondary antibody according to the manufacturer's protocol (Abcam; see S1 Table ). Antibodies references and dilutions are provided in S1 Table . The coverslips were finally incubated with Hoechst (Invitrogen, 33342) for 30 min and then mounting with mounting media (Vectashield). Chromatin was eluted from the magnetic beads with DNA isolation buffer (2% SDS, 0.1 M NaHCO 3 ) for 1h at 65°C under agitation. Extracts were reverse-crosslinked with SDS O/N at 65°C. RNAs were degraded with RNase A and proteins were finally degraded with proteinase K. Same procedure was performed for input (10 μg of DNA). DNA was finally extracted with a phenol-chloroform extraction. Quantity and quality of DNA was tested with a nanodrop (NanoDrop2000, Thermo). Samples were sequenced with the GeT core facility, Toulouse, France (http://get.genotoul.fr).

Sequencing was done HiSeq3000-HWI-J00115 according to the manufacturer's protocol. Same procedure was done for ChIPqPCR, expected that 20 μg of DNA was used with 1μg of antibody). The sequences of the primers used are available in Supplementary Table 7. For western blot analysis, extracts (Input (10% IP), No immunoprecipitated (NoIP) fraction and IP fraction were processed as ChIP extract but not incubated with the proteinase K and RNAase A. Extracts were then subjected to western blot analysis as previously described in the western blot paragraph. To compare different extracts, we loaded 2 % of Input, 0.5 % of Input, 0.5 % of NoIP fraction and 20% of IP fraction. Percentages are relative to the DNA quantity used for ChIP (100 μg).

ChIP-seq of H3K27ac, H3K4me1, H3K4me3, H3K36me3 and Pol II were done essentially as previously described [START_REF] Tolza | AP-1 Signaling by Fra-1 Directly Regulates HMGA1 Oncogene Transcription in Triple-Negative Breast Cancers[END_REF]. Briefly, after cell fixation with 1 % of paraformaldehyde at RT for 5 min, cells were incubated in cell lysis buffer (PIPES 5 mM, KCL 85 mM, NP40 0.5%, Na Butyrate 10 mM, protease inhibitors) for 10 min on ice. After mild centrifugation, nuclei were lysed in Nuclei Lysis Buffer (Tris-HCL 50 mM pH 7.5, SDS 0.125%, EDTA 10 mM, Na Butyrate 10 mM, protease inhibitors) at 4°C for 2h and, then, sonicated for 10 cycles at 4°C using BioruptorPico device from Diagenode. For immunoprecipitation of H3K4me1, H3K4me3 and H3K27ac, 150 µl of chromatin (equivalent to 4.10 6 cells) and 4.5 µg of the corresponding antibodies were used. For Pol II, 850 µl of chromatin (equivalent to 22.10 6 cells) and 20 µg of the corresponding antibody were used. Two independent replicates for each ChIP were sequenced by the MGX genomic platform (Montpellier) using the Hi-seq2500 Illumina sequencer.

Strand-specific total RNA library preparation. Total RNA was isolated using the RNAeasy midi kit (Qiagen). RNA-seq quality and quantity control were performed using a Nanodrop (NanoDrop2000, Thermo) and BioAnalyser. Library preparation and sequencing was done by GeT core facility, Toulouse, France (http://get.genotoul.fr) with the kit TruSeq Stranded total RNA according to manufacturer's institutions. Sequencing was done HiSeq3000-HWI-J00115 according to the manufacturer's protocol.

ChIP-seq data processing. The quality of the reads was estimated with FastQC (Illumina, 1.0.0). Published ChIP-seq data of H3K9me3 (GSM2258862), H3K27me3 (GSM2258850) and corresponding input (GSM2258864) in MDA-MB231 cells were downloaded from GEODATASETS (https://www.ncbi.nlm.nih.gov/geo/, GEO accession number : GSE85158) [START_REF] Franco | Enhancer transcription reveals subtype-specific gene expression programs controlling breast cancer pathogenesis[END_REF], and reanalyzed as subsequently described.

Published ChIP-seq data of BRD4 in MDA-MB231 cells were downloaded (GSM2862187) and corresponding input (GSM2862178) from GEODATASETS (https://www.ncbi.nlm.nih.gov/geo/, GEO accession number : GSE107176) [START_REF] Chan | Polycomb complexes associate with enhancers and promote oncogenic transcriptional programs in cancer through multiple mechanisms[END_REF], and reanalyzed as subsequently described. Sequenced reads were aligned to the human genome Assembly GRCh38 using STAR (2.5.1) algorithm with defaults parameters [START_REF] Dobin | STAR: ultrafast universal RNA-seq aligner[END_REF]. Details are supplied in Supplementary Table 2. Low quality reads were then filtered out using Samtools (Samtools, options -q 10 -view) [START_REF] Li | The Sequence Alignment/Map format and SAMtools[END_REF].

Conversion of BAM files to bigWig files was performed with Bamcompare tool (DeepTools utilities v3.1.3) [START_REF] Ram | deepTools2 : a next generation web server for deep-sequencing data analysis[END_REF]. Corresponding ChIP-seq data generated from genomic DNA (Input) were used as control for every bigWig files normalization (options: --normalizeUsing RPKM --operation subtract --binSize 50 bp --smoothLength 150 bp). Peaks were determined with the enrichR function of NormR package [START_REF] Helmuth | normr: Normalization and difference calling in ChIP-seq data[END_REF] (v3.8, Options: fdr = 5e-2, binsize = 500 bp for H3K27ac, H3K4me1, H3K4me3, BRD4 and binsize =3000 bp for H3K36me3, PolII, H3K9me3, H3K27me3, Ab αmH2A1.1 and Ab αmH2A1). Peaks of Ab αmH2A1.1 (and Ab αmH2A1) spaced less than 3000 bp on the linear genome were merged. All downstream analyses were mainly performed with R studio. ChIP-seq signal and peaks positions visualization were obtained with IGV [START_REF] Robinson | Integrative genomics viewer[END_REF][START_REF] Thorvaldsdóttir | Integrative Genomics Viewer (IGV): high-performance genomics data visualization and exploration[END_REF]. Boxplots were done with ggplot2 [START_REF] Wickham | ggplot2: Elegant Graphics for Data Analysis[END_REF]. Distributions of mH2A1 isoforms and H3K27me3/H3K9me3 common peaks identified at specific genomic features were calculated using ChIPseeker package with default parameters (Figs 1B and S5B) [START_REF] Yu | ChIPseeker: an R/Bioconductor package for ChIP peak annotation, comparison and visualization[END_REF]. Statistical analyses are presented in Statistics and Reproducibility paragraph.

Identification of "putative" enhancers and super-enhancers. All putative enhancers correspond to the union of H3K4me1 and H3K27ac outside TSS (+/-2 kb) to avoid TSS bias (Fig 2B) [START_REF] Blinka | Identification of Transcribed Enhancers by Genome-Wide Chromatin Immunoprecipitation Sequencing[END_REF]. TSS annotation is based on TxDb.Hsapiens.UCSC.hg38.knownGene release. Super-enhancers were determined with ROSE utility tools based on H3K27ac signal (options : stitching_distance = 12.5 kb and TSS_exclusion_zone_size : 2500 bp) (Fig 2C and2D) [START_REF] Whyte | Master Transcription Factors and Mediator Establish Super-Enhancers at Key Cell Identity Genes[END_REF][START_REF] Lovén | Selective inhibition of tumor oncogenes by disruption of super-enhancers[END_REF].

Pol II pausing index calculation. Pausing index (PI) was defined as previously [START_REF] Zhang | Attenuation of RNA polymerase II pausing mitigates BRCA1-associated R-loop accumulation and tumorigenesis[END_REF], which is the ratio of Pol II (total) density in the promoter-proximal region ([-30;300] bp centered on TSS) to the total Pol II density in the transcribed regions (TSS + 300 bp to TES). Pausing index was only calculated for actively expressed genes (n=10,198, see RNA-seq analysis). Paused genes were defined as genes that have a PI upper to 2 (n=6,821).

Venn diagrams. Intersection of peaks were determined with the function findOverlaps() from GenomicRanges package. To note that for two ChIP-seq peaks intersections, only number of overlaps is contabilised and not the number of each peaks contained per overlap. This particularity explained why number of peaks changes between venn diagrams for a same ChIP-seq (i.e RNA-seq analysis. The quality of the reads was estimated with FastQC (Illumina, 1.0.0). The reads were mapped to the human reference genome GRCh38 using the default parameters of STAR (2.5.1) [START_REF] Dobin | STAR: ultrafast universal RNA-seq aligner[END_REF]. Details are supplied in S1 Table . Low quality reads and duplicates were then filtered out using SAMtools (Samtools, options -q 10 -view ; -rmdup) [START_REF] Li | The Sequence Alignment/Map format and SAMtools[END_REF]. Unstranded normalized Bedgraph files in Read Per Millions (RPM) were obtained with STAR using Output Wiggle parameters (options : --outWigTypebedgraph; --outWigStrandUnstranded; --outWigNorm RPM) [START_REF] Dobin | STAR: ultrafast universal RNA-seq aligner[END_REF]. Gene counts were performed with htseq-count utilities with default parameters (0.8.0) [START_REF] Stempor | SeqPlots -Interactive software for exploratory data analyses, pattern discovery and visualization in genomics[END_REF].

FPKM for all genes were calculated with the formula : FPKM = (RC g x 10 6 ) / (RC p x L) where RC g corresponds to the number of reads mapped to the gene, RC p to the number of reads mapped to all protein-coding genes and L, the Length of the gene in base pairs. Differential expression analysis was performed with DESeq2 package [START_REF] Love | Moderated estimation of fold change and dispersion for RNA-seq data with DESeq2[END_REF] with cutoff |FC| > 1.5 and padj < 0.1. Corresponding volcano plot was done with EnhancedVolcano package [START_REF] Enhancedvolcano | Publication-ready volcano plots with enhanced colouring and labeling[END_REF] (Fig 4A). The mH2A1.1 KD de-regulated genes are listed in S3 Table .    GO analysis. GO analysis was performed with LIMMA package (--function goana) (3.8) [START_REF] Ritchie | Limma powers differential expression analyses for RNA-sequencing and microarray studies[END_REF] and corresponding GO terms are supplied in the S4 Table . Selection of genes related to their functions was done with biomaRt package (function getBM()) [START_REF] Durinck | BioMart and Bioconductor: a powerful link between biological databases and microarray data analysis[END_REF][START_REF] Durinck | Mapping identifiers for the integration of genomic datasets with the R/Bioconductor package biomaRt[END_REF]. We took genes related to: cytoskeleton (GO:0005856), cell adhesion (GO:0007155), Cilium (GO:0005929) and cell junction (GO:0030054) (n=6,821).

Transwell migration assay. Transwell migration assays were performed using

Transwell plates with 0.8 µm pore polycarbonate membranes (CormingTranswell, Sigma). Three days post siRNA transfection, MDA-MB231 cells (10 5 ) were seeded in the upper chamber without FBS and allowed to invade to the reverse side of the chamber under chemoattractant condition with 10% FBS medium in the lower chamber. Following incubation for 16h at 37°C, the cells were fixed with 3.7% formaldehyde for 2 min at RT. Cells permeabilization was carried out with methanol incubation for 20 min at RT. Cells were then stained with Giesma for 15 min at RT. Same final total cell number between conditions was always checked by wide field microscope to avoid proliferation bias during number of migratory cell comparison. Not migrated cells were finally removed from the upper chamber by using a cotton swab. Migrated cells adhering to the underside of the chamber were photographed

  We determined the genomic localization of mH2A1.1 in the claudin-low breast cancer cell line MDA-MB231 which expresses this histone variant at a high level compared to other types of breast cancer cell lines (S1A, S1B Fig)[START_REF] Lavigne | Increased macroH2A1.1 expression correlates with poor survival of triplenegative breast cancer patients[END_REF]. To that end, we developed a ChIP-grade polyclonal rabbit antibody that exclusively recognizes mH2A1.1 (Ab αmH2A1.1) (S1C-S1F Fig) By ChIP-seq, we determined the distribution of mH2A1.1 and of total mH2A1 using Ab αmH2A1.1 and a commercially available ChIP-grade antibody (Ab37264 (Ab αmH2A1)) (S1, S2 Tables and S1G-S1I Fig), respectively.We identified 29,112 peaks for mH2A1.1 (Ab αmH2A1.1) and 22,789 peaks for total mH2A1 (Ab αmH2A1), covering combined 13.6 % of the genome. Genomewide, ChIP-seq results obtained with these two antibodies were highly similar with a significant co-occurrence of peaks (Fisher exact test (FET): p-value < 2.2x10 -16 and Odd ratio = 27.40) and a Pearson coefficient correlation (PCC) of 0.92 (S2A Fig). Analysis of peaks detected with Ab αmH2A1 and Ab αmH2A1.1 shows that regions occupied by mH2A1.1 and mH2A1.2 do not overlap completely (Fig 1A) suggesting that some regions are preferentially bound by either one of the variants. Indeed, a significant number of sites (18,838) is preferentially recognized by Ab αmH2A1.1 and not by Ab αmH2A1 (Fig 1A). Of note, Ab αmH2A1 has lower efficacy than Ab αmH2A1.1 in detecting mH2A1.1 and may, when reducing detection thresholds, recognize a fraction of these sites (S1B, S1G, S1H Fig). Genomic regions common to both antibodies correspond to one or more peaks exclusively recognized as mH2A1.1 or both isoforms in unknown relative proportions. For further analysis, we defined three types of genomic regions: regions incorporating preferentially only mH2A1.1 (mH2A1.1 only), regions incorporating preferentially only mH2A1.2 (mH2A1.2 only) and regions incorporating either mH2A1.1 or both isoforms (mH2A1s) (Fig 1A). These ChIP-seq results raise the possibility that preferential genomic localization of the mH2A1.1 variant relate to its function. To test this possibility, we first characterized the distribution of mH2A1.1 ChIP-seq peaks with respect to selected genomic features (Fig 1B; see Materials and Methods). Nearly one third of mH2A1.1 peaks were associated with distal intergenic regions and another third of mH2A1.1 peaks with promoters (TSS+/-1kb) (Fig 1B).Next, we analyzed mH2A1.1 genomic localization with respect to the chromatin environment of the identified regions. We first integrated ENCODE ChIPseq of heterochromatin histone marks (H3K9me3 and H3K27me3)[START_REF] Franco | Enhancer transcription reveals subtype-specific gene expression programs controlling breast cancer pathogenesis[END_REF]. We found that genome-wide, mH2A1.1 ChIP-seq peaks were significantly enriched with the H3K27me3 repressive histone mark, in a manner comparable to the one observed using Ab αmH2A1 (Figs1C and S2A). We confirmed the enrichment of mH2A1.1 on heterochromatin domains in independent ChIP-qPCR experiments on WT and mH2A1.1-deficient cells using two different siRNAs directed against mH2A1.1 isoform (S3A-S3E and S4 Figs). Remarkably, even if the related PCCs were low between H3K9me3 and mH2A1 isoforms (S2A Fig), in particular for mH2A1.1, we detected that both mH2A1 isoforms overlapped significantly with H3K9me3 (Fig 1C). At heterochromatin domains (i.e. marked by H3K27me3 and/or H3K9me3), both mH2A1 ChIP-seq peaks were highly similar with a PCC of 0.94 (S2B Fig).

mH2A1. 1

 1 binds to the transcription start site of active genes. A striking genomic feature identified for mH2A1.1 was its binding to TSSs (+/-1kb) (hereafter referred to as "promoter regions") (Fig 1B) confirmed by independent ChIP-qPCR experiments in WT and mH2A1.1-deficient cells by two different siRNAs directed against mH2A1.1 isoform (S3A-S3E and S4 Figs). Therefore, we characterized the association of mH2A1 isoforms at these genomic sites in more detail. Our ChIP-seq data showed that mH2A1.1 was more frequently associated with promoter regions than mH2A1 (41% vs 35%, respectively) (Fig 3A). Interestingly, among promoter regions bound by mH2A1 isoforms, 45% of them were occupied by mH2A1.1 only, whereas 14% were occupied by mH2A1.2 only (Fig 3A). We thus wondered whether the chromatin environment of promoter regions correlates with specific binding of either isoform. At promoter regions, we found that mH2A1.1 positively correlated with active marks (H3K4me3 and Pol II) whereas mH2A1 positively correlated with repressive histone marks (H3K27me3 and H3K9me3) (S6A Fig). Interestingly, this difference in localization identified by the two antibodies was maintained when we considered the two genomic populations specific to each antibody, mH2A1.1 only and mH2A1.2 only (S6B, S6C Fig). These results strongly suggest that mH2A1.2 preferentially binds to promoters in a closed chromatin state whereas mH2A1.1 is mainly recruited to promoters in an open chromatin state. In agreement, at promoter regions, both mH2A1 ChIP-seq experiments differed from one another with a PCC of 0.41 (in comparison to a PCC of 0.92 genome-wide) (S2A, S2B Fig).

  TSSs) and downstream of TESs of expressed genes (Figs 3B, 3C and S7A-S7D). In agreement, mH2A1.1 peaks were mainly present at the TSS +/-1kb of highly transcribed genes (Fig 3D), while mH2A1 peaks were present on the TSSs +/-1kb of silent genes (Fig 3D). Similar observations were seen for active and repressive histone marks (S7E Fig). Here again, this difference in localization identified by the two antibodies was even stronger when we considered the two genomic populations specific to each antibody, mH2A1.1 only and mH2A1.2 only (Fig 3D). Notably, we observed that mH2A1.1 binding at TSSs was distributed uniformly around TSSs, peaking at maximum Pol II binding (Fig 3E). In agreement, at TSSs, 70% of mH2A1.1 peaks (two third of them being mH2A1.1 only) overlapped significantly with Pol II (Fig 3F), with a PCC of 0.48 (in comparison to a PCC of 0.07 genome-wide) (S2A and S7F Figs). Coherent with these results, mH2A1.2 only peaks were hardly detectable (3%) at Pol II TSS-bound sites (Fig 3F), and at TSSs, both mH2A1 ChIPseq are anti-correlated with a PCC of -0.07 (S2B Fig). These results highlight highly specific association of the mH2A1.1 variant with the TSSs of transcribed genes. mH2A1.1 regulates gene expression. Binding of mH2A1.1 to TSSs of transcribed genes prompted us to test whether mH2A1.1 was required for modulating gene expression. We generated RNA-seq data from MDA-MB231 cells in which mH2A1.1 protein expression was reduced by around 90 % by siRNA without significantly affecting expression of the mH2A1.2 isoform (S3A-S3C Fig). 533 genes (56.3%) were down-regulated (mH2A1.1-activated genes) and 412 genes (43.7%) were upregulated (mH2A1.1-repressed genes) in mH2A1.1 knock-down (mH2A1.1 KD) conditions compared to WT conditions (Fig 4A, S8A and S3 Table

  1 and PARP1 was shown to be involved in gene transcription, DNA repair and global cellular metabolism [14,18,29]. The vast majority of mH2A1.1 narrow peaks at TSSs was only recognized by our home-made antibody Ab αmH2A1.1 (Fig 3). We think that this is certainly due to the high affinity of Ab αmH2A1.1 towards mH2A1.1 compared to Ab αmH2A1 (S1F-S1H Fig). The quantity of mH2A1.1-restricted to the TSSs is significantly smaller than in the large domains, such as heterochromatin and SEs (Figs 1D,2C and 3C).

  Transfection of siRNAs and plasmids. At 30-50% cell confluence, transfection of siRNA (11nM) was performed using INTERFERin (Polyplus-Ozyme) according to the manufacturer's protocol. Cells in control condition were transfected with INTERFERin without any siRNA Transfection of plasmid (1μg) was done with FuGene HD (Promega) according to the manufacturer's protocol. siRNA and plasmid sequences are available in S5 Table.Two and three days post plasmid and siRNA transfection respectively, cells were recovered for experiments.

  Images were acquired with Zeiss LSM 710 big confocal microscope using an x63 PL APO oil DIC On 1.4 objective for all experiments. Images were taken in Z-stacks with a voxel size of 300 nm. A Z-stack or Standard deviation intensity projections of Zstacks are shown. Chromatin immunoprecipitation and library preparation. Cells were cross-linked in DMEM containing 1.2 % of paraformaldehyde at RT for 10 min with rotation. Cross-link was stopped by the addition of glycine to a final concentration of 0.125M for 5 min. Cell were harvested and lysed in cell lysis buffer (10 mM Tris-HCl pH 7.4, 15 mM NaCl, 60 mM KCL, 1 mM EDTA, 0.1 mM EGTA, 0.2% NP-40, 5% sucrose). After 10 min in ice, cell lysis was amplified with a 2mL dounce (Kimble Chase) to enhance the nuclei separation from cytoplasm. Cell lysis buffer containing lysed cells was deposit up to a pillow buffer (10 mM Tris-HCl pH 7.4, 15 mM NaCl, 60 mM KCL, 1 mM EDTA, 0.1 mM EGTA, 0.2% NP-40, 10% sucrose). Nuclei were then pelleted by centrifugation and wash with washing buffer (10 mM Tris-HCl pH 7.4, 15 mM NaCl, 60 mM KCL). Nuclei were then resuspended in sonication buffer (50 mM Tris-HCl pH 7.5, 150 mM KCl, 5 mM EDTA, 1% NP-40, 0.1% SDS, 0.5 % Sodium deoxycholate, Protease Inhibitor (Roche)). Chromatin was sheared using a Bioruptor (Diagenode) (30 cycles, 30 sec ON/ 30 sec OFF) in order to obtain chromatin fragments with an average size of 300-500 bp. Quality and size of chromatin fragments was monitored by ethidium-bromide stained agarose gel electrophoresis after DNA purification. Then, 100 μg of DNA was incubated with antibody O/N at 4°C on a rotation wheel. Antibodies are described in the S1 Table.

1 %

 1 3 mg of protein A magnetic dynabeads (Sigma) were added for 3h at 4°C on a rotation wheel. Immunoprecipitates were then exposed to serial washes for 5 min each on a rotation wheel at 4°C in the following buffers (two times/buffer) : WB I : 2 mM EDTA, 20 mM Tris pH 8.1, 1 % Triton 100X, 150 mM NaCl, WB II : 2 mM EDTA, 20 mM Tris pH 8.1, Triton X100, 500 mM NaCl, WB III : 1 mM EDTA, 10 mM Tris pH 8.1, 250 mM LiCl, 1 % Sodium deoxycholate, 1 % NP-40 and WB IV : 1 M EDTA, 10 mM Tris pH 8.1.

Fig 1A Ab αmH2A1. 1 (

 1 n=29,112) and Fig 1C αmH2A1.1 (n=19,867)). The area-proportional Venn diagrams were drawn based on images generated by Vennerable package. For Venn diagrams in Fig 3, 5 and S6, intersections were performed at TSS (1-bp or +/-1kb centered on TSS). Enrichment tests associated to Venn diagrams are explained in Statistics and Reproducibility paragraph. Correlation heatmaps. Correlation heatmaps were done with multiBigwigSummary (options: -bins for whole genome (S2A Fig) or -BED-file for Heterochromatin, enhancer, TSS +/-1kp; TSS 1-bp, Figs S2B, S6A and S7B ) and plotCorrelation (option: -spearman correlation heatmap) from DeepTools utilities (3.1.3) [66]. Metagenes profils. Metagene analysis profiles were performed with R Seqplot package using bigWig files (function getPlotSetArray and plotAverage) [74]. Heatmaps profiles of Pol II, Ab αmH2A1.1 and Ab αmH2A1 around TSS +/-10 kb were ranked by mH2A1.1 KD de-regulated genes (mH2A1.1-activated and repressed genes, not ranked by the deregulation level) and done with Seqplot package (function getPlotSetArray and plotHeatmap) (Fig 4B). Heatmaps profiles of Pol II, Ab αmH2A1.1 and Ab αmH2A1 around TSS +/-10 kb were ranked by pausing index and carried out with Seqplot package (function getPlotSetArray and plotHeatmap)(Fig 5A).

  using a light microscope at x200 magnification (Invitrogen EVOS Digital Color Fluorescence Microscope). Cell counting was done with ImageJ in ten different fields per condition[80]. Three independent experiments were performed for each condition. Statistics and reproducibility. All western blot, RTqPCR and Boyden Chamber assay experiments were repeated at least twice as independent biological replicates and results are presented as mean +/-sd. All statistical analyses were done with R. For Western blot, RTqPCR and Boyden Chamber, Wilcoxon tests were used to compare mean values between conditions. p-values were considered as significant when * ≤ 0.05 and highly significant when ** ≤ 0.01; *** ≤ 0.001; **** ≤0.0001. All enrichment tests were performed on a base set made up of genomic bins. Each genomic bin is defined by merging ranges of ChIP-seq data (reduce function GenomicRanges R package) ChIP-seq list : Ab αmH2A1.1, Ab αmH2A1, Pol II, H3K4me3, H3K4me1, H3K36me3, H3K27ac, H3K27me3, H3K9me3, "putative" enhancers and "putative" Super-enhancers, on whole genome or TSS +/-1kb and TSS 1-bp. d'avenir" program managed by Agence Nationale pour la Recherche (contract ANR-10-INBS-09) and by the Fondation Recherche Medical (DEQ43940 to O.C team including A.H). We thank Marc Piechaszyk for critical reading of the manuscript. We acknowledge support from the light imaging Toulouse CBI platform. The work was generously funded by the Institut National du Cancer (INCA PL-BIO-16-269) to KB. Protoc [Internet]. 2009 Aug 23 [cited 2019 Jul 23];4(8):1184-91. Available from: http://www.nature.com/articles/nprot.2009.97 80. Caroline A Schneider WSR& KWE. NIH Image to ImageJ: 25 years of image analysis. INMATEH -Agric Eng [Internet]. 2012;9(7):671-5. Available from: http://dx.doi.org/10.1038/nmeth.2089

Fig 1 .

 1 Fig 1.Histone mH2A1 isoforms are incorporated to "facultative-like"

Fig 2 .

 2 Fig 2. The mH2A1.1 isoform is enriched at super-enhancers. (A) Genome

Fig 3 .

 3 Fig 3. mH2A1.1 is recruited to the transcription start site of active genes. (A)

Fig 4 .

 4 Fig 4. mH2A1.1 actives gene expression when it binds only on their TSSs. (A)

Fig 5 .

 5 Fig 5. mH2A1.1 binds a subset of TSSs of paused genes and promotes their

Fig 6 .

 6 Fig 6. mH2A1.1 inhibits cell migration by in part favouring expression of

S1 Fig.The antibody Ab αmH2A1.1 recognizes unequivocally only the isoform mH2A1. 1 .

 1 (A) RTqPCR on MDA-MB231 and MCF7 cells showing their relative mH2A1 isoforms mRNA expression level. Relative expressions are normalized by RPLP0 mRNA. (B) Western blot on whole cell extracts of MDA-MB231 and MCF7 cells showing the better affinity of Ab αmH2A1.1 to recognize mH2A1.1 compared to Ab#37264 (Ab αmH2A1). Ab#61427 is specific to mH2A1.2. GAPDH is used as a loading control. (C) Western blot showing specific recognition of mH2A1.1 isoform by Ab αmH2A1.1 antibody. For that, HEK-293T cells were transfected with plasmids coding for Flag-mH2A1.1 (Flag1.1) or Flag-mH2A1.2 (Flag1.2) fusion overexpressedprotein. Western blot was then done with Ab αmH2A1.1, Ab#Flag and Ab#E215 (that preferentially recognizes mH2A1.2) antibodies on whole cell extracts. GAPDH is used as a loading control. (D) Immunofluorescence in HEK-293T cells showing specific recognition of mH2A1.1 isoform by Ab αmH2A1.1. (E) Western blot on ChIP extracts from HEK-293T cells overexpressing Flag1.1 or Flag1.2 showing that Ab αmH2A1.1 immunoprecipites only mH2A1.1 isoform (Methods). Different extracts were loaded : Input fraction (Input), Non immunoprecipitated fraction (NoIP) and immunoprecipitated fraction (IP). Percentages represent percentages loaded on western blot compared to quantity used for ChIP. (F) Western blot showing that Ab αmH2A1.1 is also working in ChIP in MDA-MB231 cells on the endogenous protein. (G) As in (C), but for Ab αmH2A1 (Ab#37264) antibody showing that its antibody recognizes both isoforms but it less affine for Flag1.1 than Ab αmH2A1.1. (H) As in (G), but in immunofluorescence. (I), As in (E), but for Ab αmH2A1 (Ab#37264) antibody.

Acknowledgements

We thank M. Buschbeck from JCLR Institute at Barcelona for kindly providing Flag-mH2A1.1 and Flag-mH2A1.2 expression plasmids. ChIP-seq data against mH2A1.1 and mH2A1 isoforms as well as RNA-seq data were performed in collaboration with the GeT core facility, Toulouse, France (http://get.genotoul.fr), and were supported by France Génomique National infrastructure, funded as part of "Investissement

Data availability. ChIP-seq and RNA-seq data have been deposited to GEO under accession number GSE140022. Additional data are available upon reasonable request.

Author contributions

Competing interests

The authors declare no competing interests.

S2