Note: This rebuttal was posted by the corresponding author to Review Commons. Content has not been altered except for formatting.
Learn more at Review Commons
Reply to the reviewers
We thank Review Commons and its three reviewers for your supportive and insightful responses to our manuscript. Below, we provide detailed responses to the reviewers’ individual comments and how we plan to address them during the revision.
Reviewer #1: **Major comments:**
The manuscript is very well written. The data is clearly presented. The methods are explained in sufficient detail with a few exceptions mentioned below, and statistical analysis are adequate. There are some concerns and suggestions about the experimental design and data presentation.
- Drug treatments. It is not clear whether the cells were previously grown on charcoal-stripped serum before hormone treatments. From methods, it seems they were grown in 5% FBS and directly treated with the hormones. Also, what "hormone-free medium" mean? Is it charcoal stripped Serum or not Serum at all?
For all experiments, the cells were grown in medium containing 5% FBS. Throughout the manuscript, “hormone-free” refers to medium containing 5% FBS with no dexamethasone added. Technically, this medium is not hormone-free as FBS contains low levels of cortisol. However, the levels of cortisol from the FBS in our medium seems insufficient to elicit a transcriptional response or DNA binding by GR based on experiments comparing charcoal stripped and medium containing regular 5% FBS. However, we acknowledge that it should be made clear to the reader that growth conditions technically were not hormone-free. We will make sure to include this information in both the methods and results section of a revised manuscript. In addition, we will state explicitly that our naiive cells are those that have not been exposed to a high dose of hormone.
Replicates for these data sets? The ATAC and Chip-Seq should have at least 2. The concordance of the ATAC-seq and Chip-seq replicates should be described and shown in supplemental figures.
The ChIP-seq peaks for GR are the intersect of two biological replicates. This is described in the Methods section (page 7). For the ATAC data, we used two biological replicates for the vehicle treated cells and treated two different hormones (dexamethasone and cortisol) as replicates. In a revised manuscript, we will add a supplemental figure to show the concordance between the replicates.
Fig1A - The ATAC-seq HM should be clustered to show which peaks in opening/closing and unchanged peaks also have called GR chip peaks. Showing browser shots as in Fig1B is cherry picking data and can be put in a supplementary figure as an example. This is a main point of emphasis of the manuscript so show the data. The atac peaks that do overlap with GR chip peaks should be sorted by GR peak intensity. The QPCR is then only needed to confirm the quantitative changes.
This is a good idea. As suggested by this reviewer (and also in response to a comment by one of the other reviewers), we will revise this figure panel to make the overlap between GR binding and opening and closing sites more obvious. Here are the numbers:
A549 cells:
opening sites: 49%
closing: 10%
nonchanging: 18%
U2OS cells:
opening: 54%
closing: 0.2%
nonchanging: 7%
Regarding the use of browser shots, obviously these are cherry picked examples, however in our opinion they serve a purpose beyond illustrating examples of individual loci that open or close as they also give the reader an idea of the quality of the ATAC-seq data.
To show both the ATAC sites and H3K27ac sites are specific to hormone treatment, a random set of 15K peaks not in this peak set also should be shown in HMs and should not change with the treatments. Why does the H3K27ac go down in the 6768 non changing sites with dex?
The proposed group of control peaks is essentially what we included as “non-changing” peaks. For the revision, we will refine this group and compare the H3K27ac signal between GR-occupied and non-occupied groups. Regarding reduced H3K27ac signal upon Dex treatment at non-changing sites: Notably, this comparison is based on a single ChIP-seq replicate. In our experience, ChIP-seq experiments show quite some variability between biological replicates, which limits our ability to compare signal levels quantitatively. Thus, the difference could simply reflect a difference in ChIP efficiency between the treated and untreated cells. Alternatively, it could be that there is a general redistribution of H3K27ac signal towards GR-occupied opening sites. To pin down which of these explanations is valid, we would need to perform additional experiments, e.g. using spike-ins. However, this is beyond what we can do at the moment and therefore, we will instead revise the text to make sure that the interpretation of these results is somewhat speculative.
The D & E parts of Fig1 can then be eliminated to become parts of Fig1A. Its not clear in the text that the HMs in Fig1 are all sorted in the same way.
We will revise figure 1 as requested. In our initial submission, the data was always sorted by signal intensity of the feature shown. We will revise this and sort by ATAC-signal and keep a consistent sorting order for other features shown (and stratify each group into GR-occupied or not).
- Fig. 1b (and d). The ChIP data is from 3h-hormone treatment while the ATAC-seq data is from a 20h hormone treatment. It seems a bit misleading to directly compare GR occupancy with the state of the chromatin at different time windows. Shouldn't the authors show their ATAC-seq 4h treatment data (shown in Fig S1) here instead?
We will revise the figures as suggested to show the same time point for ChIP and ATAC-seq data.
- Fig. 1f. The authors state "downregulated genes only show a modest enrichment of GR peaks". However, there is a significant enrichment of GR-peaks in repressive genes compared to non-regulated genes. It would be interesting to see how some of these peaks look in a browser shot. While the general conclusion "transcriptional repression, in general, does not require nearby GR binding", seems valid, the observation that many GR peaks appear directly bound to nearby repressed genes ought to be more emphatically recognized in the text.
This is a fair point and was also raised by the other reviewers. During the revision, we will make textual changes to acknowledge that GR binding is enriched near repressed genes, albeit less so than for activated genes. In addition, we will include genome browser shots of genes with nearby peaks that are repressed by GR.
- Concept of naïve cells (Fig. 3A). If cells are normally grown in serum-containing media, which is known to have some level of steroids, can the cells described here as "Basal expression" be truly free of a primed state? In the first part of the experimental design (+/- 4h hormone), which type of media is present here? Is it 5% FBS? A concern is that the authors may require the assumption that the (4h + 24h) period a is sufficient to erase all memory of the cells, which is exactly what they are trying to test.
See our response to the first major comment above.
It would be interesting to do a time course of the hormone-free period of the washout to determine the memory of the chromatin environment that results in the enhanced transcriptional response instead of just 24 and 48 hrs in A549 cells.
We agree that that would be interesting but this is something that we cannot include for now.
Fig 5A appears to show H3K27ac overlaying H3K27me marks near the promoter of ZBTB16 and at the GR sites within the gene locus with no reduction in H3K27me levels. This seems counterintuitive and should be explained or addressed especially since the authors use quantitative comparisons of H3K27ac levels with and without treatment in other figures.
A trivial explanation for the overlaying H3K27ac and H3K27me3 marks at the ZBTB16 locus is that the ChIP results represent a population average. From our single-cell FISH experiments, we found that only a subset of cells activates ZBTB16 expression upon hormone treatment. Thus, a potential explanation is that the cells of the population that respond are responsible for the H3K27ac signal whereas the non-responders are decorated with H3K27me3. We will include this information in a revised discussion.
Showing the changes of ZBTB16 upon 2nd stimulation via FISH is not terribly surprising and is even the most expected reason for higher RNA levels. Why does it only occur at that gene is a better question and is touched on in the discussion. It is more likely that this gene has a very low level of pre-hormone transcription compared to FKBP5 (see Fig 3e and the FISH images). ZBTB16 is in the lower 3rd of basemean RNA levels of GR responsive genes according to the RNAseq data. Selection of 1 or 2 other genes with similar basemean levels of RNA (from the RNA-Seq data) would make the data more
When compared to FKBP5, ZBTB16 indeed has very low levels of pre-hormone expression. However, this is unlikely to explain the observed “memory” for ZBTB16 given that there are other genes with similarly low pre-hormone levels that do not show more robust responses upon repeated hormone exposure (see Fig. 3B,D). For the FISH experiments, we decided to include a non-primed gene (FKBP5 as control). We agree that adding additional control genes with comparable basemean levels would be informative. For example, this would tell us if a response of only a subset of cells in the population to hormone is specific to ZBTB16. Based on single cell studies by others (PMID: 32170217), most GR target genes show a response in only a subset of cells indicating that this is unlikely a unique feature of ZBTB16 explaining the priming observed. Rather than performing additional experiments, we will revise the discussion to acknowledge the difference in basemean and the potential role of cell-to-cell variability in explaining the observed “memory” for the ZBTB16 gene.
**Minor comments:**
- In the Intro (paragraph two), the authors explain the different mechanisms by which GR might repress genes. One alternative the authors appear to have missed is the possibility of direct binding to GREs while, for example, recruiting a selective corepressor such as GRIP1 (Syed et al., 2020). There are many recent critics to the notion that transrepression via tethering is responsible for GR repressive actions at all (Escoter-Torres et al., 2020; Hudson et al., 2018; Weikum et al., 2017).
We are aware of these studies and agree that they should be included when listing the possible mechanisms by which GR can repress genes. We will revise the text accordingly.
- When the authors introduce the concept of tethering to AP-1, they go way back to the first description of tethering. However, one of the references (Ref 20) actually goes against the tethering model as they did not detect protein-protein interactions between AP-1 and GR, and also, they conclude that repression requires the DNA-binding domain.
We will pick a more appropriate reference indicative of tethering as a mechanism by which GR might repress genes.
-Figure 2. The authors state "This suggests that the few sites with persistent opening are likely a simple consequence of an incomplete hormone washout and associated residual GR binding". The authors should check the subcellular distribution of GR after their washout protocol. If the washout is not completed, GR should still be in the nuclear compartment.
The careful phrasing here was to include the possibility that GR might bind DNA even when hormone is completely washed out. However, a more likely explanation is that the washout is incomplete. The residual GR binding we find in our ChIP assays shows us that a subset of GR is indeed still chromatin-bound which implies that some GR is still in the nuclear compartment.
- The first part of the manuscript (Repression through "squelching") seems a bit disconnected from the rest of the results (reversibility in accessibility). The abstract is structured in a way that this disconnection seems much less obvious. Perhaps the authors could try to present their squelching part in the middle of the manuscript, following the flow of the abstract? This is just a suggestion.
When revising the manuscript, we will see if implementiung this suggestion is feasible.
- Figures have CAPS panel letters (A,B,C, etc) while the text calls for lower case letter (a,b,c...)
We will fix this as part of the revision.
Reviewer #2: **Major Comments**
We agree that long-term and repeated GC treatment would be very interesting to study and would yield insights that are more likely to be relevant to, for example, emerging GC-resistance during therapeutic use. We are aware of the limitations of our study and will make sure that these are acknowledged in the revised manuscript and we will point out the speculative nature of translating our findings to an in-vivo setting.
2a.) The authors show several heatmaps to indicate changes in accessibility, H3K27ac and P300 upon Dex treatment as well as GR binding patterns in Fig. 1 and S1. Those are sorted by decreasing signal strength (I assume). To make those results more comparable, I suggest to sort them all in the same way (e.g. by descending ATAC-Seq signal or fold-change).
A similar suggestion was made by reviewer 1. We agree that using the same sort order for the datasets makes it easier to link the different types of data we generated. We will present the data with a consistent sorting order and stratified by GR-occupied or not when we revise the manuscript.
2b.) In line with a.), it is unclear to the reader if those sides opening /closing are the same sides showing increased/decreased H3K27ac or P300 occupancy and if those sides bind GR. Integrating this data together with mRNA e.g as correlation plots would strengthen the author's argument that accessibility, H3K27ac and mRNA changes are indeed correlated. What about the GR binding sites that do not change accessibility or H3K27ac? What makes those different? **
Therefore, the statement "Furthermore, closing peaks, which show GC-induced loss of H3K27ac levels and lack GR occupancy (Fig. S1c-f), were enriched near repressed genes" on page 10 as well as the statement "suggesting that transcriptional repression by GR does not require nearby GR binding." in the abstract and discussion cannot be made from how the data is presented.
The first issue raised will be addressed by using the same sort order across different types of data. It might also shed light on features associated with GR binding sites that do not change accessibility or H3K27ac. Once we implement the revised sorting order, we will evaluate if the statements mentioned are indeed supported by the data.
2c.) Several recent studies have shown that GR's effect on gene expression and chromatin modification at enhancers might be locus-/context-specific ("tethering", competition, composite DNA binding) and/or recruitment of different co-regulators (see Sacta et al. 2018 (doi: 10.7554/eLife.34864), Gupte et al. 2013 (doi.org/10.1073/pnas.1309898110) and many more). Defining the GR-bound or opening/closing sides in terms of changing H3K27ac (or having H3K27ac or not) more closely would help to link those to gene expression changes e.g. in violin plots.
Furthermore, the authors could include a motif analysis to see if the different enhancer behaviours can be explained by differences in the GR motif sequence or co-occurring motifs. Thereby more closely defining the mechanism of chromatin closure a sites that lack GR binding e.g. by displacement of other transcription factors as described for p65 in macrophages (Oh et al. 2017 (doi.org/10.1016/j.immuni.2017.07.012)).
In general a more detailed analysis of the data is required before the authors could state "Instead, our data support a 'squelching model' whereby repression is driven by a redistribution of cofactors away from enhancers near repressed genes that become less accessible upon GC treatment yet lack GR occupancy." on page 10. The results might also be explained by competitive transcription factor binding, tethering or selective co-regulator recruitment (e.g. HDACs).
We will include a motif analysis comparing opening, closing and non-changing sites (stratified into GR-occupied or not) in a revised version of the manuscript. In addition, we will further investigate the redistribution of p300 upon Dex-treatment e.g. to test the correlation between p300 loss at closing sites lacking GR occupancy and transcriptional repression. We agree that the “squelching model” is just one of several explanations for repression and will provide a more comprehensive list of possible explanations beyond squelching as part of the revision.
We will discuss the difference in receptor levels between the cell lines, the different number of genomic GR binding sites and its possible implication in the observed residual binding after wash-out in U2OS-GR cells as suggested.
We agree that the coverage plots do not take the fraction of binding sites with signal into account. However, by also showing the heat maps, this information is also available to the reader. In our opinion, the coverage plots provide a straight-forward way to compare the signal for the different categories of peaks. The violin plots are an interesting alternative way to present the data, which also captures the diversity in the signal within each group. We will add violin plots to the supplementary data as requested.
We see your point. However, based on the ATAC-signal (Fig. 5D) the changes in nucleosomal occupancy upon GC treatment are the same for naiive and primed cells and revert to their base-line level after hormone withdrawal. This indicates that these loci have comparable nucleosome occupancy after wash-out. Yet, the levels for these histone modifications do not differ between primed and naiive cells indicating that these histone marks do not “mark” the promoter of primed genes after wash-out.
We are reluctant to put p-values on every chart, especially for experiments with few replicates. Importantly, we always plot the values for each individual data point, so the reader can gage if they differ between conditions. We will add p-values for figure 4 to test (support) our claim that ZBTB16 is primed whereas other GR target genes are not.
A similar suggestion was brought up by reviewer #1, here is the response we gave to this comment: When compared to FKBP5, ZBTB16 indeed has very low levels of pre-hormone expression. However, this is unlikely to explain the observed “memory” for ZBTB16 given that there are other genes with similarly low pre-hormone levels that do not show more robust responses upon repeated hormone exposure (see Fig. 3B,D). For the FISH experiments, we decided to include a non-primed gene (FKBP5 as control). We agree that adding additional control genes with comparable basemean levels would be informative. For example, this would tell us if a response of only a subset of cells in the population to hormone is specific to ZBTB16. Based on single cell studies by others (PMID: 32170217), most GR target genes show a response in only a subset of cells indicating that this is unlikely a unique feature of ZBTB16 explaining the priming observed. Rather than performing additional experiments, we will revise the discussion to acknowledge the difference in basemean and the potential role of cell-to-cell variability in explaining the observed “memory” for the ZBTB16 gene.
The fact that we do not observe elevated expression of other genes upon repeated expression could be due to the relatively short length of the hormone treatment, 4 hours, which was chosen to enrich for direct target genes of GR. These four hours might be insufficient for transcription, translation and ultimately gene regulation by the ZBTB16 protein. We have not looked at ZBTB16 protein levels.
**Minor Comments**
We will include this information in a revised version of the manuscript.
We will add the requested peak-centric view. Based on a previous study (PMID: 29385519), we expect that binding is a poor predictor of gene regulation of nearby genes, especially for repressed genes.
In our analysis, we looked at opening and closing peaks independently. If a peak is in the vicinity of multiple genes, it will only be assigned to the closest one. Thus, genes that have both and opening and a closing peak in the 50kb window will be included in both the analysis of closing sites and opening sites. We have not looked at clusters of binding sites, but agree that this would be interesting to see if the combinatorial action of multiple peaks makes regulation of the gene more likely. We will look into this during the revision process.
- The authors claim on p10 that "We could validate several examples of opening and closing sites and noticed that opening sites are often GR-occupied whereas closing sites are not occupied by GR". As most of the ChIP-Seq experiments were performed on formaldehyde-only fixed cells, the authors might miss "tethered" sides, which are mostly linked to gene repression. You might rephrase this part to most closing sites lack direct DNA binding.
Even though several studies indicate that tethered binding can be captured using formaldehyde-only fixed cells (e.g. PMID: 32619221, PMID: 15879558), we agree that the ChIP-assay might have blind spots, for instance for tethered binding, and will revise our statements as suggested.
This might be related to comment #4 given that P300 is brought to the DNA by other transcription factors whereas H3K27ac is directly DNA-bound which likely influences the cross-linking efficiency. By resorting the heat-maps, we will be able to determine the overlap between p300 recruitment and changes in H3K27ac levels (the other main enzyme that deposits this mark is CREBBP (a.k.a. CBP)).
We will include this information in a revised version of the manuscript.
We have not looked into this but a previous study by the Reddy lab (PMID: 22801371) has investigated binding sites in A549 cells that are occupied at very low Dex concentrations. They found that this is not driven by a specific GR motif but rather by the presence of binding sites for other transcription factors and chromatin accessibility.
This data for the GILZ gene is shown in Figure S2C. When we revise the manuscript we will add this information to main figures 1 and 2 as suggested.
This is shown in figure S3C and shows that expression levels of certain genes (ZBTB16 and FKBP5 but not GILZ) stay high after Dex washout (but not cortisol wash-out) consistent with persistent GR binding at a subset of GR-occupied loci for the experiments using Dex.
For both S2C and S3C, cells were treated for 4h with Dex before the wash-out. For the ZBTB16 and FKBP5 genes, the persistent GR binding after wash-out is accompanied by a preserved Dex response after wash-out. For GILZ, GR binding at one of the peaks near the GILZ gene is also preserved, yet the expression of this gene reverses to its pre-treatment levels after wash-out. A possible explanation is that the residual binding at the GILZ gene is observed for only one of several nearby GR peaks. Previous studies, where we deleted GR binding sites near the GILZ gene, have shown that the combined action of multiple GR-occupied regions is needed for robust induction of this gene (PMID: 29385519).
A trivial explanation for the overlaying H3K27ac and H3K27me3 marks at the ZBTB16 locus is that the ChIP data represents a population average. From our single-cell FISH experiments, we found that only a subset of cells activates ZBTB16expression upon hormone treatment so a potential explanation is that the cells of the population that respond are responsible for the H3K27ac signal whereas the non-responders are decorated with H3K27me3. We will include this information in a revised discussion. On a single histone, H3K27me3 and H3K27ac are mutually exclusive. However, given that a nucleosome has 2 copies of histone H3, both modifications can in principle co-exist.
We’re guessing here, but we assume the reviewer refers to the potentially slightly higher H3K27me3 levels upon Dex treatment for ChIP-seq whereas the qPCR indicates that the levels do not change? The change seen in the ChIP-seq experiment is marginal and based on a single experiment. In contrast, the qPCR data shows the results from three biological replicates and therefore is probably a more reliable source of information.
We will include this information in a revised version of the manuscript.
Cancer cell lines often have variable karyotypes and our FISH data suggests that the ZBTB16 locus is present in more than 2 copies in some of the A549 cells. Here’s the info from the ATCC website describing the karyotype of A549 cells: …” This is a hypotriploid human cell line with the modal chromosome number of 66, occurring in 24% of cells. Cells with 64 (22%), 65, and 67 chromosome counts also occurred at relatively high frequencies; the rate with higher ploidies was low at 0.4%.....”.
Upon quick inspection, we find that GR target genes are typically not marked by H3K27me3, however ZBTB16 does not appear to be the only one. When we revise the manuscript, we will look more systematically at the link between gene regulation by GR and genes marked by H3K27me3 to determine how “special” the presence of this mark is, which will also inform us about the likelihood that it is linked to the transcriptional memory observed for the ZBTB16 gene.
We are not sure if ZBTB16 regulation by GR is tissue independent. However, in contrast to most GR target genes that are regulated in a cell-type-specific manner, ZBTB16 is regulated in both cell lines we examined and has also been reported to be a GR target gene in other cell types e.g. in macrophages (PMID: 30809020).
Reviewer #3 **Major Comments:**
For sure the washout time matters and we do not doubt that the persistent changes observed upon shorter wash-out by the Hager lab are real. One of the reasons we chose the 24h period was to see if the changes observed by Lightman and Hager might persist for extended periods of time as suggested by Zaret and Yamamoto. Our findings suggest that this is not the case and that the majority of GR-induced changes are short-lived. Perhaps future studies can shed light on how long changes persist. However, given the slow dissociation between GR and Dex, we expect that it might be hard to dissect if persistent changes are indeed persisting in the absence of GR binding or reflect an incomplete hormone wash-out.
The objective of this study was to find out if persistent changes as observed in Ref33 are the exception or the rule not to test if the original observation is correct (importantly, another cell line was used in Ref33 which makes a 1:1 comparison impossible to begin with). We believe that we have convincingly shown that, for the cell lines we assayed, persistent changes are rare if occurring at al. Given that no convincing persistent changes were observed after a 24h washout, we think that it is very unlikely that such changes would be observable after even longer wash-out periods. We do not intend to include experiments using longer wash-out but will revise the discussion to emphasize that the lack of persistent changes we found might be specific to the cell lines we chose for our studies.
We agree that adding this percentage is a good idea as this would allow for a more quantitative comparison between the different groups. Here are the numbers:
A549 cells:
opening sites: 49%
closing: 10%
nonchanging: 18%
U2OS cells:
opening: 54%
closing: 0.2%
nonchanging: 7%
We will include this information in a revised version of the manuscript.
For the ATAC-seq experiments, we treated the dex-treated and cort-treated experiments as replicates to find candidate regions with persistent chromatin changes. For the ATAC-seq data, a site is 'persistent' if called (by MACS2, e.g. DEX vs EtOH) upon treatment and then again 24h after washout (DEX washout vs EtOH washout). For the ATAC-qPCR experiments, we performed 4 biological replicates and will perform a t-test to determine if the small difference we observe at some sites between the EtOH and washout is statistically significant. Given the overlapping error bars and the very small difference, don’t expect the difference to be significant even for these most promising candidates from our genome-wide analysis.
Indeed we did not find a mechanistic explanation for the ZBTB16-specific memory. Possible explanations are discussion in the following section of the results (page 14-15): “… Mirroring what we say in terms of chromatin accessibility, transcriptional responses also seem universally reversable with no indication of priming-related changes in the transcriptional response to a repeated exposure to GC for any gene with the exception of ZBTB16. Although several changes in the chromatin state occurred at the ZBTB16 locus, none of these changes persisted after hormone washout arguing against a role in transcriptional memory at this locus (Fig. 5). Similarly, the increased long-range contact frequency between the ZBTB16 promoter region and a GR-occupied enhancer does not persist after washout (Fig. 5e). Notably, our RNA FISH data showed that ZBTB16 is only transcribed in a subset of cells, hence, it is possible that persistent epigenetic changes occurring at the ZBTB16 locus also only occur in a small subset of cells and could thus be masked by bulk methods such as ChIP-seq or ATAC-seq. Another mechanism underlying the priming of the ZBTB16 gene could be a persistent global decompaction of the chromatin as was shown for the FKBP5 locus upon GR activation [35]. Likewise, sustained chromosomal rearrangements, which we may not capture by 4C-seq, could occur at the ZBTB16 locus and affect the transcriptional response to a subsequent GC exposure. Furthermore, prolonged exposure to GCs (several days) can induce stable DNA demethylation as was shown for the tyrosine aminotransferase (Tat) gene [71]. The demethylation persisted for weeks after washout and after the priming, activation of the Tat gene was both faster and more robust when cells were exposed to GCs again [71]. Interestingly, long-term (2 weeks) exposure to GCs in trabecular meshwork cells induces demethylation of the ZBTB16 locus raising the possibility that it may be involved in priming of the ZBTB16 gene [72]. However, it should be noted that our treatment time (4 hours) is much shorter. Finally, enhanced ZBTB16 activation upon a second hormone exposure might be the result of a changed protein composition in the cytoplasm following the first hormone treatment. In this scenario, increased levels of a cofactor produced in response to the first GC treatment would still be present at higher levels and facilitate a more robust activation of ZBTB16 upon a subsequent hormone exposure. Although several studies have reported gene-specific cofactor requirements [73], the 14 fact that we only observe priming for the ZBTB16 gene would make this an extreme case where only a single gene is affected by changes in cofactor levels……”.
**Minor Comments**
We will include a motif analysis for opening and closing sites in a revised version of the manuscript.
We will revise the label in a revised version of the manuscript as suggested.
We actually prefer the MA plots as they also provide information regarding the basemean counts for regulated genes. This allows one, for example, to see that other GR-regulated genes with similar basemean counts do not show a “memory” suggesting that the low expression level for ZBTB16 likely does not explain the observed priming.
We will include this information in a revised version of the figure.