Glucocorticoid stimulation induces regionalized gene responses within topologically associating domains

Transcription-factor binding to cis-regulatory regions regulates the gene expression program of a cell, but occupancy is often a poor predictor of the gene response. Here, we show that glucocorticoid stimulation led to the reorganization of transcriptional coregulators MED1 and BRD4 within topologically associating domains (TADs), resulting in active or repressive gene environments. Indeed, we observed a bias toward the activation or repression of a TAD when their activities were defined by the number of regions gaining and losing MED1 and BRD4 following dexamethasone (Dex) stimulation. Variations in Dex-responsive genes at the RNA levels were consistent with the redistribution of MED1 and BRD4 at the associated cis-regulatory regions. Interestingly, Dex-responsive genes without the differential recruitment of MED1 and BRD4 or binding by the glucocorticoid receptor were found within TADs, which gained or lost MED1 and BRD4, suggesting a role of the surrounding environment in gene regulation. However, the amplitude of the response of Dex-regulated genes was higher when the differential recruitment of the glucocorticoid receptor and transcriptional coregulators was observed, reaffirming the role of transcription factor-driven gene regulation and attributing a lesser role to the TAD environment. These results support a model where a signal-induced transcription factor induces a regionalized effect throughout the TAD, redefining the notion of direct and indirect effects of transcription factors on target genes.

Introduction

Gene regulation is controlled by the recruitment of transcriptional regulators to cis-regulatory regions and implicates an intricate interplay with the chromatin environment and the chromosome architecture (Lee and Young, 2013; Kim and Shendure, 2019; Schoenfelder and Fraser, 2019). In fact, the environment surrounding a gene is an important determinant of its transcriptional regulation as exemplified by the gene position effect (Elgin and Reuter, 2013). Furthermore, the use of transgenes in animal models highlighted that the integration site is a major determinant of expression levels. Random integration of transgenes in the mouse genome confirmed large differences in transcription levels depending on the integration site (Akhtar et al., 2013). In fact, the transgene typically adopts the chromatin environment of the integration site, which will influence gene regulation (Akhtar et al., 2013; Despang et al., 2019). However, why and how this process takes place remains poorly understood.

In higher eukaryotes, genes are clustered not only in the linear genome but also in the three-dimensional genome (Bonev and Cavalli, 2016; Rowley and Corces, 2018; Schoenfelder and Fraser, 2019; Misteli, 2020). Indeed, physical interactions between genomic regions create a multilevel structure. At high levels, compartments A and B segregate actively transcribed regions from repressed regions (Lieberman-Aiden et al., 2009; Dixon et al., 2012; Wang et al., 2016). At a smaller scale, topologically associating domains (TADs) represent self-interacting regions favoring contacts between cis-regulatory regions and genes (Lieberman-Aiden et al., 2009; Nora et al., 2012; Wang et al., 2016). In addition to physical proximity between genes, there is an expanding body of evidence suggesting that genes are transcriptionally co-regulated in response to external signals (Osborne et al., 2004; Fanucchi et al., 2013). In fact, genes belonging to the same TAD are often co-expressed (Le Dily et al., 2014; Boudaoud et al., 2017; Hsieh et al., 2020; Krietenstein et al., 2020). For instance, activated and repressed genes following steroid stimulation are often segregated into distinct TADs (Le Dily et al., 2014; Vockley et al., 2016; D’Ippolito et al., 2018). In addition to TADs, smaller structures such as sub-TADs (Rao et al., 2014; Rowley et al., 2017; Huang et al., 2019) and co-expression domains (Osborne et al., 2004; Soler-Oliva et al., 2017) have also been supporting the notion of gene co-regulation. Yet, the mechanisms driving gene co-regulation in response to external signals are still unclear.

Transcription factors modify the transcriptional program in response to external signals by converging toward sequence-specific DNA motifs within cis-regulatory regions. Accordingly, they play major roles in normal and disease development by establishing and maintaining cell states (Lee and Young, 2013). Across the genome, transcription factors bind in clusters and synergize to control the transcriptional program (Spitz and Furlong, 2012; Voss and Hager, 2015; Liu and Tjian, 2018). Transcription factors and their coregulators form condensates, creating dynamic environments surrounding active genes (Hnisz et al., 2017; Boija et al., 2018; Cho et al., 2018; Chong et al., 2018; Sabari et al., 2018; Zamudio et al., 2019). Furthermore, it is generally accepted that physical proximity between enhancer and promoter regions is responsible for transmitting the effect of distally bound transcription factors to promoter regions (Deng et al., 2012; 2014; Lee et al., 2015; Chen et al., 2018). However, transcriptional changes are often observed without evidence of direct occupancy by the induced transcription factors or proximity of the cis-regulatory regions with the gene promoter (Alexander et al., 2019; Benabdallah et al., 2019; Schoenfelder and Fraser, 2019). These observations suggest that transcription factors create a regionalized impact on a gene domain through modulation of coregulators and/or the chromosome architecture.

Signal-induced transcription factor stimulations, including steroid nuclear receptors, are associated with a global redistribution of transcriptional coregulators, leading to the activation or repression of specific genes (López-Maury et al., 2008; Voss and Hager, 2015). Among them, the glucocorticoid receptor (GR, NR3C1) leads to the rapid activation and repression of target genes (Reddy et al., 2009). Interestingly, gene activation and repression are not always correlated with GR DNA occupancy (Thormann et al., 2018), suggesting the implication of alternative mechanisms involved in gene regulation. Considering that the GR interacts with a large number of transcriptional coregulators, activation and repression mechanisms were suggested to involve competition mechanisms with other transcription factors (Schmidt et al., 2016). Since, GR-responsive genes are spatially connected through pre-stimulation interactions (D’Ippolito et al., 2018; Portuguez et al., 2022), an interesting possibility is that signal-induced transcription factors, when activating and repressing, modulate the microenvironment within a gene domain leading to co-regulated transcriptional changes.

Here, we show that glucocorticoids, through the GR, induce a reorganization of transcriptional coregulators MED1 and BRD4. These regionalized changes in gene regulation are secluded within TADs, showing gains and losses in MED1 and BRD4. The presence of a gene within a changing TAD is an important determinant of the gene response. Our model proposes that regionalized gene regulation within TADs is a direct consequence of the modulation of the levels of transcriptional coregulators in response to glucocorticoids.

Materials and methods

Cell culture

A549 (ATCC, CCL-185) cells were grown in F12K medium (Gibco, 21127022). The culture medium was supplemented with 10% fetal bovine serum (Invitrogen, qualified 12483020), 100 μM MEM non-essential amino acids (Cellgro, 25-0250), 2 mM L-glutamine (Gibco, 25030-081), 100 U/mL penicillin, and 100 μg/mL streptomycin (Gibco, 15170-063). For Dex (Sigma, D1756) stimulations, cells were kept for 3 days in phenol red-free DMEM (Corning, #17-205-CV) supplemented with 5% charcoal/dextran-stripped fetal bovine serum (Fisher, #SH3006803), 100 μM MEM non-essential amino acids, 2 mM L-glutamine, 100 U/mL penicillin, and 100 μg/mL streptomycin. A Dex concentration of 100 nM was used for 60 min.

ChIP-seq

ChIP-seq experiments were performed as described previously (Bilodeau et al., 2009; Kagey et al., 2010; Fournier et al., 2016; Boudaoud et al., 2017). In brief, 50 million cells were cross-linked for 10 min with 1% formaldehyde and quenched with 125 mM glycine for 5 min. Cells were then washed with PBS, pelleted, flash frozen, and stored at −80°C. Sonicated DNA fragments were immunoprecipitated with antibodies directed against MED1 (Bethyl Laboratories, A300-793A) and BRD4 (Bethyl Laboratories, A301-985A50). Library preparation and high-throughput sequencing were performed at the McGill University and Génome Québec Innovation Centre (MUGQIC), Montréal, Canada. Previously published GR ChIP-seq raw data were retrieved from the ENCODE portal. Analysis of raw sequencing reads was performed using the standard analysis pipelines for quality control, enrichment quantification, and visualization from the Canadian Center for Computational Genomics (Bourgey et al., 2019) for ChIP-seq analysis (version 4.1.2). In brief, reads were trimmed using Trimmomatic (Bolger et al., 2014). High-quality reads were aligned to the human reference genome (hg38) using the BWA (Li and Durbin, 2009) aligner. PCR duplicates were removed using Picard’s MarkDuplicates (http://broadinstitute.github.io/picard/). Narrow peaks were called using MACS2 (Zhang et al., 2008) callpeak with the following options: --nomodel --gsize 2479938032.8 and supplying the sequenced corresponding input DNA as the background control. Peaks overlapping with ENCODE DAC exclusion list regions (Amemiya et al., 2019) (accession number ENCSR636HFF) and belonging to non-standard chromosomes were removed. Analysis of read coverages per bins was performed using multiBamSummary from deepTools 3.5.1 (Ramírez et al., 2016). To identify differentially occupied regions, we used the R package DiffBind (Ross-Innes et al., 2012) [using the DESeq2 method (Love et al., 2014)] with the following parameters: summits = TRUE, minOverlaps = 2, and adjusted p-value (false discovery rate) < 5%. The volcano plot was generated using the R package EnhancedVolcano. Heatmaps were generated using computeMatrix and plotHeatmap from deepTools 3.5.1. To generate genomic visualizations, reads from BAM files were normalized using reads per kilobase per million mapped reads (RPKM), extended to 225 bp using the bamCoverage function from deepTools 3.5.1 (Ramírez et al., 2016), and uploaded to the University of California, Santa Cruz (UCSC) Genome Browser (Kent et al., 2002). To identify enriched DNA-binding motifs, the CentriMo tool (version 5.5.3) (Bailey and MacHanick, 2012) from MEME suite was used with --local --score 5.0 --ethresh 10.0 parameters as the input and the JASPAR database for non-redundant transcription factor binding sites in eukaryotes (Castro-Mondragon et al., 2022) as reference.

Genomic region–gene associations were performed using the Genomic Regions Enrichment of Annotations Tool implemented in the R package rGREAT (McLean et al., 2010; Gu and Hübschmann, 2023) using “basal plus extension” as a gene regulatory domain definition and a maximum extension of 10 kb. To assign distal cis-regulatory regions to genes, genomic regions were integrated with chromatin interactions from Hi-C data (described in the following section). The ChromHMM 18-state model dataset from human A549 cells following 1 h Dex stimulation was retrieved from the ENCODE consortium website (accession number ENCSR931PHX). The percentages were computed as the proportion of nucleotides overlapping with each chromatin state. The 18 chromatin states are as follows: active transcription start site (TssA), flanking transcription start site (TssFlnk), upstream flanking transcription start site (TssFlnkU), downstream flanking transcription start site (TssFlnkD), strong transcription (Tx), weak transcription (TxWk), genic enhancer 1 (EnhG1), genic enhancer 2 (EnhG2), active enhancer 1 (EnhA1), active enhancer 2 (EnhA2), weak enhancer (EnhWk), ZNF genes and repeats (ZNF/Rpts), heterochromatin (Het), bivalent/poised TSS (TssBiv), bivalent enhancer (EnhBiv), repressed Polycomb (ReprPC), weak repressed Polycomb (ReprPCWk), and Quiescent/low (Quies).

Chromosome architecture

A549 Hi-C raw data (D’Ippolito et al., 2018; McDowell et al., 2018) were retrieved from the ENCODE portal (Davis et al., 2018) and processed using the HiC-Pro pipeline version 3.1.0 (Servant et al., 2015). In brief, paired-end reads were aligned to the hg38 reference genome using Bowtie 2 (Langmead and Salzberg, 2012), and default parameters were used to remove duplicate and low-map quality reads and assign reads to MboI restriction fragments. Hi-C interaction matrices were generated at a resolution of 50, 10, and 5 kb. Significant chromatin interactions were identified using FitHic2 version 2.0.8 (Kaul et al., 2020) at 5 kb and 10 kb resolution of the interaction matrix in control and Dex (1 h) conditions. The GenomicInteractions R package (Harmston et al., 2015) was used for manipulating chromatin interaction data. TADs were identified at a resolution of 50 kb using Armatus version 2.3 (Filippova et al., 2014) with a gamma parameter of 0.8. For each TAD, the number of regions gaining and losing MED1 and BRD4 was computed to calculate the TAD score using the following formula: gain/(gain + loss).

Transcriptomic

Previously published Dex-stimulated A549 RNA-seq raw data were retrieved from the ENCODE portal and processed using the MUGQIC RNA-Seq pipeline version 4.1.2 (Bourgey et al., 2019). In brief, reads were trimmed for adaptor sequences using Trimmomatic (Bolger et al., 2014). High-quality reads were aligned to the human reference genome (hg38) using the STAR aligner (Dobin et al., 2013). PCR duplicates were removed using Picard’s MarkDuplicates (http://broadinstitute.github.io/picard/). Gene counts were determined using featureCounts (version 2.0.1) (Liao et al., 2014) with the genomic annotation Ensembl release 104. Samples considered as outliers were removed after the visual inspection of the PCA plot and assessment of the distance between samples. Differentially expressed genes were identified using DESeq2 (Love et al., 2014) and called significant when the Benjamini–Hochberg-corrected p-values were under 0.05. Upregulated genes were selected at a minimum log2 fold change > 0.75 and downregulated genes at a minimum log2 fold change < −0.75. A gene was considered activated or repressed if it was selected at least once between 0 and 6 h. A gene was defined as active if having at least one read in 60% of the samples between 0 and 6 h. Genes defined as inactive were not considered in downstream analysis.

Results

Glucocorticoid response implicates redistribution of MED1 and BRD4

To explore the generation of regional effects on gene regulation, we examined the recruitment of the bromodomain-containing protein 4 (BRD4) and the mediator complex subunit MED1, which are functional coregulators of the GR (Chen and Roeder, 2007; Pradhan et al., 2016) and are found within mobile phase separation droplets (Sabari et al., 2018). Experiments were carried out in A549 cells following a 60-min treatment with Dex. Between 8,189 and 27,310, regions were identified in the different datasets, and the consensus was determined by the overlap between replicates (Supplementary Table S1). As expected, a large proportion of regions occupied by MED1 were shared with BRD4 in control (84.4%) and Dex-stimulated cells (93.7%) (Figure 1A). Accordingly, the read coverage across the genome was highly correlated between the two transcriptional coregulators in control (Pearson’s r = 0.95 and p < 2.2e-16) and Dex-stimulated (Pearson’s r = 0.92 and p < 2.2e-16) cells (Figure 1B). Concurrent gains and losses in MED1 and BRD4 were frequently observed throughout the genome. For example, regulatory regions of ANGPTL4 and IL11 genes, which are positively and negatively Dex-regulated genes, respectively, (S1A-B Fig), showed increased and decreased levels of BRD4 and MED1 at regions bound by the GR (Figure 1C). Overall, concurrent increases and decreases at 5,442 and 1,097 regions, respectively, were identified combining the MED1 and BRD4 replicates (Figure 1D; Supplementary Table S2). These regions will be referred to as regions gaining and losing MED1 and BRD4 thereafter. To confirm that regions integrating replicates for BRD4 and MED1 properly identified differential regions for both coregulators, we generated ChIP-seq density heatmaps (Figure 1E). Increases and decreases in both MED1 and BRD4 were observed in regions defined as gaining and losing transcriptional coregulators. At these regions, the variations in MED1 and BRD4 signals following the Dex stimulation were correlated (for gains, Pearson’s r = 0.94 and p < 2.2e-16; for losses, Pearson’s r = 0.91 and p < 2.2e-16, Supplementary Figure S1C). Therefore, our results confirm that Dex stimulation reorganizes MED1 and BRD4 genome-wide.

FIGURE 1

FIGURE 1. Redistribution of MED1 and BRD4 is correlated across the genome following glucocorticoid stimulation. (A) Venn diagram representing MED1 and BRD4 consensus occupied regions in A549 cells in control (Ctrl) and Dex-stimulated conditions. (B) Correlations of genome-wide read coverage between MED1 and BRD4 ChIP-seq samples in control and Dex-stimulated conditions. Coverage calculation was performed for consecutive bins of equal size (10 kb). A log10 transformation was applied after adding a pseudocount of 1 to the coregulator read count. Pearson’s correlation was computed using the read count raw values. (C) ChIP-seq occupancy profiles of GR, MED1, and BRD4 at ANGPTL4 and IL11 genes, which are known to be activated and repressed by glucocorticoids, respectively. ChIP-seq profiles are displayed in RPKM. Gene depictions are presented below the ChIP-seq tracks. The 18-model chromatin states were retrieved from ENCODE and gathered in categories as follows: active TSS (TssA), flanking TSS (TssFlnk, TssFlnkU, and TssFlnkD), strong transcription (Tx), weak transcription (TxWk), genic enhancer (EnhG1 and EnhG2), active enhancer (EnhA1 and EnhA2), weak enhancer (EnhWk), and quiescent/low (Quies). (D) Volcano plot of differentially occupied regions combining ChIP-seq replicates for MED1 and BRD4 after a 1-h Dex stimulation. Regions gaining (n = 5,442) and losing (n = 1,097) MED1 and BRD4 are represented by red and green dots, respectively. Regions associated with ANGPTL4 and IL11 are displayed. (E) Top—density heatmaps representing MED1 and BRD4 ChIP-seq intensities at regions gaining (n = 5,442) and losing (n = 1,097) transcriptional coregulators following Dex stimulation in A549 cells. Regions were sorted in descending order based on the mean value per region. Bottom—average ChIP-seq signal for regions gaining and losing MED1 and BRD4. A region of 5 kb centered on the occupied region is displayed. (F) Overlaps of regions gaining and losing MED1 and BRD4 with the 18-model chromatin states.

To validate that the GR was directly implicated in the rearrangement of transcriptional coregulators, regions gaining and losing BRD4 and MED1 were overlaid with the GR genome-wide occupancy using the available data (D’Ippolito et al., 2018; McDowell et al., 2018). A large subset of regions gaining (97.8%) and losing (40.7%) MED1 and BRD4 was occupied by the GR (Supplementary Table S2), consistent with a global reorganization of transcriptional coregulators following the direct action of the GR. To functionally assess the genomic regions with the differential recruitment of MED1 and BRD4, we integrated information from the 18-model chromatin states (Ernst and Kellis, 2012) (Figure 1F). Most regions gaining MED1 and BRD4 were associated with enhancers (EnhA1 and EnhA2), while regions losing coregulators were primarily found at promoter (Tss and TssFlnkU) regions. As expected, gains were associated with the GR DNA-binding motif, while losses were enriched for the DNA-binding motif of the AP-1 family of transcription factors (Supplementary Figure S1D). These findings confirm the GR-associated modulation of cis-regulatory regions following Dex stimulation.

Gains and losses in MED1 and BRD4 correlate with gene level variations

To establish whether the differential recruitment of MED1 and BRD4 was associated with variations at the gene level, we used an available transcriptomic dataset of A549 cells stimulated with Dex (D’Ippolito et al., 2018; McDowell et al., 2018). A total of 1,759 differentially expressed genes (911 activated and 848 repressed) were identified between 0 and 6 h of Dex stimulation (Supplementary Figure S2A; Supplementary Table S3). Regions gaining and losing MED1 and BRD4 following Dex treatment were assigned linearly to genes using rGREAT (McLean et al., 2010; Gu and Hübschmann, 2023) or connected to a gene via chromatin interactions using publicly available Hi-C data (D’Ippolito et al., 2018) (see Materials and Methods). As expected, gains in MED1 and BRD4 were enriched in Dex-activated genes (Fisher’s test, OR = 9.3, and p-value = 3.58e-164). On the other hand, losses in transcriptional coregulators were enriched in Dex-repressed genes (Fisher’s test, OR = 4.5, and p-value = 4.2e-40). To determine if the differential recruitment of MED1 and BRD4 was associated with the amplitude of RNA-level variations, we compared fold changes. Variations in MED1 and BRD4 were correlated with variations at the RNA levels at all Dex stimulation time points (Figure 2A; Supplementary Figure S2B). In addition to differential recruitment, the number of regions gaining and losing MED1 and BRD4 could be a determining factor in the gene output. When genes were associated with two or more regions gaining or losing transcriptional coregulators MED1 and BRD4, the mean change in RNA levels was greater (Figure 2B). The presence of two or more regions with differential recruitment was associated with an increase between 1.23- and 2.12-folds and a decrease between 1.49- and 1.91-folds at the RNA level. Therefore, our results support that the differential recruitment of transcriptional coregulators MED1 and BRD4 at cis-regulatory regions is associated with variations at the gene levels.

FIGURE 2

FIGURE 2. Differential recruitment of MED1 and BRD4 is associated to gene level variations. (A) Scatter plot showing the correlation between variations in RNA levels of Dex-regulated genes and the recruitment of MED1 and BRD4 at their cis-regulatory regions after 1 h of Dex stimulation. For RNA levels, the log2 of the fold change is represented. For recruitment of transcriptional coregulators, the log2 of the read density at each differential region is represented. Regions gaining and losing MED1 and BRD4 were associated to genes based on their linear proximity or the presence on 3D chromatin interactions. Correlations were assessed using Pearson’s correlation method. (B) Time course analysis of the variations at the RNA level of genes associated with 1 or multiple (2+) regions with the differential recruitment of transcriptional coregulators. For RNA levels, the average of the log2 of the fold change for each gene within the group and the standard error of the mean (error bars) are represented. Mann–Whitney U tests were used to compare means between genes associated with 1 or multiple regions at 0.5, 1, 2, 3, 4, 5, and 6 h. The Benjamini–Hochberg procedure was applied on the empirical p-values to correct for multiple testing. ns, not significant; *p < 0.05 and **p < 0.01.

Regional gains and losses in MED1 and BRD4 in response to Dex

The chromosomal architecture is known to delimit different transcriptional activities within the nucleus and influence the range of action of enhancer regions (Dowen et al., 2014; Schoenfelder and Fraser, 2019). Among key factors for chromosome organization, TADs could represent physical frontiers to subdivide gene domains. Furthermore, TADs are maintained following hormonal stimulation and associated with the creation of coordinated regulatory units (Le Dily et al., 2014; D’Ippolito et al., 2018). To directly test the relationship of regions gaining and losing MED1 and BRD4 with TADs, we determined TAD boundaries at 1 h after Dex stimulation in A549 cells by processing available Hi-C data using Armatus (Filippova et al., 2014; D’Ippolito et al., 2018; McDowell et al., 2018) (Supplementary Table S4). Each region differentially occupied by MED1 and BRD4 was assigned to its corresponding TAD, and ratios of gains and losses were calculated [score = (gain)/(gain + loss)] (Figure 3A). Of the 4,957 TADs, 37.3% (1,851 TADs) were associated with at least one region with the differential recruitment of MED1 and BRD4. TAD scores were distributed following a bimodal distribution (Hartigan’s dip test and p < 2.2e-16), highlighting a bias toward homogeneity for the differential recruitment of transcriptional coregulators (Figure 3B). Indeed, TADs were biased toward either gaining (75% with a score ≥ 0.7, referred to as up) or losing (16% with a score ≤ 0.3, referred to as down) MED1 and BRD4 (Figure 3C). On average, 3.4 regions with the differential recruitment of MED1 and BRD4 were found among TADs scored as up, down, or balanced. These results support the regionalized recruitment of transcriptional coregulators following stimulation of the GR.

FIGURE 3

FIGURE 3. TADs are biased toward gaining or losing MED1 and BRD4 in response to a glucocorticoid stimulation. (A) The ratios of regions gaining and losing MED1 and BRD4 within each TAD are biased. The TAD score [= number of gains/(number of gains + losses)] was calculated for each of the 1,851 TADs with at least one region with the differential recruitment of MED1 and BRD4. Top—ranking of TADs based on the score. Values of 0 and 1 are, respectively, equal to a TAD with only regions losing or gaining MED1 and BRD4. Middle—the number of regions gaining and losing MED1 and BRD4 in each TAD is displayed. Bottom—the size of each TAD (Mb) is displayed. (B) Density plot representing the bimodal distribution of TAD scores. (C) Quantification of the number of TADs defined as up (score ≥ 0.7), balanced (0.3 < score < 0.7), and down (score ≤ 0.3).

The TAD environment surrounding a gene influences the response to Dex

While regions gaining and losing MED1 and BRD4 were strongly associated with variations in gene levels, not all Dex-regulated genes were associated with the differential recruitment of transcriptional coregulators. Indeed, the fact that 60% and 78.9% of activated and repressed genes following Dex treatment were not associated with differential coregulators raised questions about the mechanisms involved. We hypothesized that the activity of the TAD environment of a gene was an important determinant of gene level variations in addition to the differential recruitment of transcriptional regulators by the GR. The majority of 1,759 Dex-regulated genes (76%) were enriched within a TAD, gaining or losing MED1 and BRD4 (permutation test, n = 10,000, and p = 9.99e-5) (Supplementary Table S3). As expected, activated and repressed genes were enriched within TADs scored as up (Fisher’s test, OR = 5.02, and p = 7.56e-113) and down (Fisher’s test, OR = 2.07, and p = 7.79e-14). If the TAD environment surrounding a gene is an important determinant in the response to Dex, differentially expressed genes should be found within up and down TADs, matching the changes independently from the differential recruitment of transcriptional coregulators. Globally, 39.8% of Dex-regulated genes assigned to a responsive TAD were associated with a region gaining or losing MED1 and BRD4, compared to 60.2% that were not (Supplementary Table S3). Interestingly, whether the differentially expressed genes were associated to regions gaining or losing MED1 and BRD4 or not, they were enriched within a TAD matching their activity (Figures 4A, B). Indeed, chi-squared tests revealed a significant association between the direction of the gene response and the activity of the TAD, whether the gene was occupied (χ² = 214.42 and p = 2.98e-45) or not (χ² = 81.58 and p = 8.07e-17) by regions gaining or losing MED1 and BRD4 (Figures 4B, C). Analysis of Pearson residuals using a critical absolute value of 2 showed that Dex-regulated genes associated with the differential recruitment of MED1 and BRD4 or not were enriched within a TAD category matching the gene activity (Figure 4C). To validate our observations, we subdivided Dex-regulated genes into bound or not by the GR. Statistically significant relationships between the direction of the gene response and the activity of the TAD, whether the gene was bound by the GR (χ² = 163.69 and p = 2.36e-34) or not (χ² = 71.24 and p = 1.24e-14), were confirmed (Supplementary Figure S3). These results suggest that, in the absence of a measurable differential recruitment of transcriptional coregulators and GR binding, the presence of a gene within a specific TAD is an important determinant of the response.

FIGURE 4

FIGURE 4. Dex-regulated genes are enriched within TADs matching their activity independently from the differential recruitment of MED1 and BRD4. (A) Dex-regulated genes associated with the differential recruitment of MED1 and BRD4 or not were found within TADs matching their activity. Dex-activated and -repressed genes were assigned TADs and subdivided whether they were associated to a region with the differential recruitment of transcriptional coregulators (n = 532) or not (n = 805). Data are represented as a percentage of the total number of genes. (B) Representation of the 1,337 Dex-regulated genes found within a responsive TAD. RNA levels in the fold change (log2) are displayed from 1 to 6 h after Dex stimulation. Hierarchical clustering (using the Euclidean distance) was applied to the gene expression matrix and is represented by the dendogram. Regions gaining or losing MED1 and BRD4 are represented by red and green lines, respectively. TAD scores were calculated and are represented as mentioned before. (C) Association plot illustrating the dependence between changes at the RNA level for Dex-regulated genes and the category of the TAD (up, balanced, and down). The height of each bar is proportional to the Pearson residual, while the width is proportional to the square root of the expected frequency so that the area of the rectangle is proportional to the difference between observed and expected frequencies. Residual values are colored if greater than 2 (enrichment, blue) or less than −2 (depletion, red). Top—Dex-regulated genes associated with the differential recruitment of MED1 and BRD4. Bottom—Dex-regulated genes not associated with the differential recruitment of MED1 and BRD4.

The TAD environment has less impact than the recruitment of the GR and transcriptional coregulators

Our results support that the differential amount of MED1 and BRD4 is associated with changes at the gene level (Figure 2). If true, differentially expressed genes following Dex stimulation that are influenced by the environment of the TAD should be less affected compared to those associated with regions gaining and losing MED1 and BRD4. To determine if the amplitude of the differential gene regulation was equivalent between genes associated or not with regions gaining or losing MED1 and BRD4, we evaluated the trajectories of gene expression changes throughout the 6-h time course for each category of TAD (up, balanced, and down) (Figure 5). As expected, global gene activity was correlated with the TAD classification, and gene activation or repression was observed in TADs gaining or losing MED1 and BRD4, respectively. The contact with regions gaining or losing transcriptional coregulators was associated with stronger activation (between 2.73- and 6.51-fold increases) and repression (between 1.49- and 2.16-fold decreases). Similar results were obtained when differentially occupied genes were subdivided in bound or not by the GR with the exception that no differences were observed between repressed genes (Supplementary Figure S4). Taken together, these results show that while primary DNA binding by a transcription factor and the differential recruitment of transcriptional coregulators produce a stronger gene response, the surrounding TAD environment also influences the gene output.

FIGURE 5

FIGURE 5. Differential recruitment of MED1 and BRD4 translates into a stronger gene response. Gene level variations match the transcriptional activity of the TAD following Dex treatment. Genes associated to a region with the differential recruitment of MED1 and BRD4 (blue) or not (yellow) are displayed up to 6 h after treatment. For each TAD category, mean variations in RNA levels through time and the standard error of the mean (error bars) are represented. Mann–Whitney U tests were used to compare RNA level variations between genes with and without the differential recruitment of transcriptional coregulators at 0.5, 1, 2, 3, 4, 5, and 6 h. The Benjamini–Hochberg procedure was applied on the empirical p-values to correct for multiple testing. ns, not significant; *p < 0.05, **p < 0.01, and ***p < 0.001.

Discussion

Model for regionalized transcriptional regulation

Our results establish that following ligand activation, the GR elicits regionalized gains and losses in transcriptional coregulators MED1 and BRD4, creating favorable or unfavorable environments for transcriptional regulation. The primary consequence of the GR recruitment is to modulate levels of transcriptional coregulators, which not only leads directly to the activation or repression of the bound target genes but also influences other genes in the vicinity (Figures 1, 4). While we cannot exclude the possibility that genes not bound by the GR and not associated with regions gaining or losing MED1 and BRD4 are false negatives due to technical reasons, we envision two possibilities to explain how a regionalized effect on transcriptional regulation is achievable. On one hand, phase separation droplets created by the accumulation of transcription regulators could diffuse and be used as cargo between genes (Silveira and Bilodeau, 2018). Accordingly, modulation of the local concentration of transcriptional regulators could facilitate or hamper their diffusion to neighboring genes and affect transcription. This model is supported by the GR implication in the formation of biological condensates during gene activation (Stortz et al., 2020; Frank et al., 2021). On the other hand, modifying the recruitment of transcriptional regulators could impact the dynamic of the regulatory regions themselves and the frequency of contact with promoters within a TAD. This model is supported by the existence of subdomains or chromatin modules within TADs characterized by short-range enhancer–promoter and promoter–promoter interactions and the increased dynamics of regulatory regions after hormonal stimulation (Le Dily and Beato, 2018; Hsieh et al., 2020; Krietenstein et al., 2020; van Mierlo et al., 2023). The low resolution of Hi-C approaches can explain why we failed to detect physical interactions between unbound differentially expressed genes and GR-occupied regions. Be that as it may, our results globally support regional transcriptional consequences for genes in proximity to regulatory regions responding to a transcription factor, whether binding or physical interactions are detected.

Implication of the model for gene positioning effects in mammalian cells

The proposed model represents an important shift in the conception of direct and indirect effects elicited by transcription factors. Current methods to determine the primary effects of a transcription factor are based on the integration of multi-omics data. For example, chromosome architecture data are coupled to genome-wide localization of the transcription factor to determine transcription factor-bound regions and assign them to genes. Following the integration of gene expression datasets, differentially expressed genes that are unbound are typically discarded and labeled as indirect effects. Here, we are proposing that recruitment of a transcription factor and its associated coregulators has two types of primary effects on transcription: 1) a direct effect through DNA binding of the transcription factor at cis-regulatory regions and 2) a domain-dependent effect based on the position in a specific environment. It is important to note that what we are describing as domain-dependent could represent highly dynamic transcription factors or physical contacts with cis-regulatory regions making detection of enhancer–promoter contacts more difficult. The distinction between these models and the precise definition of the boundaries in which they are acting will require further investigations. Nonetheless, our interpretation is an important distinction from the current definition of indirect effects. Based on our assessment, inducing a transcription factor will redistribute transcriptional coregulators, leading to primary direct and indirect effects. With this nomenclature, gene expression changes independent of direct DNA binding or position effects would be considered as secondary.

Distinction between gene activation and repression

Transcription factors activate and repress genes by modulating recruitment of coregulators (Lambert et al., 2018). Our results support that gains in transcriptional coregulators MED1 and BRD4 are dependent on the presence of the GR on chromatin. However, for most regions where losses in MED1 and BRD4 were observed, the GR was not detected (Supplementary Table S2). While the presence of the GR at repressed regions could be harder to detect, the result also suggests that the recruitment of the GR at chromatin is not necessary to remove transcriptional coregulators. The GR was shown to titrate and sequester transcriptional coregulators, suggesting a passive mechanism of gene repression referred to as squelching (Schmidt et al., 2016; Bothe et al., 2021; Portuguez et al., 2022). Considering that genes repressed by the GR are active in basal conditions, repression can be driven by opportunity. When the concentration of the GR increases in the nucleus, interference with the transcriptional program in place would provide the transcriptional coregulators required for subsequent transcriptional activation. This interpretation is supported by the fact that higher levels of MED1 and BRD4 are found at Dex-repressed genes compared to Dex-activated genes in normal conditions (Portuguez et al., 2022). This type of passive mechanism would result in repression being distributed throughout active genes rather than being targeted. This model would explain why a limited number of regions with a significant loss in MED1 and BRD4 were identified (Supplementary Table S2). Interestingly, GR-responsive genes share spatial domains specialized in activation or repression (Portuguez et al., 2022). Whether activation domains require spatial association to a repressed domain is an open question. Interestingly, while GR has been extensively studied molecularly, no mutant with the ability to activate transcription without also repressing has previously been reported to our knowledge (the ability of mutants to repress without activating is frequent) (Beck et al., 2011). It will be interesting to determine if repression is mandatory for the ability of the GR to activate transcription.

Conclusion

The influence of the genomic microenvironment has been associated with the gene position effect in mammalian cells. We are proposing that TADs are being modulated by the redistribution of transcriptional coregulators. Therefore, when a transgene is added to different genomic contexts, access to mobile transcriptional regulators or dynamic cis-regulatory regions will differ. Furthermore, this model explains the fast kinetics of differential expression observed for genes not bound by an induced transcription factor but responding to stimulation. It remains to be determined if, during the biological response of a cell, activation and repression mechanisms are molecularly linked.

Data availability statement

The data generated for this publication is available on GEO (Edgar et al., 2002) under accession number GSE226487. The first replicate for the whole cell extract (GSM2040031) and MED1 (GSM2040033) in control conditions matching their Dex-stimulated counterparts are already available under accession number GSE76893. The ChIP-seq data are viewable on the UCSC genome browser here: https://genome.ucsc.edu/s/ckntav/A549_Dex_response. All publicly available sequencing datasets used in the manuscript are listed in Supplementary Table S5.

Author contributions

CT, ÉF, MF, and SB conceptualized the study and experiments. MF, FK, MC, MS, and F-AB-S executed the experiments. CT, ÉF, and AB performed the in silico experiments. CT, ÉF, AB, GB, AD, and SB analyzed the data. CT and SB wrote the manuscript. All authors contributed to the article and approved the submitted version.

Funding

This work was supported by funds from the Canada Research Chair in Transcriptional Genomics (Grant #950-228321 to SB), the Natural Sciences and Engineering Research Council of Canada (Grants #436266-2013 and #2019-06490 to SB), the Canadian Institutes of Health Research (Grants #286769 and #451568 to SB), doctoral research scholarships from the Fonds de Recherche du Québec (to CT and MS), and the merit scholarship program for foreign students (to MS).

Acknowledgments

The authors thank all lab members for insightful discussions and feedback.

Conflict of interest

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Publisher’s note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors, and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

Supplementary material

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fgene.2023.1237092/full#supplementary-material

SUPPLEMENTARY FIGURE S1 | Dex-dependent changes at the RNA and chromatin levels. (A) Gene level variations for ANGPTL4 and IL11 genes following Dex stimulation. Variations in RNA levels through time were determined after smoothing the trajectories using the local regression (LOESS) method. (B) RNA-seq profiles of ANGPTL4 and IL11 genes following Dex stimulation. Tracks are displayed in reads per kilobase per million mapped reads (RPKM). Gene depictions are presented below the RNA-seq tracks. (C) Scatter plot displaying the difference in normalized count values between Dex and control conditions for MED1 and BRD4. Each dot represents a genomic region defined as differentially recruiting transcriptional coregulators. Correlations were assessed using Pearson’s correlation method. (D) Motif enrichment analysis of genomic regions gaining (top) or losing (bottom) MED1 and BRD4. The top three DNA-binding motifs recovered using the JASPAR database are displayed.

SUPPLEMENTARY FIGURE S2 | Differential recruitment of MED1 and BRD4 is associated with the gene response to Dex. (A) Quantification of the number of Dex-regulated genes after 0.5, 1, 2, 3, 4, 5, and 6 h of stimulation. Activated (red) and repressed (green) genes are represented at each time point. (B) Scatter plots showing the correlation between variations in RNA levels of Dex-regulated genes and the recruitment of transcriptional coregulators at their cis-regulatory regions after 0.5, 2, 3, 4, 5, and 6 h. For RNA levels, the log2 of the fold change is represented. For recruitment of MED1 and BRD4, the log2 of the read density at each differential region is represented. Regions gaining and losing MED1 and BRD4 were associated to genes based on their linear proximity or the presence on 3D chromatin interactions. Correlations were assessed using Pearson’s correlation method.

SUPPLEMENTARY FIGURE S3 | Dex-regulated genes are enriched within TADs matching their activity independently from GR binding. (A) Dex-regulated genes bound by the GR or not were found within TADs matching their activity. Dex-activated and -repressed genes were assigned TADs and subdivided whether they were associated to a region with the differential recruitment of the GR (n = 857) or not (n = 480). Data are represented as a percentage of the total number of genes. (B) Representation of the 1,337 Dex-regulated genes found within a responsive TAD. RNA levels in fold change (log2) are displayed from 1 to 6 h after Dex stimulation. Hierarchical clustering (using the Euclidean distance) was applied to fold changes and is represented by the dendogram. Regions gaining or losing MED1 and BRD4 are represented by red and green lines, respectively. TAD scores were calculated and are represented as before. (C) Association plot illustrating the dependence between changes at the RNA level for Dex-regulated genes and the category of the TAD (up, balanced, and down). The height of each bar is proportional to the Pearson residual, while the width is proportional to the square root of the expected frequency so that the area of the rectangle is proportional to the difference between observed and expected frequencies. Residual values are colored if greater than 2 (enrichment, blue) or less than -2 (depletion, red). Top—Dex-regulated genes with GR binding. Bottom—Dex-regulated genes without GR binding.

SUPPLEMENTARY FIGURE S4 | Differential binding of the GR translates into stronger gene activation. Genes associated to GR binding (blue) or not (yellow) are displayed up to 6 h after treatment. For each TAD category, mean variation in RNA levels through time and the standard error of the mean (error bars) are represented. Mann–Whitney U tests were used to compare RNA level variations between genes with and without GR binding at 0.5, 1, 2, 3, 4, 5, and 6 h. The Benjamini–Hochberg procedure was applied on the empirical p-values to correct for multiple testing. ns, not significant; *p < 0.05 and ***p < 0.001.

SUPPLEMENTARY TABLE S1 | List of genomic regions for MED1 and BRD4 ChIP-seq datasets in control and Dex-stimulated conditions.

SUPPLEMENTARY TABLE S2 | List of genomic regions with the differential recruitment of MED1 and BRD4 following Dex stimulation.

SUPPLEMENTARY TABLE S3 | List of Dex-regulated genes.

SUPPLEMENTARY TABLE S4 | List of genomic coordinates for TADs.

SUPPLEMENTARY TABLE S5 | List of publicly available datasets used in this study.

References

Akhtar, W., Jong, J. D., Pindyurin, A. V., Pagie, L., Meuleman, W., Ridder, J. D., et al. (2013). Chromatin position effects assayed by thousands of reporters integrated in parallel. Cell. 154, 914–927. doi:10.1016/j.cell.2013.07.018

PubMed Abstract | CrossRef Full Text | Google Scholar

Alexander, J. M., Guan, J., Li, B., Maliskova, L., Song, M., Shen, Y., et al. (2019). Live-cell imaging reveals enhancer-dependent sox2 transcription in the absence of enhancer proximity. Elife 8, e41769. doi:10.7554/eLife.41769

PubMed Abstract | CrossRef Full Text | Google Scholar

Amemiya, H. M., Kundaje, A., and Boyle, A. P. (2019). The ENCODE blacklist: Identification of problematic regions of the genome. Sci. Rep. 9, 9354–9355. doi:10.1038/s41598-019-45839-z

PubMed Abstract | CrossRef Full Text | Google Scholar

Bailey, T. L., and MacHanick, P. (2012). Inferring direct DNA binding from ChIP-seq. Nucleic Acids Res. 40, e128. doi:10.1093/nar/gks433

PubMed Abstract | CrossRef Full Text | Google Scholar

Beck, I. M., Bosscher, K. D., and Haegeman, G. (2011). Glucocorticoid receptor mutants: Man-made tools for functional research. Trends Endocrinol. Metab. 22, 295–310. doi:10.1016/j.tem.2011.03.009

PubMed Abstract | CrossRef Full Text | Google Scholar

Benabdallah, N. S., Williamson, I., Illingworth, R. S., Grimes, G. R., Therizols, P., Bickmore, W. A., et al. (2019). Decreased enhancer-promoter proximity accompanying enhancer activation. Mol. Cell. 76, 473–484.e7. doi:10.1016/j.molcel.2019.07.038

PubMed Abstract | CrossRef Full Text | Google Scholar

Bilodeau, S., Kagey, M. H., Frampton, G. M., Rahl, P. B., and Young, R. A. (2009). SetDB1 contributes to repression of genes encoding developmental regulators and maintenance of ES cell state. Genes. Dev. 23, 2484–2489. doi:10.1101/gad.1837309

PubMed Abstract | CrossRef Full Text | Google Scholar

Boija, A., Klein, I. A., Sabari, B. R., Dall’Agnese, A., Coffey, E. L., Zamudio, A. V., et al. (2018). Transcription factors activate genes through the phase-separation capacity of their activation domains. Cell. 175, 1842–1855.e16. doi:10.1016/j.cell.2018.10.042

PubMed Abstract | CrossRef Full Text | Google Scholar

Bolger, A. M. A. M., Lohse, M., and Usadel, B. (2014). Trimmomatic: A flexible trimmer for illumina sequence data. Bioinformatics 30, 2114–2120. doi:10.1093/bioinformatics/btu170

PubMed Abstract | CrossRef Full Text | Google Scholar

Bonev, B., and Cavalli, G. (2016). Organization and function of the 3D genome. Nat. Rev. Genet. 17, 661–678. doi:10.1038/nrg.2016.112

PubMed Abstract | CrossRef Full Text | Google Scholar

Bothe, M., Buschow, R., and Meijsing, S. H. (2021). Glucocorticoid signaling induces transcriptional memory and universally reversible chromatin changes. Life Sci. Alliance 4, e202101080. doi:10.26508/LSA.202101080

PubMed Abstract | CrossRef Full Text | Google Scholar

Boudaoud, I., Fournier, É., Baguette, A., Vallée, M., Lamaze, F. C., Droit, A., et al. (2017). Connected gene communities underlie transcriptional changes in Cornelia de Lange Syndrome. Genetics 207, 139–151. doi:10.1534/genetics.117.202291

PubMed Abstract | CrossRef Full Text | Google Scholar

Bourgey, M., Dali, R., Eveleigh, R., Chen, C., Letourneau, L., Fillon, J., et al. (2019). GenPipes: An open-source framework for distributed and scalable genomic analyses. Gigascience 8, giz037–11. doi:10.1093/gigascience/giz037

PubMed Abstract | CrossRef Full Text | Google Scholar

Castro-Mondragon, J. A., Riudavets-Puig, R., Rauluseviciute, I., Berhanu Lemma, R., Turchi, L., Blanc-Mathieu, R., et al. (2022). Jaspar 2022: The 9th release of the open-access database of transcription factor binding profiles. Nucleic Acids Res. 50, D165–D173. doi:10.1093/nar/gkab1113

PubMed Abstract | CrossRef Full Text | Google Scholar

Chen, W., and Roeder, R. G. (2007). The Mediator subunit MED1/TRAP220 is required for optimal glucocorticoid receptor-mediated transcription activation. Nucleic Acids Res. 35, 6161–6169. doi:10.1093/nar/gkm661

PubMed Abstract | CrossRef Full Text | Google Scholar

Chen, H., Levo, M., Barinov, L., Fujioka, M., Jaynes, J. B., and Gregor, T. (2018). Dynamic interplay between enhancer–promoter topology and gene activity. Nat. Genet. 50, 1296–1303. doi:10.1038/s41588-018-0175-z

PubMed Abstract | CrossRef Full Text | Google Scholar

Cho, W.-K., Spille, J.-H., Hecht, M., Lee, C., Li, C., Grube, V., et al. (2018). Mediator and RNA polymerase II clusters associate in transcription-dependent condensates. Science 361, 412–415. doi:10.1126/science.aar4199

PubMed Abstract | CrossRef Full Text | Google Scholar

Chong, S., Dugast-Darzacq, C., Liu, Z., Dong, P., Dailey, G. M., Cattoglio, C., et al. (2018). Imaging dynamic and selective low-complexity domain interactions that control gene transcription. Science 361, eaar2555–9. doi:10.1126/science.aar2555

PubMed Abstract | CrossRef Full Text | Google Scholar

D’Ippolito, A. M., McDowell, I. C., Barrera, A., Hong, L. K., Leichter, S. M., Bartelt, L. C., et al. (2018). Pre-established chromatin interactions mediate the genomic response to glucocorticoids. Cell. Syst. 7, 146–160.e7. doi:10.1016/j.cels.2018.06.007

PubMed Abstract | CrossRef Full Text | Google Scholar

Davis, C. A., Hitz, B. C., Sloan, C. A., Chan, E. T., Davidson, J. M., Gabdank, I., et al. (2018). The Encyclopedia of DNA elements (ENCODE): Data portal update. Nucleic Acids Res. 46, D794–D801. doi:10.1093/nar/gkx1081

PubMed Abstract | CrossRef Full Text | Google Scholar

Deng, W., Lee, J., Wang, H., Miller, J., Reik, A., Gregory, P. D., et al. (2012). Controlling long-range genomic interactions at a native locus by targeted tethering of a looping factor. Cell. 149, 1233–1244. doi:10.1016/j.cell.2012.03.051

PubMed Abstract | CrossRef Full Text | Google Scholar

Deng, W., Rupon, J. W., Krivega, I., Breda, L., Motta, I., Jahn, K. S., et al. (2014). Reactivation of developmentally silenced globin genes by forced chromatin looping. Cell. 158, 849–860. doi:10.1016/j.cell.2014.05.050

PubMed Abstract | CrossRef Full Text | Google Scholar

Despang, A., Schöpflin, R., Franke, M., Ali, S., Jerković, I., Paliou, C., et al. (2019). Functional dissection of the Sox9–Kcnj2 locus identifies nonessential and instructive roles of TAD architecture. Nat. Genet. 51, 1263–1271. doi:10.1038/s41588-019-0466-z

PubMed Abstract | CrossRef Full Text | Google Scholar

Dixon, J. R., Selvaraj, S., Yue, F., Kim, A., Li, Y., Shen, Y., et al. (2012). Topological domains in mammalian genomes identified by analysis of chromatin interactions. Nature 485, 376–380. doi:10.1038/nature11082

PubMed Abstract | CrossRef Full Text | Google Scholar

Dobin, A., Davis, C. A., Schlesinger, F., Drenkow, J., Zaleski, C., Jha, S., et al. (2013). Star: Ultrafast universal RNA-seq aligner. Bioinformatics 29, 15–21. doi:10.1093/bioinformatics/bts635

PubMed Abstract | CrossRef Full Text | Google Scholar

Dowen, J. M., Fan, Z. P., Hnisz, D., Ren, G., Abraham, B. J., Zhang, L. N., et al. (2014). Control of cell identity genes occurs in insulated neighborhoods in mammalian chromosomes. Cell. 159, 374–387. doi:10.1016/j.cell.2014.09.030

PubMed Abstract | CrossRef Full Text | Google Scholar

Edgar, R., Domrachev, M., and Lash, A. E. (2002). Gene expression omnibus: NCBI gene expression and hybridization array data repository. Nucleic Acids Res. 30, 207–210. doi:10.1093/nar/30.1.207

PubMed Abstract | CrossRef Full Text | Google Scholar

Elgin, S. C. R., and Reuter, G. (2013). Position-effect variegation, Heterochromatin formation, and gene silencing in Drosophila. Cold Spring Harb. Perspect. Biol. 5, a017780. doi:10.1101/cshperspect.a017780

PubMed Abstract | CrossRef Full Text | Google Scholar

Ernst, J., and Kellis, M. (2012). ChromHMM: Automating chromatin-state discovery and characterization. Nat. Methods 9, 215–216. doi:10.1038/nmeth.1906

PubMed Abstract | CrossRef Full Text | Google Scholar

Fanucchi, S., Shibayama, Y., Burd, S., Weinberg, M. S., and Mhlanga, M. M. (2013). Chromosomal contact permits transcription between coregulated genes. Cell. 155, 606–620. doi:10.1016/j.cell.2013.09.051

PubMed Abstract | CrossRef Full Text | Google Scholar

Filippova, D., Patro, R., Duggal, G., and Kingsford, C. (2014). Identification of alternative topological domains in chromatin. Algorithms Mol. Biol. 9, 14–11. doi:10.1186/1748-7188-9-14

PubMed Abstract | CrossRef Full Text | Google Scholar

Fournier, M., Bourriquen, G., Lamaze, F. C., Côté, M. C., Fournier, É., Joly-Beauparlant, C., et al. (2016). FOXA and master transcription factors recruit Mediator and Cohesin to the core transcriptional regulatory circuitry of cancer cells. Sci. Rep. 6, 34962. doi:10.1038/srep34962

PubMed Abstract | CrossRef Full Text | Google Scholar

Frank, F., Liu, X., and Ortlund, E. A. (2021). Glucocorticoid receptor condensates link DNA-dependent receptor dimerization and transcriptional transactivation. Proc. Natl. Acad. Sci. U. S. A. 118, e2024685118. doi:10.1073/pnas.2024685118

PubMed Abstract | CrossRef Full Text | Google Scholar

Gu, Z., and Hübschmann, D. (2023). rGREAT: an R/bioconductor package for functional enrichment on genomic regions. Bioinformatics 39, btac745–19. doi:10.1093/bioinformatics/btac745

PubMed Abstract | CrossRef Full Text | Google Scholar

Harmston, N., Ing-Simmons, E., Perry, M., Barešić, A., and Lenhard, B. (2015). GenomicInteractions: An R/Bioconductor package for manipulating and investigating chromatin interaction data. BMC Genomics 16, 963–969. doi:10.1186/s12864-015-2140-x

PubMed Abstract | CrossRef Full Text | Google Scholar

Hnisz, D., Shrinivas, K., Young, R. A., Chakraborty, A. K., and Sharp, P. A. (2017). A phase separation model for transcriptional control. Cell. 169, 13–23. doi:10.1016/j.cell.2017.02.007

PubMed Abstract | CrossRef Full Text | Google Scholar

Hsieh, T. S., Cattoglio, C., Slobodyanyuk, E., Hansen, A. S., Rando, O. J., Tjian, R., et al. (2020). Resolving the 3D landscape of transcription-linked mammalian chromatin folding. Mol. Cell. 78, 539–553.e8. doi:10.1016/j.molcel.2020.03.002

PubMed Abstract | CrossRef Full Text | Google Scholar

Huang, H., Chen, S. T., Titus, K. R., Emerson, D. J., Bassett, D. S., and Phillips-Cremins, J. E. (2019). A subset of topologically associating domains fold into mesoscale core-periphery networks. Sci. Rep. 9, 9526. doi:10.1038/s41598-019-45457-9

PubMed Abstract | CrossRef Full Text | Google Scholar

Kagey, M. H., Newman, J. J., Bilodeau, S., Zhan, Y., Orlando, D. A., van Berkum, N. L., et al. (2010). Mediator and cohesin connect gene expression and chromatin architecture. Nature 467, 430–435. doi:10.1038/nature09380

PubMed Abstract | CrossRef Full Text | Google Scholar

Kaul, A., Bhattacharyya, S., and Ay, F. (2020). Identifying statistically significant chromatin contacts from Hi-C data with FitHiC2. Nat. Protoc. 15, 991–1012. doi:10.1038/s41596-019-0273-0

PubMed Abstract | CrossRef Full Text | Google Scholar

Kent, W. J., Sugnet, C. W., Furey, T. S., Roskin, K. M., Pringle, T. H., Zahler, A. M., et al. (2002). The human genome browser at UCSC. Genome Res. 12, 996–1006. doi:10.1101/gr.229102

PubMed Abstract | CrossRef Full Text | Google Scholar

Kim, S., and Shendure, J. (2019). Mechanisms of interplay between transcription factors and the 3D genome. Mol. Cell. 76, 306–319. doi:10.1016/j.molcel.2019.08.010

PubMed Abstract | CrossRef Full Text | Google Scholar

Krietenstein, N., Abraham, S., Venev, S. V., Mirny, L. A., Dekker, J., Rando, O. J., et al. (2020). Ultrastructural details of mammalian chromosome architecture. Mol. Cell. 78, 554–565.e7. doi:10.1016/j.molcel.2020.03.003

PubMed Abstract | CrossRef Full Text | Google Scholar

Lambert, S. A., Jolma, A., Campitelli, L. F., Das, P. K., Yin, Y., Albu, M., et al. (2018). The human transcription factors. Cell. 172, 650–665. doi:10.1016/j.cell.2018.01.029

PubMed Abstract | CrossRef Full Text | Google Scholar

Langmead, B., and Salzberg, S. L. (2012). Fast gapped-read alignment with Bowtie 2. Nat. Methods 9, 357–359. doi:10.1038/nmeth.1923

PubMed Abstract | CrossRef Full Text | Google Scholar

Le Dily, F., and Beato, M. (2018). Signaling by steroid hormones in the 3D nuclear space. Int. J. Mol. Sci. 19, 306. doi:10.3390/ijms19020306

PubMed Abstract | CrossRef Full Text | Google Scholar

Le Dily, F., Baù, D., Pohl, A., Vicent, G. P., Serra, F., Soronellas, D., et al. (2014). Distinct structural transitions of chromatin topological domains correlate with coordinated hormone-induced gene regulation. Genes. Dev. 28, 2151–2162. doi:10.1101/gad.241422.114

PubMed Abstract | CrossRef Full Text | Google Scholar

Lee, T. I., and Young, R. A. (2013). Transcriptional regulation and its misregulation in disease. Cell. 152, 1237–1251. doi:10.1016/j.cell.2013.02.014

PubMed Abstract | CrossRef Full Text | Google Scholar

Lee, K., Hsiung, C. C., Huang, P., Raj, A., and Blobel, G. (2015). Dynamic enhancer-gene body contacts during transcription elongation. Genes. Dev. 29, 1992–1997. doi:10.1101/gad.255265.114

PubMed Abstract | CrossRef Full Text | Google Scholar

Li, H., and Durbin, R. (2009). Fast and accurate short read alignment with Burrows – wheeler transform. Bioinformatics 25, 1754–1760. doi:10.1093/bioinformatics/btp324

PubMed Abstract | CrossRef Full Text | Google Scholar

Liao, Y., Smyth, G. K., and Shi, W. (2014). FeatureCounts: An efficient general purpose program for assigning sequence reads to genomic features. Bioinformatics 30, 923–930. doi:10.1093/bioinformatics/btt656

PubMed Abstract | CrossRef Full Text | Google Scholar

Lieberman-Aiden, E., van Berkum, N. L., Williams, L., Imakaev, M., Ragoczy, T., Telling, A., et al. (2009). Comprehensive mapping of long-range interactions reveals folding principles of the human genome. Science 326, 289–293. doi:10.1126/science.1181369

PubMed Abstract | CrossRef Full Text | Google Scholar

Liu, Z., and Tjian, R. (2018). Visualizing transcription factor dynamics in living cells. J. Cell. Biol. 217, 1181–1191. doi:10.1083/jcb.201710038

PubMed Abstract | CrossRef Full Text | Google Scholar

López-Maury, L., Marguerat, S., and Bähler, J. (2008). Tuning gene expression to changing environments: From rapid responses to evolutionary adaptation. Nat. Rev. Genet. 9, 583–593. doi:10.1038/nrg2398

PubMed Abstract | CrossRef Full Text | Google Scholar

Love, M. I., Huber, W., and Anders, S. (2014). Moderated estimation of fold change and dispersion for RNA-seq data with DESeq2. Genome Biol. 15, 550. doi:10.1186/s13059-014-0550-8

PubMed Abstract | CrossRef Full Text | Google Scholar

McDowell, I. C., Barrera, A., D’Ippolito, A. M., Vockley, C. M., Hong, L. K., Leichter, S. M., et al. (2018). Glucocorticoid receptor recruits to enhancers and drives activation by motif-directed binding. Genome Res. 28, 1272–1284. doi:10.1101/gr.233346.117

PubMed Abstract | CrossRef Full Text | Google Scholar

McLean, C. Y., Bristor, D., Hiller, M., Clarke, S. L., Schaar, B. T., Lowe, C. B., et al. (2010). GREAT improves functional interpretation of cis-regulatory regions. Nat. Biotechnol. 28, 495–501. doi:10.1038/nbt.1630

PubMed Abstract | CrossRef Full Text | Google Scholar

Misteli, T. (2020). The self-organizing genome: Principles of genome architecture and function. Cell. 183, 28–45. doi:10.1016/j.cell.2020.09.014

PubMed Abstract | CrossRef Full Text | Google Scholar

Nora, E. P., Lajoie, B. R., Schulz, E. G., Giorgetti, L., Okamoto, I., Servant, N., et al. (2012). Spatial partitioning of the regulatory landscape of the X-inactivation centre. Nature 485, 381–385. doi:10.1038/nature11049

PubMed Abstract | CrossRef Full Text | Google Scholar

Osborne, C. S., Chakalova, L., Brown, K. E., Carter, D., Horton, A., Debrand, E., et al. (2004). Active genes dynamically colocalize to shared sites of ongoing transcription. Nat. Genet. 36, 1065–1071. doi:10.1038/ng1423

PubMed Abstract | CrossRef Full Text | Google Scholar

Portuguez, A. S., Grbesa, I., Tal, M., Deitch, R., Raz, D., Kliker, L., et al. (2022). Ep300 sequestration to functionally distinct glucocorticoid receptor binding loci underlie rapid gene activation and repression. Nucleic Acids Res. 50, 6702–6714. doi:10.1093/nar/gkac488

PubMed Abstract | CrossRef Full Text | Google Scholar

Pradhan, M. A., Blackford, J. A., Devaiah, B. N., Thompson, P. S., Chow, C. C., Singer, D. S., et al. (2016). Kinetically defined mechanisms and positions of action of two new modulators of glucocorticoid receptor-regulated gene induction. J. Biol. Chem. 291, 342–354. doi:10.1074/jbc.M115.683722

PubMed Abstract | CrossRef Full Text | Google Scholar

Ramírez, F., Ryan, D. P., Grüning, B., Bhardwaj, V., Kilpert, F., Richter, A. S., et al. (2016). deepTools2: a next generation web server for deep-sequencing data analysis. Nucleic Acids Res. 44, 160–165. doi:10.1093/nar/gkw257

CrossRef Full Text | Google Scholar

Rao, S. S. P., Huntley, M. H., Durand, N. C., Stamenova, E. K., Bochkov, I. D., Robinson, J. T., et al. (2014). A 3D map of the human genome at kilobase resolution reveals principles of chromatin looping. Cell. 159, 1665–1680. doi:10.1016/j.cell.2014.11.021

PubMed Abstract | CrossRef Full Text | Google Scholar

Reddy, T. E., Pauli, F., Sprouse, R. O., Neff, N. F., Newberry, K. M., Garabedian, M. J., et al. (2009). Genomic determination of the glucocorticoid response reveals unexpected mechanisms of gene regulation. Genome Res. 19, 2163–2171. doi:10.1101/gr.097022.109

PubMed Abstract | CrossRef Full Text | Google Scholar

Ross-Innes, C. S., Stark, R., Teschendorff, A. E., Holmes, K. a., Ali, H. R., Dunning, M. J., et al. (2012). Differential oestrogen receptor binding is associated with clinical outcome in breast cancer. Nature 481, 389–393. doi:10.1038/nature10730

PubMed Abstract | CrossRef Full Text | Google Scholar

Rowley, M. J., and Corces, V. G. (2018). Organizational principles of 3D genome architecture. Nat. Rev. Genet. 19, 789–800. doi:10.1038/s41576-018-0060-8

PubMed Abstract | CrossRef Full Text | Google Scholar

Rowley, M. J., Nichols, M. H., Lyu, X., Ando-Kuri, M., Rivera, I. S. M., Hermetz, K., et al. (2017). Evolutionarily conserved principles predict 3D chromatin organization. Mol. Cell. 67, 837–852.e7. doi:10.1016/j.molcel.2017.07.022

PubMed Abstract | CrossRef Full Text | Google Scholar

Sabari, B. R., Agnese, A. D., Boija, A., Klein, I. A., Coffey, E. L., Abraham, B. J., et al. (2018). Coactivator condensation at super-enhancers links phase separation and gene control. Science 361, 1–11. doi:10.1126/science.aar3958

CrossRef Full Text | Google Scholar

Schmidt, S. F., Larsen, B. D., Loft, A., and Mandrup, S. (2016). Cofactor squelching: Artifact or fact? BioEssays 38, 618–626. doi:10.1002/bies.201600034

PubMed Abstract | CrossRef Full Text | Google Scholar

Schoenfelder, S., and Fraser, P. (2019). Long-range enhancer–promoter contacts in gene expression control. Nat. Rev. Genet. 20, 437–455. doi:10.1038/s41576-019-0128-0

PubMed Abstract | CrossRef Full Text | Google Scholar

Servant, N., Varoquaux, N., Lajoie, B. R., Viara, E., Chen, C. J., Vert, J. P., et al. (2015). HiC-Pro: an optimized and flexible pipeline for Hi-C data processing. Genome Biol. 16, 259. doi:10.1186/s13059-015-0831-x

PubMed Abstract | CrossRef Full Text | Google Scholar

Silveira, M. A., and Bilodeau, S. (2018). Defining the transcriptional ecosystem. Mol. Cell. 72, 920–924. doi:10.1016/j.molcel.2018.11.022

PubMed Abstract | CrossRef Full Text | Google Scholar

Soler-Oliva, M. E., Guerrero-Martínez, J. A., Bachetti, V., and Reyes, J. C. (2017). Analysis of the relationship between coexpression domains and chromatin 3D organization. PLoS Comput. Biol. 13, e1005708–e1005725. doi:10.1371/journal.pcbi.1005708

PubMed Abstract | CrossRef Full Text | Google Scholar

Spitz, F., and Furlong, E. E. M. (2012). Transcription factors: From enhancer binding to developmental control. Nat. Rev. Genet. 13, 613–626. doi:10.1038/nrg3207

PubMed Abstract | CrossRef Full Text | Google Scholar

Stortz, M., Pecci, A., Presman, D. M., and Levi, V. (2020). Unraveling the molecular interactions involved in phase separation of glucocorticoid receptor. BMC Biol. 18, 59–20. doi:10.1186/s12915-020-00788-2

PubMed Abstract | CrossRef Full Text | Google Scholar

Thormann, V., Rothkegel, M. C., Sch, R., Glaser, L. V., Djuric, P., Li, N., et al. (2018). Genomic dissection of enhancers uncovers principles of combinatorial regulation and cell type-specific wiring of enhancer – promoter contacts. Nucleic Acids Res. 46, 2868–2882. doi:10.1093/nar/gky051

PubMed Abstract | CrossRef Full Text | Google Scholar

van Mierlo, G., Pushkarev, O., Kribelbauer, J. F., and Deplancke, B. (2023). Chromatin modules and their implication in genomic organization and gene regulation. Trends Genet. 39, 140–153. doi:10.1016/j.tig.2022.11.003

PubMed Abstract | CrossRef Full Text | Google Scholar

Vockley, C. M., D’Ippolito, A. M., McDowell, I. C., Majoros, W. H., Safi, A., Song, L., et al. (2016). Direct GR binding sites potentiate clusters of TF binding across the human genome. Cell. 166, 1269–1281.e19. doi:10.1016/j.cell.2016.07.049

PubMed Abstract | CrossRef Full Text | Google Scholar

Voss, T. C., and Hager, G. L. (2015). Dynamic regulation of transcriptional states by chromatin and transcription factors. Nat. Rev. Genet. 15, 69–81. doi:10.1038/nrg3623

CrossRef Full Text | Google Scholar

Wang, S., Su, J., Beliveau, B. J., Bintu, B., Moffitt, J. R., Wu, C., et al. (2016). Spatial organization of chromatin domains and compartments in single chromosomes. Science 353, 598–602. doi:10.1126/science.aaf8084

PubMed Abstract | CrossRef Full Text | Google Scholar

Zamudio, A. V., Dall’Agnese, A., Henninger, J. E., Manteiga, J. C., Afeyan, L. K., Hannett, N. M., et al. (2019). Mediator condensates localize signaling factors to key cell identity genes. Mol. Cell. 76, 753–766.e6. doi:10.1016/j.molcel.2019.08.016

PubMed Abstract | CrossRef Full Text | Google Scholar

Zhang, Y., Liu, T., Meyer, C. A., Eeckhoute, J., Johnson, D. S., Bernstein, B. E., et al. (2008). Model-based analysis of ChIP-seq (MACS). Genome Biol. 9, R137. doi:10.1186/gb-2008-9-9-r137

PubMed Abstract | CrossRef Full Text | Google Scholar