Active transcription and epigenetic reactions synergistically regulate meso-scale genomic organization – Nature.com

Numerical simulations capture experimentally observed features of chromatin organization

We have developed a mathematical model to capture dynamic chromatin organization in the nucleus, in terms of its compaction into the heterochromatic phase or decompaction into the euchromatic phase (Fig.1a). We treat the meso-scale genomic organization as a dynamic, far-from-equilibrium process, governed by the energetics of phase-separation in conjunction with the kinetics of epigenetic reactions and the formation of chromatin loops aided by supercoiled DNA extrusion through cohesin due to RNAPII-mediated transcription. The model ingredients are depicted schematically in Fig.1a. We begin by defining the energetics of the chromatin distribution in terms of the entropic-enthalpic balance of chromatin-chromatin interactions, the chromatin-lamina interactions as well as the penalty on the formation of phase boundaries via Eq. (6) (refer Methods, and Supplementary Section S1.2 in the SI). The gradients in the free-energy landscape, defined as the chemical potential (refer Supplementary Eq. (S3)), drive the dynamic evolution of chromatin towards the two energy wells corresponding to the euchromatin and heterochromatin phases via Eq. (7a, b) (refer Methods, Supplementary Section (S1.4) in the SI). Interconversion of the two phases of chromatin can occur via (a) epigenetic regulation of histone acetylation and methylation (Fig.1b), and (b) supercoiling-driven extrusion of chromatin loops from heterochromatin into euchromatin along the phase boundaries (Eq. (7b)) as shown in Fig.1c.

a Schematic of a portion of a nucleus showing the multiple mechanisms involved in chromatin organization such as chromatin-chromatin interactions, the chromatin-lamina interactions and epigenetic regulation. Additionally, extrusion of chromatin loops due to DNA supercoiling which is increased by transcriptional activity also plays a role in meso-scale genomic organization. While this may occur within either chromatin phases (red circle), we further explore the role of chromatin loop extrusion at the heterochromatin-euchromatin interface (black circle). b The model captures the chromatin-chromatin interaction energetics via a double well free energy description as shown in the contour plot. The two wells correspond to the heterochromatin (red circle) and euchromatin phases (blue circle). Any initial configuration (light blue circle) spontaneously decomposes into these wells at steady state. The dynamics of this transition are governed by diffusion and reaction kinetics comprising of epigenetic regulation and kinetics of supercoiling-driven chromatin extrusion (red box inset). c Loading of cohesin assisted by NIPBL/MAU2 initiates the formation of chromatin loops. Cohesin can also be dynamically unloaded via unloading factors viz. WAPL/PDS5. Active processes such as RNAPII mediated transcription further drive the extrusion of trapped DNA, supercoiling it into chromatin loops.

The process of phase separation is initiated by adding a random perturbation to the initially uniform chromatin configuration (as shown in Fig.2a, left panel) which captures the intrinsic intranuclear heterogeneities. As the simulation progresses heterochromatin domains (in red, center panel of Fig.2a) spontaneously nucleate and grow. The evolution ultimately stabilizes resulting in a steady state (right panel of Fig.2a) with a quasi-periodic distribution of stable domains of heterochromatin rich phase (({phi }_{h}={phi }_{h}^{max })) in red and euchromatin rich phase (left({phi }_{h}=0right)) in blue. Each of these domains are nearly circular (see Supplementary Section S2 of SI for a discussion on non-circular lamellar domains) with characteristic sizes. Concomitantly, heterochromatin domains localized to the nuclear lamina (called LADs) of comparable sizes appear in our simulations (Fig.2a).

a Visualization of the chromatin organization obtained from the simulations. The initial chromatin organization is a homogenous distribution with a small perturbation added, resulting in nucleation of heterochromatin domains (center panel) which grow into heterochromatin domains of characteristic sizes at a steady state. b Super-resolution visualizations of chromatin organization observed in-vivo via STORM imaging of HeLa nuclei (left panel, scale bar 3m, data previously reported in ref. 19, n=19 nuclei) and ChromSTEM imaging of BJ fibroblast nuclei (right panel, scale bar 1m, n=1 nucleus) show that chromatin organization in nucleus is characterized by interspersed heterochromatic domains of comparable sizes. c The smooth boundaries of the chromatin packing domains as seen in ChromSTEM observations are captured by the model. d Numerically predicted trend of sizes of heterochromatin domains as the transcription-mediated chromatin extrusion rate increases. e Schematic diagrams of the step-by step events (events i through vi) involved in the nucleation, growth and stabilization of heterochromatin domains at a steady state. f Plot of theoretically evaluated growth rate of heterochromatin domains with (red) and without (blue) reactions. Reactions give rise to a stable domain radius. In the absence of reactions, no stable heterochromatin domain length scales are observed. g The evaluation of stable radius (blue) and stable LAD thickness (red) as transcription mediated surface reactions are changed. Here, the relative radius is defined as the steady state radius relative to its value when transcription is zero, i.e., relative radius = ({widetilde{R}}_{d}^{{SS}}/{widetilde{R}}_{d}^{{SS}}|{Gamma }_{a}=0). The relative LAD thickness is similarly defined.

The meso-scale distribution of chromatin throughout the nucleus predicted by the mathematical model presents a striking qualitative similarity with the experimentally observed distribution of DNA in the nucleus using ChromSTEM, and STORM as reported previously19 (Fig.2b). Domains of compacted chromatin with a characteristic size are observed via a high histone density distinguished from regions of low histone density (Fig.2b). Lastly, the preferential accumulation of heterochromatin domains along the nuclear periphery seen via STORM imaging (Fig.2b), again with similar size scale, is also in excellent agreement with the experiments.

When defining the free energy density of chromatin organization in the nucleus (see Supplementary Eq. (S1) in SI), we penalized the formation of sharp interfaces via an interface penalty (eta), defined as the energy cost associated with the formation of the interfaces between heterochromatin and euchromatin phases. As we show in the SI (Supplementary Section S1.5), the energy penalty (eta) results in the formation of a smooth rather than a sharp interface between the heterochromatin and the euchromatin phases. Numerical simulations of chromatin organization exhibit such smooth interfaces around chromatin domains, as shown in the zoomed in image in Fig.2c (right panel). The width of the interface (delta) is controlled by the competition between the interfacial and bulk energy contributions (refer Supplementary Section S1.5).

Smooth chromatin phase boundaries are indeed observed in-vivo via Chrom-STEM imaging (Supplementary Section S1.11). We characterized the 3D chromatin density around individual heterochromatin domains in a BJ fibroblast nucleus using Chrom-STEM (Fig.2c, left panel; Supplementary Fig.S5). We estimated the average chromatin density within concentric circles emerging from the center of individual domains to the periphery (Fig.2c, Supplementary Fig.S5). The chromatin density was highest at the core of the domain and dropped slowly from the center of the domain to the periphery. The smooth decrease in radial density indicates that the chromatin domain boundaries are not abrupt (Fig.2c), in agreement with the numerical simulations.

We next investigate how the size scaling of the heterochromatin domains is regulated by the epigenetic reactions acetylation and methylation of histones and supercoiling-driven chromatin extrusion which together can lead to interconversion between heterochromatin and euchromatin. First, we see that in the absence of the epigenetic reactions and chromatin extrusion multiple domains of a characteristic size are not obtained as shown in Supplementary Fig.S10 (detailed discussion in Supplementary Section S5). In this case, although nucleation of multiple heterochromatin domains occurs even without reactions (Supplementary Fig.S10a), all of them merge into a single large cluster driven by Ostwald ripening so as to minimize the interface formation.

The model also predicts that the size of the heterochromatin domains in the interior and periphery can be regulated by the epigenetic reaction rates of acetylation and methylation as shown in Supplementary Fig.S6 (Supplementary Section S2). We see that as methylation increases the size of the interior domains increases too. On the other hand, increase in acetylation results in the formation of smaller heterochromatin domains. The trends followed by the domains towards the interior of the nucleus are replicated by the LADs as well. Lastly, we identify that the size scales of the domains the domain radii in the interior of the nucleus and the LAD thickness along its periphery depend on the level of transcription governed supercoiling-driven chromatin extrusion rate ({widetilde{Gamma }}_{a}) (Fig.2d, Supplementary Fig.S6). We note that, as the transcription (({widetilde{Gamma }}_{a})) is increased, the sizes of the heterochromatin domains decrease, both in the interior as well as at the periphery. At the same time, we also note that as chromatin extrusion rate is increased, the average volume fraction of heterochromatin (left({bar{phi }}_{h}right)) in the nucleus decreases, while that of euchromatin (left({bar{phi }}_{e}right)) increases.

Next, we theoretically predict an explicit dependence of the sizes of interior heterochromatic domains and LADs on epigenetic and transcription reactions and the diffusion kinetics of the epigenetic marks.

Intuitively, in the presence of more repressive methylation the overall heterochromatin content in the nucleus should increase, while in higher histone acetylation conditions the overall euchromatin content will increase. Thus, the epigenetic reactions can independently determine the average volume fractions of each form of chromatin, thereby breaking the detailed balance condition where the free energies of each phase determine their relative abundance in a thermodynamic equilibrium. A mathematical relation between the average volume fraction of each chromatin phase and the epigenetic reaction parameters can be determined by averaging the chromatin evolution equation (Eq. (7b)) at a steady state (i.e. (frac{partial {{{{{{rm{phi }}}}}}}_{{{{{{rm{d}}}}}}}}{partial widetilde{t}}=0)). In the absence of transcription driven chromatin extrusion (i.e. ({widetilde{Gamma }}_{a}=0)), we see that the epigenetic kinetics regulates the average heterochromatin content of the nucleus as, ({bar{phi }}_{h}approx frac{{widetilde{Gamma }}_{{me}},left(1-{bar{phi }}_{n}right)}{{widetilde{Gamma }}_{{me}}+1}) (Supplementary Eq. (S23), refer Supplementary Section S3 for more details).

The presence of transcription-mediated loop extrusion kinetics (i.e., ({widetilde{Gamma }}_{a}, ne , 0) in Eq. (7b)) further augments the deviation from thermodynamic equilibrium (i.e., the breaking of detail balance) via surface reactions that actively extrude DNA at the interface of heterochromatic domains. In the presence of transcription, the average heterochromatin (and euchromatin) content in the nucleus becomes (refer Supplementary Eq. (S22)),

$$begin{array}{c}{bar{phi }}_{h}, approx , frac{{widetilde{Gamma }}_{{me}},left(1-{bar{phi }}_{n}right)}{{widetilde{Gamma }}_{{me}}+1+kappa {widetilde{Gamma }}_{a}},,{bar{phi }}_{e}, approx , frac{left(1+kappa {widetilde{Gamma }}_{a}right),left(1-{bar{phi }}_{n}right)}{{widetilde{Gamma }}_{{me}}+1+kappa {widetilde{Gamma }}_{a}},end{array}$$

(1)

where (kappa) is a function of ({{{{{{rm{phi }}}}}}}_{{{{{{rm{h}}}}}}}^{max }), volume fraction change across the interface (Delta phi), and the length of the interface between the two chromatin phases (refer Supplementary Section S3 for derivation). Since supercoiling-mediated chromatin extrusion converts the tightly packed heterochromatin into low density transcriptionally active euchromatin phase, as extrusion rate ({widetilde{Gamma }}_{a}) increases, the average heterochromatin content decreases.

Thus, the overall mean chromatin composition of the nucleus (left({bar{phi }}_{h},{bar{phi }}_{e}right)) is determined by the reaction kinetics of epigenetic regulation along with transcription. The reaction kinetics alone would drive a homogenous chromatin organization with (left({bar{phi }}_{h},{bar{phi }}_{e}right)). On the (left({phi }_{d},{phi }_{n}right)) phase space we see that the average composition (shown as a light blue circle in Fig.1b) determined by reactions is energetically unfavorable it does not lie in the energy wells and hence must evolve in time.

Next, we show that the average composition of the two chromatin phases, shown in Fig.2e(i), plays a key role in the emergence of the characteristic sizes of the heterochromatin domains. To illustrate this, we first observe that the mean chromatin composition (left({bar{phi }}_{h},{bar{phi }}_{e}right)) lies in neither of the energy wells as shown in Fig.1b (light blue circle) and is thus energetically unfavorable. The need to reduce the total free energy in the nucleus drives the system to phase separate by nucleating heterochromatin domains (Fig.2e(iii)) corresponding to the red energy well labeled heterochromatin in Fig.1b surrounded by euchromatin domains corresponding to the dark blue energy well labeled euchromatin. The events entailing the individual steps in the nucleation and growth of a single droplet of heterochromatin due to phase separation, as shown in Fig.2e, are as follows:

Due to phase separation, the heterochromatin volume fraction immediately outside the droplet is ({phi }_{h}=0) corresponding to the euchromatic energy well. Far away from the droplet, the mean composition (left({bar{phi }}_{h},{bar{phi }}_{e}right)) remains undisturbed. The resulting spatial gradient in the chromatin composition (blue curve in Fig.2e(iv)) sets up a diffusive flux of heterochromatin into the droplet, allowing it to grow.

On the other hand, within the heterochromatin droplet (with ({phi }_{h}={phi }_{h}^{max })) histone acetylation reactions will allow conversion of heterochromatin inside the droplet into euchromatin outside. Active supercoiling-mediated chromatin loop extrusion further adds to the heterochromatin outflux. Together loop extrusion and acetylation oppose the diffusive influx of heterochromatin and thereby reduce the size of the droplet (Fig.2e).

Based on the above observations, the rate at which the nucleated heterochromatin droplet grows can be written in terms of the balance of reaction-diffusion gradient driven influx and acetylation and transcription driven outflux of heterochromatin as (refer Supplementary Section S4, Supplementary Eq. (S25)),

$$4pi {widetilde{R}}_{d}^{2}frac{d{widetilde{R}}_{d}}{dwidetilde{t}}=,underbrace{{4pi {widetilde{R}}_{d}{bar{phi }}_{h}}}_{{{{{{mathrm{inwards}}}}}} , {{{{{mathrm{diffusion}}}}}}} - underbrace{{frac{4}{3}pi {widetilde{R}}_{d}^{3}{phi }_{h}^{max }}}_{begin{array}{c}{{{{{mathrm{Acetylation}}}}}} , {{{{{mathrm{working}}}}}}\ {{{{{mathrm{against}}}}}} , {{{{{mathrm{inwards}}}}}}\ {{{{{mathrm{diffusion}}}}}}end{array}}-underbrace{{4pi {widetilde{R}}_{d}^{2}frac{delta }{2}{widetilde{Gamma }}_{a}{phi }_{h}^{max }}}_{begin{array}{c}{{{{{mathrm{Chromatin}}}}}} , {{{{{mathrm{extrusion}}}}}} , {{{{{mathrm{working}}}}}}\ {{{{{mathrm{against}}}}}} , {{{{{mathrm{inwards}}}}}}\ {{{{{mathrm{diffusion}}}}}}end{array}}$$

(2)

where (delta) is the rescaled width of the interface, which is in turn related to the length scale obtained via the competition between the interfacial energy and chromatin-chromatin interaction (refer Supplementary Section S1.5). The resulting evolution of the droplet growth rate (left(d{widetilde{R}}_{d}/dwidetilde{t}right)) as the radius of the droplet increases is shown in Fig.2e. Notice the two fixed points (Fig.2f, labeled critical and stable radius) where (d{widetilde{R}}_{d}/dwidetilde{t}=0). Beyond the critical radius the domains grow in size.

The second fixed point (stable radius) corresponds to the rescaled steady state (i.e., (d{widetilde{R}}_{d}/dwidetilde{t}=0)) heterochromatin domain size as determined by the active epigenetic and the transcriptional regulation in tandem with passive diffusion, and can be written as (derivation shown in Supplementary Section S4, Supplementary Eq. (S27)),

$$begin{array}{c}{widetilde{R}}_{d}^{{ss}}=-frac{3{widetilde{Gamma }}_{a}delta }{4}+sqrt{{left(frac{3{widetilde{Gamma }}_{a}delta }{4}right)}^{2}+frac{3}{{phi }_{h}^{max }}frac{{widetilde{Gamma }}_{{me}}left(1-{bar{phi }}_{n}right)}{1+{widetilde{Gamma }}_{{me}}+kappa {widetilde{Gamma }}_{a},}}.end{array}$$

(3)

From Eq. (3), we observe that the steady state droplet radius (left({widetilde{R}}_{d}^{{ss}}right)) depends on both diffusion and reaction kinetics. With increase in methylation, ({widetilde{R}}_{d}^{{ss}}) increases implying bigger heterochromatin domains. On the other hand, with increase in either the acetylation or transcription-mediated loop extrusion the steady state radius decreases. The quantitative dependence of the steady state radius on transcriptional kinetics is shown in Fig.2g (blue solid line). Note that the steady state radius shown in Fig.2g is normalized relative to the steady state radius with no transcription. Thus, our theory predicts an increase in the sizes of compacted chromatin domains in the interior of the nucleus upon inhibition of transcription.

The size dependence of chromatin domains along the nuclear periphery can be similarly determined by the balance of reaction, transcription, and diffusion kinetics for the LADs. The affinity of chromatin to the nuclear periphery due to the chromatin-lamina interactions in Eq. (6) induces a preferential nucleation of LADs. A schematic representation of heterochromatin compaction along the nuclear periphery resulting in LAD growth is shown in Fig.2e. As with the interior heterochromatin droplet, phase-separation drives the heterochromatin compaction (left({phi }_{h}={phi }_{h}^{max }right)) within the LADs, while the chromatin immediately outside corresponds to the euchromatin energy minimal well (left({phi }_{h}=0right)). Far away from the peripheral LAD nucleation sites, the chromatin composition remains undisturbed at the average composition of (left({bar{phi }}_{h},{bar{phi }}_{e}right)). The variation of chromatin composition with distance from nuclear periphery is shown in Fig.2e (blue line). Like in the case of the interior heterochromatin droplets, the heterochromatin composition gradient driven diffusive influx is balanced by the epigenetic and transcriptional regulated heterochromatin outflux, which determines the rescaled steady-state thickness of the LADs (refer to the Supplementary Section S7, Supplementary Eq. (S34)),

$${widetilde{x}}_{t}^{{ss}}=frac{{widetilde{Gamma }}_{{me}}left(1-{bar{phi }}_{n}right)}{{phi }_{h}^{max }left(1+{widetilde{Gamma }}_{{me}}+kappa {widetilde{Gamma }}_{a}right),}-frac{delta {widetilde{Gamma }}_{a}}{2}$$

(4)

As with the interior domains, we observe that the LADs become thicker with increase in methylation, while they become thinner with increasing acetylation or chromatin extrusion rates. A quantitative dependence of steady state LAD thickness on transcription rate based on Eq. (4) is plotted in Fig.2g (red dashed line). Our theory predicts an increase in the sizes of LADs along the nuclear periphery upon inhibition of transcription. While the theoretical analysis helps develop a fundamental biophysical understanding of the role of energetics and kinetics in chromatin phase separation, a nucleus-wide chromatin organization and its dynamic evolution can only be obtained numerically.

Next, we use the in-silico model to make testable quantitative predictions of the meso-scale chromatin organization in the nucleus. We also report the in-vivo nuclear chromatin reorganization upon transcription inhibition using complimentary STORM19 and ChromSTEM on nuclei from multiple cell lines. The choice of the parameters for rates of acetylation ({widetilde{Gamma }}_{{ac}}), methylation ({widetilde{Gamma }}_{{me}}), and the strength of chromatin-lamina interactions ({widetilde{V}}_{L}), were held constant for all the following simulations, and the choice of the level of spatial noise is discussed in the Supplementary Section S8. We calibrate the active chromatin supercoiling-driven loop extrusion rate ({Gamma }_{a}) to obtain an in-silico change in the interior domain sizes quantitatively comparable to that observed upon transcriptional inhibition. The calibrated model is then used to predict the change in LAD thickness due to inhibition of transcription, which upon comparison with experimental images serves to validate the model. A schematic for the workflow utilized to calibrate and cross-validate the model predictions in the interior and along periphery of the nucleus is shown in Supplementary Fig. (S14) (Supplementary Section S8).

ChromSTEM was used to obtain super-resolution images in terms of statistical descriptions of chromatin packing domains for BJ fibroblasts. ChromSTEM allows the quantification of 3D chromatin conformation with high resolution22. ChromSTEM mass density tomograms were collected for BJ fibroblasts treated with Actinomycin D (ActD) (Fig.3a, center) and compared to DMSO treated mock controls (Fig.3a, left) to evaluate the average size and density of chromatin packing domains. We have previously demonstrated that chromatin forms spatially welldefined higherorder packing domains and that, within these domains, chromatin exhibits a polymeric power-law scaling behavior with radially decreasing mass density moving outwards from the center of the domain23. As the ChromSTEM intensity in the reconstructed tomogram is proportional to the chromatin mass density, we estimated the size of the domains based on where the chromatin mass scaling and the radial chromatin density deviate from their predicted behavior (discussed in Supplementary Section S1.11). Based on the statistical analysis of individual packing domains, in a single tomograph shown in Fig.3a, we observed 71 domains in DMSO and 48 domains in the ActD-treated nucleus. Of the identified domains, the average domain radius ((pm) S.E) of BJ cells treated with DMSO and ActD was estimated to be 103.5 (pm) 4.73nm and 129.7 (pm) 6.78nm, respectively (Fig.3a, right panel), representing a 20.2% increase in size. Overall, fewer domain centers, and larger chromatin packing domains were experimentally observed upon ActD treatment compared to the control.

a ChromSTEM tomogram reconstructions for DMSO (left panel) and ActD treated (center panel) BJ fibroblasts. The domains radii for BJ cells treated ActD (right panel, n=48 domains) show 1.25 times (unpaired two tail t-test, p=0.002) increase compared to control (n=71 domains). b Representative live-cell PWS images (1-hour ActD treatment). Scale bars=5m. Box plots compare the domain sizes between DMSO control and ActD treated cells. Sample size HCT116: n=63 nuclei (control), 65 (ActD), p=0.05; A549: n=102 (control), 84 (ActD), p=1e7; U2OS: n=116 (control), 75 (ActD), p=1e12; n=103 (control), 150 (ActD), p=0.04. c Heatmap density of DNA super-resolution images in DMSO control (left panel, n=19 nuclei) and ActD (right panel, n=20 nuclei) treated HeLa nuclei. All scale bars 3m. d Loss of chromatin loop extrusion due to absence of RNAPII results in increased heterochromatin domain size (in red, nucleosomes not shown for clarity). e Numerical prediction of chromatin organization in DMSO control and ActD treated nucleus. f Zoomed in views of DMSO and ActD treated nuclei localized to the nucleus interior (top panels) and the periphery (bottom panels). Red and blue boxes shown in c are zoomed into. All scale bars 1m. g Left: Simulations predict domains in ActD nuclei are on average 1.63 times larger than in DMSO nuclei (n=127 (DMSO), 77 (ActD) unpaired two tail t-test, p=0) while LADs are 1.37 times thicker (n=38 (DMSO), 15 (ActD); unpaired two tail t-test, p=0). Right: Domain radii observed experimentally in ActD treated nuclei (n=3584 loci, 20 nuclei) are 1.61 times (unpaired two tail t-test, p=0) larger than in DMSO nuclei (n=5830 loci,19 nuclei), while LADs are 1.3 times thicker (n=1082 loci (DMSO), 1015 loci (ActD), unpaired two tail t-test, p=0.0006). All boxplots show the mean (cross), median (horizontal line), upper and bottom quartiles (box outlines) and the maximum and minimum non-outlier data points (whiskers). All source data are provided as a source data file.

In addition to evaluating domain properties using ChromSTEM, we utilized live-cell partial wave spectroscopy (PWS) imaging to observe the change in chromatin organization after transcription inhibition in various cell lines (Fig.3b). The PWS images demonstrate a significant reduction in average chromatin packing scaling upon ActD treatment in live cells across four different cell types. Next, the size of the domains is quantitatively approximated via polymer scaling relationships discussed in Supplementary Section S1.1322,24. The quantification of the domain sizes (boxplots in Fig.3b) shows that, for all cell types studied, packing domains are larger for upon transcription inhibition with ActD treatment in agreement with the ChromSTEM results on BJ fibroblasts.

Additionally, we have previously used STORM imaging to observe the nucleus wide changes in chromatin organization caused by transcription abrogation in HeLa nuclei after ActD treatment19. Heatmaps of chromatin density obtained via Voronoi tessellation-based color-coding of STORM images (see19 for analysis) are shown in Fig.3c. The zoomed in images of heatmaps of the chromatin cluster density (Fig.3f) clearly show the increasing heterochromatin domain sizes when RNAPII activity is inhibited, in agreement with our theoretical and numerical predictions (Fig.2d, e). Importantly, we see that the changes in chromatin organization occur not only in the interior domains of the nucleus but also along its periphery (Fig.3f, g).

Altogether these complementary imaging techniques establish that nucleus wide increase in sizes of compacted chromatin domains occurs upon the loss of transcription in a wide range of cell lines.

The chromatin cluster density maps obtained from STORM imaging were further analyzed to quantify the sizes of heterochromatin domains after DMSO and ActD treatment. A density-based threshold was used to isolate the high-density heterochromatin regions, which were then clustered via a density based spatial clustering algorithm (see Supplementary Section S1.8) and further sub-classified into LADs and interior domains depending on the distance from nuclear periphery (Supplementary Section S1.9). The quantitatively extracted distribution of interior heterochromatin domain radii for DMSO and ActD treated nuclei shows that their mean radius after transcription inhibition was nearly 1.61 times that in DMSO controls (Fig.3g).

Indeed, our model (Eq 3-4, Fig.2d, g) predicts that loss of transcription results in increased heterochromatin domain size. This is because under control conditions, extrusion of heterochromatin phase into euchromatin occurs. We assume, based on previous experimental findings19, that the presence of RNAPII activity drives the supercoiling of the DNA loop, thereby extruding it from the heterochromatin phase into the euchromatin phase at the phase boundaries (Fig.3c, left panel). However, when RNAPII is inhibited with ActD treatment (Fig.3c, right panel), the absence of this driving force for supercoiling-mediated loop extrusion keeps more DNA in the heterochromatin phase thereby increasing the domain sizes. The in-silico chromatin distribution predicted under control (left panel) and transcription inhibited (({Gamma }_{a}=0), right panel) conditions is shown in Fig.3e. The phase separated heterochromatin domains (left({phi }_{h}={phi }_{h}^{max }right)) are shown in red in a loosely compacted euchromatin background (blue, ({phi }_{h}=0)). We quantify the change in the sizes of the heterochromatin domains predicted by the model as the active extrusion rate ({Gamma }_{a}) is parametrically varied. The value of ({Gamma }_{a}) under control conditions is chosen (Supplementary TableS2) such that the change in the interior domain sizes with respect to transcription inhibition (with ({Gamma }_{a}=0)) is quantitatively the same as observed experimentally.

Next, we quantitatively validate the choice of ({Gamma }_{a}) under control conditions by comparing the predicted change in LAD thickness against that quantified from the STORM images. Our theoretical predictions (Eq. (4)) show that the reduction in transcription increases the thickness of the LADs reflecting the behavior predicted in the interior of the nucleus (Fig.2d, g). Our simulations of chromatin distribution in the nucleus (Fig.3e) show that inhibition of transcription (({Gamma }_{a}=0)) results in thicker LADs. Of note, the chromatin-lamina interaction strength (left({V}_{L}right)) stays unchanged between the two simulations. Yet, we see a higher association of chromatin with the periphery. Upon quantitative comparison (Fig.3g, left panel) we see that the LADs grow approximately 1.37 times thicker upon loss of transcription.

To validate this prediction, we compare the predicted change in LAD thickness with that quantified from in-vivo STORM imaging. (Fig.3g, refer to Supplementary Sections S1.8 and S1.9 for procedure). The quantified comparison of LAD thickness between DMSO and ActD nuclei (Fig.3g) shows nearly 1.3 times increase upon ActD treatment, in close quantitative agreement with the model prediction. Overall, with both model predictions and cellular observations, our results suggest that impairment of transcription plays a significant role in determining the size scaling of the interior heterochromatin domains and LADs.

We next enquire how, in addition to altering the size of the compacted domains, abrogation of transcription changes the extent of DNA packing. For this we analyzed the chromatin distribution in HeLa nuclei under DMSO and ActD treatments from STORM images previously generated19. Under control conditions the distribution of DNA is qualitatively more homogenous while ActD treated nuclei exhibit more isolated distinct domains of compacted chromatin surrounded by region of very low chromatin density (Fig.4a). For quantification, we plot the chromatin intensity along a horizontal line chosen to run across two heterochromatin domains with euchromatin between them (see zoomed images in Fig.4b, blue and red horizontal line). The chromatin intensity, plotted in Fig.4c (in blue) shows that even in the euchromatin region, the DNA presence is substantial. On the other hand, chromatin intensity across a horizontal line chosen across a heterochromatin domain in ActD nucleus (Fig.4b, c; in red) shows a much steeper gradient outside the domain.

a Grayscale heatmap density rendering of super-resolution images of DNA in control (DMSO, left panel, n=19 nuclei) and actinomycin D (ActD, right panel, n=20 nuclei) treated HeLa nuclei. All scale bars 3m. b Zoomed in views of DMSO and ActD treated nuclei. Boxes shown in a are zoomed into. All scale bars - 1m. c Along the blue (DMSO) and red (ActD) line segments, we plot the chromatin heatmap intensity (corresponding to the total DNA content) for the DMSO-treated control nucleus (in blue) and ActD-treated nucleus (in red). The DMSO-treated nucleus shows a wider distribution of small heterochromatin domains, while the ActD treated nucleus shows a greater compaction with isolated large heterochromatin domains. d Numerical prediction of distribution of total DNA (in grayscale) in a nucleus with (DMSO) and without (ActD) transcription mediated chromatin extrusion. e Distribution of total DNA content along the blue (red) line in d under DMSO (ActD) treatment. The black dashed line shows the level of total DNA predicted in the euchromatin phase of DMSO and ActD treated nuclei.

The increased presence of DNA in the euchromatic phase in presence of transcription as observed experimentally is captured by the simulations. The in-silico distribution of DNA (measured as the sum of volume fractions of the chromatin phases, ({phi }_{e}+{phi }_{h})) in a nuclear region far from LADs is plotted in Fig.4d for control and transcription inhibited in-silico nuclei. We see that the euchromatic phase (outside white circles) is darker when transcription is inhibited, indicating the presence of much lesser DNA than in control euchromatin. A quantification of the total DNA along cut-lines chosen in the control and ActD in-silico nuclei confirm the observations (Fig.4e).

Since the lack of transcription inhibits supercoiling-mediated chromatin loop extrusion from heterochromatin into euchromatin, we see a reduced density of DNA in the euchromatin phase of the nucleus under ActD conditions. Further, due to the lack of chromatin extrusion out of the heterochromatin domains when transcription is inhibited, we also observe that they are larger in size. Thus, transcription, via chromatin loop extrusion, results in removal of DNA from compacted heterochromatin region by converting it into active euchromatin form.

Taken together, our results suggest that transcription not only affects the scaling of the lengths (radius or thickness) of the heterochromatin domains, but also significantly changes the relative amounts of DNA in the euchromatin and heterochromatin phases.

We have established that change in transcription activity affects the global chromatin organization of the nucleus via altered supercoiling mediated loop extrusion. In turn, chromatin loop extrusion is initiated by the loading of cohesin onto DNA via a balance between cohesin loaders such as NIPBL and cohesin unloaders like WAPL (Fig.1c2,12,25,26). If the chromatin loop extrusion is responsible for the global chromatin reorganization, altering the cohesin loading/unloading balance must also result in chromatin reorganization. Thus, next, we study the chromatin arrangement in WAPL-deficient (WAPL) nuclei marked by increased levels of loaded cohesin.

In vivo, WAPL depletion causes an accumulation of large amounts of cohesin on chromatin27. This results in a much more homogenous distribution of DNA, which was previously termed blending due to excessive extrusion of chromatin loops, as shown schematically in Fig.5a19. In our mathematical model, WAPL deficiency is simulated as an increase in the rate of chromatin extrusion (({Gamma }_{a})). Based on the theoretical size scaling of the interior heterochromatin domains and LADs, as seen from Eq. (3) and Fig.2g, our model predicts that increase in ({Gamma }_{a}) would result in a decrease in the radius of the steady state heterochromatin domains (Fig.5b).

a Schematic representation of chromatin loop extrusion. WAPL-depletion results in increased cohesin loading and excessive transcription-driven chromatin loop extrusion. Note that nucleosomes, despite being present, are not shown to improve clarity. b Numerical prediction of distribution of heterochromatin domains in the interior and the LADs along the periphery (all domains in red) in a nucleus without (Cas9) and with (WAPL) cohesin unloading disruption. c Heatmap density of DNA super-resolution images in d control (Cas9, left panel) and WAPL knock-out (WAPL) treated HeLa nuclei. All scale bars - 3m. d Left: Zoomed in views of Cas9 and WAPL treated nuclei focusing on the interior heterochromatin domains. White solid boxes shown in c are zoomed into. All scale bars - 1m. Right: Zoomed in views of Cas9 and WAPL treated nuclei along the nuclear periphery. White dashed boxes shown in c are zoomed into. All scale bars - 1m. e Quantification of heterochromatin domain radius in the interior of Cas9- and WAPL - treated nuclei. (n=2386 loci in 6 nuclei for Cas9-treatment and 2416 loci in 7 nuclei for WAPL treatment). WAPL treated nuclei exhibit a significantly lower ((sim) 0.86 times) mean heterochromatin radius (unpaired two-tailed t-test, p=6e10). Quantification of LAD thickness along the periphery of Cas9- and WAPL - treated nuclei. (n=219 loci in 6 nuclei for Cas9-treatment and 169 loci in 7 nuclei for WAPL treatment). WAPL treated nuclei exhibit a significantly lower ((sim) 0.43 times) mean LAD thickness (unpaired two-tailed t-test, p=1e13). f Boxplot in left panel shows the distribution of domain radii predicted numerically. WAPL nuclei have a mean domain radius 0.8 times that of Cas9-treated nuclei (unpaired two-tailed t-test, p=0). Boxplot in right panels shows the distribution of LAD thicknesses predicted numerically. WAPL nuclei have a mean LAD thickness 0.82 times that of Cas9-treated nuclei. All boxplots show the mean (cross), median (horizontal line), upper and bottom quartiles (box outlines) and the maximum and minimum non-outlier data points (whiskers) of the plotted distribution. All source data are provided as a source data file.

STORM images of HeLa nuclei without (labeled Cas9) and with WAPL-deficiency previously revealed genome-wide changes in the chromatin organization induced by excessive loading of cohesin (Fig.5c, d)19. A visual comparison between representative zoomed-in regions (white boxes in Fig.5c) demonstrates the reduction of heterochromatin domain sizes in the interior of the nuclei in WAPL nuclei (Fig.5d). Using clustering analysis (refer Supplementary Section S1.8 and S1.9), we quantify the altered chromatin domain sizes in control and WAPL HeLa cell nuclei. We observe that WAPL nuclei with increased chromatin blending have heterochromatin domains with a mean radius approximately 15% smaller than control nuclei (Fig.5e).

In-silico, we parametrically vary the active chromatin extrusion rate ({Gamma }_{a}) above the control level (SupplementaryTable S2, determined forcontrol treatment). The value of ({Gamma }_{a}) for WAPL nuclei is chosen (Supplementary TableS2) such that the decrease in the size of interior heterochromatin domains reduces by 15% (Fig.5f) to agree with the experimental observation (Fig.5e).

As discussed previously (Fig.2g), the model predicts that the effects of chromatin extrusion observed in the interior domains of the nucleus are replicated along the nuclear periphery. Simulation of nuclear chromatin organization (Fig.5b) reveals that by changing only the rate of chromatin extrusion ({Gamma }_{a}), keeping all other parameters including chromatin-lamina interaction potential ({V}_{L}) constant, we see a reduction in the association of chromatin with the lamina. Specifically, a 2.5-fold increase in ({Gamma }_{a}) calibrated to occur due to WAPL-deficiency predicts a 51.2% decrease in the average LAD thickness, as shown in Fig.5f.

The predicted change in LAD thickness is consistent with previous experimental observations and was further quantitatively validated by measuring the thickness of LADs in STORM images of control and WAPL nuclei (Fig.5e)19. A reduction in the sizes of domains, as seen in the nucleus interior, can also be observed at the nuclear periphery, as shown in a representative zoomed in region (white dashed boxes in Fig.5c) in Fig.5d. The mean thickness of the LADs at the nuclear periphery is approximately 20% smaller for WAPL nuclei (Fig.5h) as compared to the control-treated nuclei.

Together, these results confirm that the meso-scale spatial chromatin organization is strongly regulated by the chromatin loop formation, and this effect can be modulated not only by the transcription activity, but also by altering the extent of loading or unloading of cohesin rings on the DNA. These results provide further evidence for the link between transcriptional regulation and nucleus-wide chromatin distribution via transcription-driven supercoiling mediated chromatin loop extrusion.

Since we have established, via both quantitative analysis of experimental data and simulations, that extrusion of chromatin loops is governed by both cohesin loading/unloading balance and RNAPII mediated transcription, a question of their tandem role emerges.

To simulate the individual effects of cohesin loading and transcriptional activity, we decompose the overall active chromatin extrusion rate into its distinct constitutive steps. The individual steps involved in the process of supercoiling mediated chromatin loop extrusion from heterochromatin into euchromatin (as discussed previously in Section Introduction) are shown in Fig.6a. As a first step, a balance between the loading of cohesin via NIPBL/MAU225 on chromatin occurring at a rate ({Gamma }_{l}) and its unloading via by WAPL/PDS52,12,26 occurring at a rate ({Gamma }_{{ul}}) results in the association of cohesin rings with chromatin at an overall rate ({Gamma }_{{coh}}={Gamma }_{l}-{Gamma }_{{ul}}). In other words, ({Gamma }_{{coh}}) denotes the overall rate of cohesin loading on DNA. The entrapment of DNA by cohesin is followed by the extrusion of supercoiled loops of chromatin via DNA supercoiling by the RNAPII mediated transcription, at a rate denoted by ({Gamma }_{{tr}}). Thus, as shown in Fig.6a, by assuming a first-order reaction kinetics for both steps, the overall rate of active chromatin extrusion ({Gamma }_{a}) at the interface of heterochromatin and euchromatin is proposed to be multiplicatively decomposed as,

$$begin{array}{c}{Gamma }_{a}={Gamma }_{{tr}}{Gamma }_{{coh}}={Gamma }_{{tr}}left({Gamma }_{l}-{Gamma }_{{ul}}right)end{array}$$

(5)

a Schematic showing the associative sub-steps of chromatin extrusion incorporating cohesin loading v/s unloading balance and active transcriptional work done by RNAPII. The rate of active extrusion of chromatin loops (left({Gamma }_{a}right)) is determined by both sub-steps. Note that nucleosomes, despite being present, are not represented in this schematic to better display the chromatin loops. b Numerical prediction of distribution of heterochromatin domains in the interior and the LADs along the periphery (all domains in red) in a nucleus in control (Cas9-DMSO treatment, top-left panel), transcription inhibited (Cas9-ActD, top right), WAPL knock-out treated (WAPL-DMSO, bottom left) and simultaneous WAPL knock-out along with transcription inhibition treated (WAPL-ActD, bottom right). c Heatmap density rendering of super-resolution images of DNA in control (Cas9-DMSO treatment, left panel), transcription inhibited (Cas9-ActD, center left), WAPL knock-out treated (WAPL-DMSO, center right) and simultaneous WAPL knock-out along with transcription inhibition treated (WAPL-ActD) HeLa nuclei. All scale bars 3m. d Quantification of heterochromatin domain radius in the interior (plain colored boxes) as well as the LAD thickness along the nuclear periphery (hatched boxes) of Cas9-DMSO (3328 loci in 13 nuclei), Cas9-ActD (4042 loci in 11 nuclei), WAPL-DMSO (1548 loci in 10 nulcei) and WAPL-ActD (1926 loci in 11 nuclei) treated nuclei. As previously, ActD treated nuclei exhibited a significantly increased domain size (unpaired two-tailed t-test, p=0) while WAPL treated nuclei exhibit a significantly lower mean heterochromatin radius (unpaired tw-tailed t-test, p=0). However, the differences between Cas9-ActD treated and WAPL-ActD treated nuclei was insignificant (unpaired two-tailed t-test, p (sim) 0.9). All boxplots show the mean (cross), median (horizontal line), upper and bottom quartiles (box outlines) and the maximum and minimum non-outlier data points (whiskers) of the plotted distribution. All source data are provided as a source data file.

In addition to the extrusion of loops via RNAPII mediated DNA supercoiling activity12,13,19,28,29,30, in vitro experiments proposed that cohesin once transiently loaded onto DNA, could independently drive the formation of loops via its ATPase machinery9,11,31,32,33. Cell based experiments demonstrated that in WAPL cells, clusters of cohesin in WAPL cells assemble together into vermicelli-like structures and these structures disappear upon transcription inhibition, but not upon partial loss of cohesin19. These results, taken together, present strong evidence for the important role of transcription in powering cohesin mediated loop extrusion. While the relative role of cohesins motor activity and transcription in loop extrusion inside cells remains to be determined, here we focus on the latter given the previous in vivo experimental findings. We indeed show that a kinetic model captured by Eq. (5) sufficiently explains the effect of extrusion of the specific chromatin loops extending from transcriptionally silenced heterochromatin into genetically active euchromatin on determining the meso-scale chromatin domain sizes.

The chromatin organization is simulated in a nucleus under control and transcription inhibition treatments for nuclei with and without WAPL deficiency. The chromatin organization in a control nucleus (labeled Cas9-DMSO), simulated via parameters listed in Supplementary TableS1 is shown in Fig.6b, top-left panel. The individual inhibition of transcriptional activity without affecting the cohesin loading (Cas9-ActD) results in a chromatin organization with increased heterochromatin domains sizes and LAD thickness, as shown in Fig.6b, top-right panel. On the other hand, the simulation of chromatin distribution in nucleus with depleted cohesin unloading, without disturbing the transcriptional activity, (WAPL-DMSO) is shown in Fig.6b, bottom-left panel. Finally, the chromatin distribution predicted in a WAPL nucleus with inhibited transcription (WAPL-DMSO-treatment) is shown in Fig.6b, bottom-right panel. As shown in Fig. 3e and Fig. 3g, ActD (mathematically, ({Gamma }_{{tr}}=0) in Eq. (5)) results in larger heterochromatin domains and thicker LADs, while WAPL nuclei (increased cohesin loading; mathematically, ({Gamma }_{{ul}}/{Gamma }_{l}) increases in Eq. (5)) show the opposite effect with smaller heterochromatin domains and LADs. For a WAPL nuclei in which transcription is inhibited (WAPL ActD; mathematically, ({Gamma }_{{tr}}=0) and ({Gamma }_{{ul}}/{Gamma }_{l}) increases in Eq. (3)), the model predicts that inhibition of transcription returns the chromatin organization to the control (Cas9-ActD) levels. Transcription inhibition thus blocks the reduction in chromatin domain sizes induced due to WAPL deficiency due to lack of impetus for chromatin supercoiling.

To quantitatively validate the model predictions, we investigate the in-vivo chromatin organization under individual and tandem changes in transcription and cohesin unloading by re-analyzing previously reported super-resolution images shown as heatmap density plots in Fig.6c19. Visual inspection of this data agrees with the model predictions that transcriptional inhibition counteracts the chromatin blending observed in DMSO treated WAPL nuclei, which was also previously reported19. We thus focused on extracting the radius of heterochromatin domains and LAD thickness to further validate the model results quantitatively (Fig.6d). Cas9 ActD treated nuclei show an increased heterochromatin domain radius compared to control while WAPL nuclei show a significant reduction in domain radius and LAD thickness (Fig.6d). However, WAPL ActD treated nuclei show no significant difference in comparison to Cas9 ActD treated nuclei (Fig.6d), in quantitative agreement with the numerical predictions.

These results further confirm that the effect of transcription on global chromatin distribution occurs via supercoiling mediated chromatin loop extrusion, especially at the interface of heterochromatin and euchromatin phases. Furthermore, these results also present a significant validation of the mathematical phase-field model of chromatin organization in the nucleus.

Read this article:
Active transcription and epigenetic reactions synergistically regulate meso-scale genomic organization - Nature.com