Nucleosome repositioning during differentiation of a human myeloid leukemia cell line.

Cell differentiation is associated with changes in chromatin organization and gene expression. In this study, we examine chromatin structure following differentiation of the human myeloid leukemia cell line (HL-60/S4) into granulocytes with retinoic acid (RA) or into macrophage with phorbol ester (TPA). We performed ChIP-seq of histone H3 and its modifications, analyzing changes in nucleosome occupancy, nucleosome repeat length, eu-/heterochromatin redistribution and properties of epichromatin (surface chromatin adjacent to the nuclear envelope). Nucleosome positions changed genome-wide, exhibiting a specific class of alterations involving nucleosome loss in extended (∼1kb) regions, pronounced in enhancers and promoters. Genes that lost nucleosomes at their promoters showed a tendency to be upregulated. On the other hand, nucleosome gain did not show simple effects on transcript levels. The average genome-wide nucleosome repeat length (NRL) did not change significantly with differentiation. However, we detected an approximate 10 bp NRL decrease around the haematopoietic transcription factor (TF) PU.1 and the architectural protein CTCF, suggesting an effect on NRL proximal to TF binding sites. Nucleosome occupancy changed in regions associated with active promoters in differentiated cells, compared with untreated HL-60/S4 cells. Epichromatin regions revealed an increased GC content and high nucleosome density compared with surrounding chromatin. Epichromatin showed depletion of major histone modifications and revealed enrichment with PML body-associated genes. In general, chromatin changes during HL-60/S4 differentiation appeared to be more localized to regulatory regions, compared with genome-wide changes among diverse cell types studied elsewhere.

Introduction

Cell differentiation is known to be accompanied by chromatin changes and, in particular, nucleosome repositioning, since only about 37.5% of human nucleosomes have DNA sequence-determined positions. The magnitude of chromatin changes and their effect on gene expression depend upon the particular biologic system. Here, we consider 2 differentiation pathways of the human cell line HL-60/S4.

The immortalized myeloid leukemia cell line (HL-60), isolated from a female patient, was originally described to be an acute promyelocytic leukemia (APL); but was subsequently reclassified as an acute myeloblastic leukemia (AML) with maturation. HL-60 cells are multipotential and can be readily differentiated into granulocytes with DMSO or retinoic acid (RA), into monocytes with vitamin D3, or into macrophage with phorbol ester (TPA). The multipotential character of HL-60 cells to differentiate along various myeloid directions was summarized in an early review by one of the original discoverers of this important cell line. The mutagenically-derived subline HL-60/S4 available through ATCC as #CRL-3306, exhibits distinctly faster differentiation than the parent HL-60 cell line. HL-60/S4 cells develop granulocytic nuclear segmentation in 4 d, whereas the parent HL-60 line requires at least 6 d for the same level of differentiation. We have previously studied HL-60/S4, exploring nuclear shape, chromatin structure and cytoskeletal changes during differentiation induced by RA, TPA and vitamin D3. Furthermore, we have recently introduced the concept of “epichromatin,” which involves a small percentage (∼5%) of the HL-60/S4 genome, is highly enriched with retrotransposon Alu, and is reproducibly located adjacent to the nuclear envelope in undifferentiated, RA and TPA treated HL-60/S4 cells. While being adjacent to the nuclear envelope, this type of regions should be distinguished from lamina-associated domains (LADs) defined elsewhere, because epichromatin regions are identified after short formaldehyde fixation times (“snapshots”); whereas LADs are defined after longer transfection times, yielding 30–40% of the genome. Most recently, we have also described the transcriptomes of the undifferentiated, granulocytic and macrophage forms of HL-60/S4.

Differential gene expression, the basic mechanism underlying cell differentiation, depends upon a complex interplay between nucleosome positioning, histone and DNA covalent modifications and chromatin binding proteins. Integrating these various changes at individual promoters represents a current frontier in understanding the molecular biology of cell differentiation and cancer formation. In recent studies, epigenetic profiling has been performed in several blood cells types to study the differentiation-induced changes of the epigenome. For example, the BLUEPRINT consortium studied the differentiation of purified circulating monocytes from healthy volunteers into macrophages, and quantified this process in terms of the correlation of changes of H3K4me1, H3K4me3, H3K27ac, and DNase I accessibility with changes of gene expression. A more recent study has followed up this analysis, including 3 additional histone modifications and DNA methylation data, as well as performing NOME-seq to study the relation of nucleosome positioning with DNA methylation. However, the change of the chromatin density per se was not systematically analyzed in these works. Perhaps, one of the reasons for this is the intrinsic variability of nucleosome positions between healthy individuals – this problem can be overcome by using well-defined cell lines rather than primary cells. For example, in single-cell studies it is possible to distinguish quite nicely between HL-60 and GM12878 cells just by the characteristic differences of the chromatin accessibilities

In the present study, we have performed ChIP-seq with antibodies against histone H3 and 6 histone H3 modifications, combining these data with the RNA-seq to characterize the differentiation-specific chromatin changes in the HL-60/S4 cell system.

Results

Nucleosome repositioning at promoters and enhancers

We have determined nucleosome occupancy profiles using Micrococcal Nuclease (MNase)-assisted ChIP-seq with H3 antibody in untreated, RA-treated and TPA-treated HL-60/S4 cells. Complex scenarios of nucleosome gain or loss and nucleosome repositioning take place upon HL-60/S4 cells differentiation into granulocyte (RA-treated) or macrophage (TPA-treated) cell forms. Using a “sliding window” approach with the help of NucTools software, we defined 1000-bp genomic regions which are characterized by significant loss or gain of nucleosome occupancy upon differentiation with RA or TPA (see Methods). We have identified 22,024 regions with nucleosome loss and 224,140 regions with nucleosome gain upon RA treatment. Similarly, there were 28,732 regions with nucleosome loss and 237,850 regions with nucleosome gain upon TPA treatment. This is different from the previous NOME-seq study in primary monocyte cells differentiating into macrophages, which claimed preferential nucleosome depletions upon differentiation. These data are not directly comparable since we were calling fixed-size 1000-bp windows of differential histone H3 occupancy, while the NOME-seq study was calling single-nucleosome peaks, but the difference is interesting anyway.

Figure 1A and B present the statistics of differential gain/loss of nucleosomes in 3 types of HL-60/S4 genomic regions (promoters, enhancers and retrotransposon Alu elements) as a consequence of RA or TPA treatment. Promoters and enhancers were enriched within regions that lost nucleosomes, but not enriched within regions that gained nucleosomes. Alu elements showed no significant gain or loss of nucleosomes. Promoters of upregulated genes that lost nucleosomes upon RA treatment were enriched within a very large number of KEGG pathways including among others “B-cell receptor signaling pathway” (P = 9.5e-5) and “T cell receptor signaling” (P = 1.4e-3). Promoters of downregulated genes that lost nucleosomes upon RA treatment were enriched only with “DNA replication” (P = 2.7e-3) and “mismatch repair” pathways (P = 6.6e-3). Promoters of upregulated genes that lost nucleosomes upon TPA treatment were most significantly enriched for “cGMP-PKG signaling” (P = 4.4e-5) and “cAMP signaling” pathways (P = 2.8e-3). Promoters of downregulated genes that lost nucleosomes upon TPA treatment were enriched with “propanoate metabolism” (P = 1.3e-4) and “butanoate metabolism” pathways (P = 6.8e-4).

Figure 1.

Statistics of the distribution of regions that gained or lost nucleosomes upon HL-60/S4 differentiation by RA (A) or TPA (B), among different functional elements (enhancers, promoters, Alu repeats). Bar color code: light blue, nucleosomes gained, relative to 0 (undifferentiated) cells; yellow, nucleosomes lost, relative to 0. Fold enrichment equal to 1.0 signifies that a given genomic feature exhibits no net change in nucleosome numbers. Following RA or TPA treatment, more nucleosomes are lost, than are gained, in promoters and enhancers; but this is not observed with Alu repeat sequences. (C-D) Relation between changes of gene expression and nucleosome repositioning upon RA treatment (C) and TPA treatment (D). “Number of genes up/down,” the number of up or downregulated genes for 3 groups: “All genes,” all significantly expressed genes; “Nuc lost,” genes that lost nucleosomes at their promoters; “Nuc gained,” genes that gained nucleosomes at their promoters. The group of genes that lost nucleosomes at their promoters contains a larger proportion of upregulated genes, than the group of genes that gained nucleosomes at the promoter.

A typical parameter of chromatin structural changes occurring during cell differentiation is “nucleosome occupancy,” which can be interpreted as the fraction of cells from the population in which a given region of DNA is occupied by a histone octamer. Figure 1C and D show that when comparing the numbers of up/downregulated genes with lost nucleosomes, the genes which have lost nucleosomes were preferentially upregulated (i.e., increased transcript levels) compared with those which gained nucleosomes or to all genes. In particular, the ratio of the number of upregulated versus downregulated genes was around 1.14, both for “all genes” and for genes that gained nucleosomes at their promoters. On the other hand, this ratio (up/down gene regulation) increased to 1.39 and 1.5 for genes that lost nucleosomes at their promoters upon treatment by RA and TPA, respectively. Fisher's exact test confirmed that nucleosome loss at promoters upon RA treatment correlates with a significant increase in the number of upregulated genes (χ2-test, P = 1.7e-5); but the difference for genes that gained nucleosomes was not significant (χ2-test, P = 0.64). The same effect was observed for TPA-treated genes.

Figure 2 presents nucleosome occupancy changes for specific genes belonging to haematopoietic differentiation pathways that have lost nucleosomes in promoter regions, comparing undifferentiated (O) to differentiated (RA and TPA) cell forms. This set of genes generally exhibited large increases in transcript levels after both RA and TPA treatment (with the exception of CD38, which showed a decrease after TPA). At our H3 ChIP-seq resolution, we mostly observe distinct changes that can be localized to regions enclosing 1–5 nucleosomes. The loss of nucleosomes and corresponding increase of gene expression is particularly apparent in the CD38 (RA), IRF7 (RA and TPA), EGR2 (RA) and EGR3 (RA and TPA) gene examples. On the other hand, nucleosome gain does not seem to correlate with changes of gene expression, as exemplified in Supplementary Figure S1. Thus, nucleosome loss in promoter regions, but not nucleosome gain, appears to be a particular functional mechanism during HL-60/S4 differentiation by RA and TPA.

Figure 2.

Examples of nucleosome repositioning in genes that exhibit significant changes in mRNA levels in the differentiated state, compared with the undifferentiated HL-60/S4 cell form. A blue-green bar under the X-axis shows the gene body and the direction of transcription is indicated by an arrow. Promoters are defined as (+/−1000 bp) centered around the TSS (Transcription Start Site). Plot line colors: O (undifferentiated), black; RA treated, red; TPA treated, blue. Images available in color online.

Nucleosome repeat length during cell differentiation

Global nucleosome repeat length (NRL) was calculated for undifferentiated, RA and TPA-treated HL-60/S4 cells, following the calculation procedure described previously. Fig. 3A demonstrates that there are no significant differences in global NRL (192 +/−1.5 bp) among the 3 cell states. This apparent constancy of NRL within the HL-60/S4 cell states is interesting, because, for example, the genome-wide NRL values measured for mouse embryonic stem cells differentiation varied by 5–7 bp.

Figure 3.

Nucleosome repeat length (NRL) in different cell states and genomic regions. (A) Genome-wide average NRL is 192 bp without significant difference upon HL-60/S4 differentiation by RA or TPA. (B) Average nucleosome occupancy profile around bound CTCF. The corresponding NRL values are indicated on the figure. (C) Average nucleosome occupancy around bound PU.1 transcription factor. (D) The same as (C), but only the density of “+” strand reads is shown. The NRL determined for nucleosomes periodically arranged near PU.1 sites is indicated on this panel. The locations of CTCF and PU.1 binding sites were obtained from ENCODE ChIP-seq data on undifferentiated HL-60 cells; nucleosome positions were calculated from present data on undifferentiated and differentiated HL-60/S4 cells. Plot line colors: O (undifferentiated), black; RA treated, red; TPA treated, blue.

As pointed out in the earlier publications, specific genomic regions may display local variations and marked differences from the global NRL values. One such example is shown in Fig. 3B, which illustrates the global nucleosome occupancy surrounding CTCF binding sites in HL-60/S4 cells (undifferentiated, 180.8 +/−2.0 bp; RA, 181.8 +/−2.6 bp; TPA, 181.2 +/−2.4 bp), which is ∼10 bp smaller than the global average NRL. The center of the CTCF binding site exhibits a clear deficiency of bound nucleosomes, compared with the adjacent chromatin environment. Another specific example is the nucleosome occupancy around bound PU.1 transcription factors in the 3 HL-60/S4 cell states (Fig. 3C and D). The NRL values change by ∼2 bp during cell differentiation (undifferentiated, 177.4 +/−2.5 bp; RA, 179.4 +/−3.2 bp; TPA, 179.7 +/−2.0 bp). The chromatin region around bound PU.1 has ∼12 bp shorter NRL than the global NRL, and is also slightly smaller than that for CTCF. It remains to be demonstrated how general is the change in NRL around transcription factor (TF) binding sites, and whether such changes depend upon TF binding. In the examples above, we used CTCF and PU.1 binding sites determined by ChIP-seq in HL-60 cells by the ENCODE consortium.

Nucleosome occupancy changes associated with histone modifications

We have determined histone modifications using ChIP-seq in untreated, RA-treated and TPA-treated HL-60/S4 cells. As shown previously for other differentiation systems, nucleosome occupancy signatures of different histone modification domains vary significantly. The differentiation of HL-60/S4 is no exception. Figure 4 illustrates the nucleosome occupancy profiles of chromatin domains of 4 different histone modifications that are related to chromatin functions, comparing nucleosome occupancy in undifferentiated, RA- and TPA-differentiated cell states. Sites that were enriched with H3K4me3 (a marker for active promoters) in untreated cells reveal reduced (but still apparent) nucleosome occupancy in the center of their domain following treatment with RA or TPA (Fig. 4A). Sites that were enriched with H3K4me1, H3K9ac and H3K27ac in untreated cells reveal less dramatic changes in nucleosome occupancy following treatment (Fig. 4B, C and D). In contrast to activating marks, the sites enriched with H3K4me1 (which can be associated with repression) in untreated cells appear to show an increased nucleosome occupancy following either treatment. The finding that major decrease of the nucleosome density peak occurs at H3K4me3 sites in undifferentiated HL-60 cells is consistent with our statistics shown in Fig. 1A and B, which suggests significant nucleosome loss at some promoters upon HL-60/S4 differentiation.

Figure 4.

Nucleosome occupancies at different histone modification domains exhibit changes upon HL-60/S4 differentiation with RA or TPA. Panels: (A) H3K4me3; (B) H3K4me1; (C) H3K9ac; (D) H3K27ac. The largest changes occur in H3K4me3 domains, indicating major nucleosome loss at promoters upon HL-60/S4 differentiation (consistent with Fig. 1A and B). Coordinates of H3K4me1, H3K9ac and H3K27ac domains are obtained from the corresponding ChIP-seq experiments performed in this study. Coordinates of H3K4me3 domains are based upon ENCODE ChIP-seq data on undifferentiated HL-60 cells; nucleosome occupancy was calculated from present data on undifferentiated and differentiated HL-60/S4 cells. Plot line colors: O (undifferentiated), black; RA treated, red; TPA treated, blue.

Nucleosome features in HL-60/S4 epichromatin

Using the new data on nucleosome positioning and histone modifications during HL-60/S4 differentiation obtained in this study, we attempt to deepen our understanding of “epichromatin” regions, previously mapped in undifferentiated, RA and TPA treated HL-60/S4 cells. The concept of epichromatin (i.e., the surface of chromatin adjacent to the interphase nuclear envelope and at the “outer” surface of mitotic chromosomes) is based upon specifically localized immunostaining experiments. Two mouse monoclonal antibodies of divergent origin, PL2–638 and 1H6, yield identical staining patterns on fixed and permeabilized cells of considerable phylogenetic diversity, spanning from human to plant cells, implying the existence of a conserved chromatin epitope. Recently, 2 of us (ALO and DEO) enriched and characterized epichromatin from fixed and sonicated HL-60/S4 cells (0, RA and TPA cell states) by using a modified immunoprecipitation method (“xxChIP-seq”) using PL2–6 and 1H6. Some of the previously described characteristics of HL-60/S4 epichromatin regions include: 1) Representation of only ∼5% of the human genome. 2) Enrichment of retrotransposon Alu, ∼10-fold more concentrated than average for the human genome. 3) Epichromatin regions averaged ∼1 kb in size. 4) Epichromatin-containing DNA sequences, derived from 0, RA and TPA cells reacted with PL2–6 and 1H6 revealed a discontinuous distribution along each chromosome, which appeared to be relatively constant in chromosome location, comparing the 3 cell states. Still, the question of what exactly is epichromatin remains unanswered, and therefore we investigated the relationship between epichromatin regions and nucleosome occupancy changes during HL-60/S4 differentiation.

Since ChIP-seq is inherently a very noisy technique, we computed a list of common epichromatic regions, defined as the intersection of all replicate experiments between untreated and RA-treated cells. The intersection of 2 genomic regions is defined in the following way: 2 regions intersect if they share at least one base pair (based on their genomic coordinates). As input for these calculations, we took the coordinates of epichromatic regions for each replicate experiment reported earlier. Defined in this way, replicate PL2–6 immunoprecipitates of undifferentiated HL-60/S4 cells exhibited 40,288 epichromatin regions; replicate RA-treated cells, 77,897 regions. Defining only common regions based on all replicates, both antibodies (PL2–6 and 1H6) and 2 cell states (undifferentiated and RA-treated) resulted in 6541 common epichromatin regions.

We have found the following 3 properties of the common HL-60/S4 epichromatin regions: (1) Alu enrichment, compared with average genomic background levels, is about 2-fold at the center of epichromatin, falling to background at +/−2.5 kb from the center of the distribution (Fig. 5A). (2) The average nucleotide frequency reveals increased GC nucleotide content at the center of epichromatin (Fig. 5B). (3) Nucleosome occupancy (i.e., histone H3 read density) in common epichromatin regions exhibits a similar spatial distribution (comparing O, RA and TPA), showing highest occupancy in the center, declining to background levels at +/− ∼2 kb (Fig. 5C). This latter analysis is consistent with the notion that the epichromatin epitope is highly exposed within the epichromatin regions. In addition, ∼12% of the epichromatin data set from PL2–6 and 1H6 antibody experiments on untreated cells intersect with broad DNase I sensitive regions (identified in ENCODE for HL-60 cells). These results imply a significant enrichment of “open” chromatin within epichromatin regions. On the other hand, the calculated NRL, based on all candidate epichromatin regions, was not different from the genome-wide NRL, while the number of confirmed common epichromatin regions was too small to quantify NRL.

Figure 5.

Average feature distribution of common HL-60/S4 epichromatin regions as a function of the distance from the center of these regions. (A) Enrichment of retrotransposon Alu repeats. (B) Nucleotide frequencies. (C) Nucleosome occupancy profiles. Plot line colors for panel C: O (undifferentiated), black; RA treated, red; TPA treated, blue.

Signatures of histone post-translational modifications within common epichromatin regions are presented in Fig. 6. In all cases, modifications show a reduced presence in the middle of epichromatin: H3K9me3, a marker for constitutive heterochromatin (Fig. 6A); H3K4me3 (Fig. 6B) and H3K9ac (Fig. 6C), both markers for active promoters (although H3K9ac appears to be less depressed in TPA-treated cells); H3K27ac (Fig. 6D) and H3K4me1 (Fig. 6E), both markers for active enhancers; and H3K36me3 (Fig. 6F), a marker for actively transcribing genes. In summary, from the perspective of histone post-translational modifications, epichromatin is not clearly transcriptionally repressed, nor transcriptionally active.

Figure 6.

Signatures of histone modification marks around common epichromatin regions. Panels: (A) H3K9me3; (B) H3K4me3; (C) H3K9ac; (D) H3K27ac; (E) H3K4me1; (F) H3K36me3. Plot line colors: O (undifferentiated), black; RA treated, red; TPA treated, blue.

Finally, our analysis of the common epichromatin regions revealed a significant increase in gene promoters, compared with the average promoter density in the total human genome. Six percent of epichromatin regions contain promoters, corresponding to a 2.3-fold enrichment with respect to the genome average. Gene ontology (GO) analysis of these enriched promoters, unexpectedly, revealed enrichment for genes related to PML bodies (Fig. 7A and Supplementary Table 1). It is not clear what is the functional significance of having so many (i.e., 41) PML-associated genes in proximity to the interphase nuclear envelope (by virtue of their close continuity to epichromatin). The changes in transcript levels, compared with undifferentiated cells, for the specific PML-associated genes are frequently not the same for RA and TPA treated cells, implying that the PML bodies likely have different compositions in the 3 cell states. Figure 7B and C present nucleosome occupancy and transcript level changes for 2 important proteins of PML bodies (PML and SP100), also implying differential gene regulation between the 3 cell states. Because PML bodies appear to have multiple connections to stress responses the compositional heterogeneity of PML bodies suggests functional heterogeneity in the different cell states.

Figure 7.

(A) Gene Ontology (GO) analysis of common epichromatin promoter regions, indicating that these domains are enriched in genomic regions containing genes related to PML bodies. (B and C) Illustrating that PML bodies may have variable composition within the differentiated states: (B) Nucleosome repositioning in the PML gene, comparing 0, RA and TPA cell states; (C) Nucleosome repositioning in the SP100 gene, comparing 0, RA and TPA cell states. Plot line colors: O (undifferentiated), black; RA treated, red; TPA treated, blue.

Discussion

The present study used antibodies against histone H3 and several histone modifications in a ChIP-seq study of various nuclear elements (promoters, enhancers, TF binding sites and epichromatin) on undifferentiated (0) and differentiated (RA and TPA treated) forms of the myeloid leukemic cell line HL-60/S4. This yielded information about the relation of nucleosome loss to gene expression changes, the effects on nucleosome repeat length (NRL), nucleosome positions and nucleosome occupancy in and around CTCF and PU.1 binding sites and within common epichromatin regions. We have found changes in nucleosome occupancy surrounding various histone modification domains and revealed the distribution of histone modifications within and surrounding common epichromatin regions.

Our previous work on HL-60/S4 cells since 1998 was facilitated by their robust growth ability in the undifferentiated state and their rapid and reproducible cell differentiation in response to RA or TPA treatment. Furthermore, we recently obtained the transcriptomes of these 3 cell states (0, RA and TPA), and demonstrated that these cells have a very stable (abnormal) karyotype. Transcription factor binding measurements in HL-60/S4 are still missing, but we can extrapolate, to some extent, from the ENCODE ChIP-seq data for a related HL-60 cell line. All these data complemented the experimental data obtained in the current work, allowing us to perform the computational analyses.

Some of the most interesting conclusions from this study are the following: 1) During RA and TPA induced differentiation, regions that lost (but not regions that gained) nucleosomes are enriched at functional genomic elements such as promoters and enhancers (Fig. 1A and B). The number of genes upregulated upon nucleosome loss at the corresponding promoters was significantly elevated in comparison with genes that gained nucleosomes or in comparison with all genes (Fig. 1C and D). Documenting this observation, we presented a set of specific genes that show increased transcript levels and decreased promoter nucleosome occupancy in the differentiated state (Fig. 2). This observation has an intriguing parallel in the aforementioned transcriptome study, which demonstrated that following differentiation, more genes contributed more mRNA transcripts than genes that contributed decreased transcripts. 2) The global nucleosome repeat length (NRL) of HL-60/S4 cells appears essentially unchanged during cell differentiation, in contrast to other cell differentiation systems, such as stem cells differentiation into neuronal precursor cells and embryonic fibroblasts. 3) The NRLs around CTCF and PU.1 binding sites in undifferentiated and differentiated HL-60/S4 cells are about 10 bp shorter than the global NRL, arguing that transcription factor binding can lead to local changes of NRL. 4) Nucleosome occupancies can also change during HL-60/S4 differentiation in regions of specific histone modifications; this is, for example, well revealed around H3K4me3 domains. In the latter case, we observed a significant loss of nucleosomes following RA or TPA treatment (consistent with conclusion #1 above). 5) Epichromatin regions in undifferentiated and differentiated HL-60/S4 cells maintain their unusual characteristics, consistent with earlier observations. Alu and nucleosomes exhibit their peak densities at the “center” of the common epichromatin regions, which also shows a “spike” in G/C richness, compared with the surrounding regions. Common epichromatin regions also reveal a depletion of histone post-translational modifications considered indicative of either active or repressed transcription. In general, epichromatin is not within the gene body, but rather, is away from genes. If we assume that gene regulatory regions can extend up to 1000 kb from the promoter (unless it intersects another gene; the default parameter in the software GREAT used in this analysis), then it appears that the regulatory regions responsible for PML-associated genes (Supplementary Table 1) are enriched in epichromatin. At least 2 of these “close” promoters (i.e., for SP100 and DAXX) regulate transcription for proteins that are key components of ALT-associated PML bodies, believed to play a role in telomere elongation independent of telomerase. The large number of common chromosomal locations of epichromatin, comparing the 3 cell states, suggests stability of nuclear architecture during the lineage differentiation of HL-60/S4 cells. It is important to emphasize that epichromatin is identified by a conserved epitope (involving histones H2A, H2B and DNA) that remains to be defined, but likely represents a highly conserved structural feature of nucleosomes. A recent study presents evidence that etoposide-induced DNA damage of an embryonic teratocarcinoma cell line results in a disappearance of the peripheral epichromatin epitope in damaged nuclei, but not in undamaged nuclei. The epitope disappearance may represent histone-DNA loss, histone modification or masking by other proteins, an intriguing probe for changes at the nuclear periphery.

Finally, the present manuscript argues that in HL-60/S4 cells, there is a statistically significant association between the loss of nucleosomes within gene promoters (in differentiated cells) and increased levels of gene transcripts. Nucleosome depletion at promoters in many systems is quantitatively linked to the level of gene expression. However, it is also clear that there are promoters that require the presence of nucleosomes to initiate the binding of transcription factors, leading to subsequent transcription. Our observation that, during HL-60/S4 differentiation, nucleosome loss and nucleosome gain affect gene expression asymmetrically (i.e., a clear association with nucleosome loss and a more complicated response to nucleosome gain) agrees with the concept that there are several different mechanisms by which nucleosome positioning can affect gene expression. Interestingly, these changes in nucleosome occupancy involve extended regions of one to several kb, which is different from previously considered single-nucleosome repositioning scenarios and parallels a similar effect recently found at the V(D)J recombination domains The molecular mechanisms of the change of chromatin density in extended regions that we report here may be linked to chromatin compartmentalization in domains such as TADs (topologically associated domains), which have been recently shown to be correlated with regions of differential chromatin accessibility. An exploration of this idea may provide a direction for new studies.

Materials and methods

Cell culture and cell preparation

HL-60/S4 cells were maintained in RPMI 1640 medium, plus 10% fetal calf serum and 1% Pen/Strep/Glutamine. This cell line is available from ATCC (www.atcc.org, #CRL-3306). For the current experiments, cells were maintained in duplicate T-75 flasks, each flask containing 30 ml. Undifferentiated HL-60/S4 cells were seeded at a concentration of 0.83×105. RA treated cells were seeded at a concentration of 1.67×105, with RA (Sigma-Aldrich R2625) added to 1 µM. TPA treated cells were seeded at a concentration of 4 × 105, with TPA (Sigma-Aldrich P1585) added to 16 nM. After 4 d, approximately 107 were removed from each flask.

Harvested cells were centrifuged at 300xg for 5 min, resuspended in PBS and centrifuged again. The cell pellets were suspended in 20 ml PBS, made 1% HCHO (fresh ampule of 16% HCHO, methanol-free) and incubated 10 min at RT on a rotator. Fixation was stopped by addition of 2.5 M glycine to a final concentration of 125 mM and incubated 5 min at RT on a rotator. Fixed cells were centrifuged at 1000xg for 5 min and washed twice in PBS plus 50mM PMSF. Cell pellets were suspended in 10 ml of ice-cold “swelling buffer” (10 mM KCl, 1 mM MgCl2, 0.1% NP-40, 1 mM DTT, 25 mM Hepes pH 7.8, made 0.5 mM PMSF and containing Roche COMPLETE), incubated 10 min on ice, then centrifuged at 1000xg followed by plunging the cell pellets in liquid N2.

ChIP-seq protocol

The frozen cell pellets were resuspended in MNase (Micrococcal Nuclease) buffer (25 mM KCl, 4 mM MgCl2, 1 mM CaCl2, 50 mM Tris/HCl pH 7.4) and 10 U MNase per 1×106 cells were added. After 15 min, incubation at 37°C MNase was stopped by adding 10x Covaris buffer (100 mM Tris pH 8.0, 2 M NaCl, 10 mM EDTA, 5% N-lauroylsarcosine, 1% Na-deoxycholate, supplemented with protease inhibitors). The samples were sonicated for 15 min with the following parameters with a Covaris S2 system: burst 200, cycle 20%, and intensity 8. Following centrifugation the supernatant was collected and directly used for IP.

After IgG preclearance the sheared chromatin was incubated overnight with protein G magnetic beads (Cell signaling, 9006) and 4 μg of anti-pan H3 (Abcam, ab1791), H3K4me1 (Abcam, ab8895), H3K4me3 (Abcam, ab8580), H3K9ac (Active Motif, 39917), H3K9me3 (Abcam, ab8898), H3K27ac (Abcam, ab4729) or H3K27me3 (Abcam, ab6002), H3K36me3 (Abcam, ab9050). After washes with 1x Covaris buffer (10 mM Tris–HCl, pH 8.0, 200 mM NaCl, 1 mM EDTA, 0.5% N-lauroylsarcosine, 0.1% Na–deoxycholate), high-salt-buffer (50 mM HEPES pH 7.9, 500 mM NaCl, 1mM EDTA, 1% Triton X-100, 0.1% Na–deoxycholate, 0.1% SDS), lithium buffer (20 mM Tris–HCl pH 8.0, 1 mM EDTA, 250 mM LiCl, 0.5% NP-40, 0.5% Na–deoxycholate) and 10 mM Tris–HCl, chromatin was eluted from the magnetic beads (elution buffer: 50 mM Tris pH 8.0, 1 mM EDTA, 1% SDS, 50 mM NaHCO3) and the crosslink was reversed overnight. After RNase A and proteinase K digestion, DNA was purified and cloned in a barcoded sequencing library for the Illumina sequencing platform. In brief, after DNA repair and A-addition NEBNext adapters (NEB, E7335) were ligated and digested with the USER enzyme. Barcodes (NEB, E7335) were introduced via PCR with a maximum of 14 cycles by the NEBNext polymerase (NEB, M0541). Size selection for mononucleosomal insert fragments was done with Ampure XP beads (Agencourt, A63880). Each ChIP-seq library was sequenced on the Illumina HiSeq 2000 with 50bp single-end in 2 replicates, whereas the H3-chip samples were sequenced with 50 bp paired-ends in 2 replicates.

Analysis of H3 ChIP-seq

H3 ChIP-seq reads were mapped using Bowtie allowing up to 1 mismatch and counting only unique hits, obtaining the following number of mapped 100-bp paired-end reads: 248 million for undifferentiated HL-60/S4 cells, 319 million for RA-treated HL-60/S4, and 338 million reads for TPA-treated HL-60/S4. The average DNA fragments lengths were correspondingly 160 bp, 159 bp and 159 bp, reflecting a moderate chromatin digestion (in a strongly digested chromatin DNA fragments are closer to the 147bp nucleosome DNA length). Mapped reads were processed using NucTools to generate genome-wide nucleosome occupancy landscapes, extract individual genomic regions and calculate NRL, as described previously. Aggregate nucleosome occupancy profiles around genomic features were calculated using HOMER and visualized using OriginPro 2016 (OriginLab Corporation). For the analysis performed in this manuscripts, we used the corresponding RefSeq definitions of promoters (extended +/−1,000 from TSS), FANTOM definitions of enhancers and Dfam definitions of ALU repeats. The coordinates of transcription factor binding sites and DNaseI-sensitive regions in HL-60 cells determined by the ENCODE consortium were obtained from the following GEO entries: GSM749688 (CTCF), GSM1010843 (PU.1) and GSM736626 (DNaseI hypersensitivity).

Analysis of ChIP-seq of histone modifications

The basic data processing was performed similarly to H3 ChIp-seq, with a difference that these were single-end reads, and the sequencing was to a smaller depth (30–40 million reads for each of 6 histone modification for each of the 3 cell conditions). Peak calling was performed using HOMER and repeated with MACS and SICER to confirm that the number of peaks did not differ significantly for different peak calling procedures. In the case of H3K4me3, the data from the ENCODE consortium obtained for HL-60 cells (GSM945222) had higher genomic coverage than our H3K4me3 ChIP-seq in HL-60/S4, but the aggregate profiles around genomic features were qualitatively similar.

Finding genomic regions which lost/gained nucleosomes

1000-bp regions which lost or gained nucleosome were determined genome-wide using NucTools with the following criteria applied to the pairs of cell states (e.g. RA vs 0): If the average nucleosome occupancy in RA-treated cells (PRA) in a given genomic window was higher than that in the control state (P0), i.e., if 2x (PRA – P0)/(PRA + P0) > 0.99 this was assigned as a nucleosome gain region. If 2x (PRA – P0)/(PRA + P0) < −0.99 the corresponding region was considered to have a nucleosome loss.

Analysis of RNA-seq data

We used RNA-seq data for undifferentiated HL-60/S4 and RA- and TPA-differentiated counterparts (4 replicates for each state) reported in the accompanying paper (Mark Welch et al, Nucleus, in press; deposited at the Short Read Archive (SRA), bioproject PRJNA303179). Because different methods can produce inconsistent results we measured gene-level differential expression using DESeq2, EBSeq, EdgeR, and limma-voom. For our data, EBSeq proved the most conservative, with the fewest genes called as significantly differentially expressed and the fewest genes called, that were not called by other methods (Supplemental Table 2). Therefore, only genes that EBSeq determined to have posterior probability of differential expression (PPDE) >0.95 were considered in the analysis (Supplemental Table 3).

Epichromatin coordinates

Coordinates of epichromatin domains in undifferentiated HL-60/S4 and RA-differentiated cells were determined in our previous publications using modified immunoprecipitation method (“xxChIP-seq”) with PL2–6 and 1H6 antibodies. Here we have defined a subset of these regions that was reproduced in all replicate experiments with both antibodies, both in undifferentiated HL-60/S4 and RA-differentiated cells. Defined in this way, replicate PL2–6 immunoprecipitates of undifferentiated HL-60/S4 cells exhibited 40,288 epichromatin regions; replicate RA-treated cells, 77,897 regions. Defining only common regions based on all replicates, both antibodies (PL2–6 and 1H6) and both cell states (undifferentiated and RA-treated) resulted in 6,541 epichromatin regions.

Enrichment analysis

Enrichment of genomic regions which lost or gained nucleosomes at functional genomic elements (promoters, enhancers, Alu repeats) was calculated as a ratio of several regions falling to a given class of functional elements vs. the number of random regions (generated for a data set of the same and with the same region length distribution) falling to a given class of functional elements. Random regions were generated using BedTools. Gene Ontology enrichment at promoters with differential nucleosome occupancy was calculated using DAVID (2016 update) with default parameters. The enrichment of epichromatin regions in functional regulatory elements was calculated using GREAT with the following parameters: Each gene was assigned a basal regulatory domain of 5kb upstream and 1kb downstream of the TSS (regardless of other nearby genes). The gene regulatory domain is extended in both directions to the nearest gene's basal domain but no more than 1000kb in one direction.

Permission

No human patients or experimental animals were used. Approval by the Ethics Committee of Human Experimentation was not required.