Epigenetic Transitions and Knotted Solitons in Stretched Chromatin

The spreading and regulation of epigenetic marks on chromosomes is crucial to establish and maintain cellular identity. Nonetheless, the dynamic mechanism leading to the establishment and maintenance of tissue-specific, epigenetic pattern is still poorly understood. In this work we propose, and investigate in silico, a possible experimental strategy to illuminate the interplay between 3D chromatin structure and epigenetic dynamics. We consider a set-up where a reconstituted chromatin fibre is stretched at its two ends (e.g., by laser tweezers), while epigenetic enzymes (writers) and chromatin-binding proteins (readers) are flooded into the system. We show that, by tuning the stretching force and the binding affinity of the readers for chromatin, the fibre undergoes a sharp transition between a stretched, epigenetically disordered, state and a crumpled, epigenetically coherent, one. We further investigate the case in which a knot is tied along the chromatin fibre, and find that the knotted segment enhances local epigenetic order, giving rise to “epigenetic solitons” which travel and diffuse along chromatin. Our results point to an intriguing coupling between 3D chromatin topology and epigenetic dynamics, which may be investigated via single molecule experiments.


INTRODUCTION
Cell-line-specific features in multi-cellular organisms are achieved by regulation of "epigenetic marks", biochemical modifications of DNA and histone octamers which do not affect the underlying genomic sequence. Thus, the average pattern of epigenetic marks in a given-cell line well correlates with the pattern of genes which are transcriptionally active or inactive within that cell-line [1][2][3][4][5][6]. Dissecting the biophysical mechanisms leading to the de novo establishment, spreading and maintainance of epigenetic marks is consequently a key step towards a better understanding of the dynamic organisation of genomes and of chromosomal re-arrangement throughout the cell cycle [7], in cellular ageing [8] and pluripotency [9]. Notwithstanding their pivotal role, these mechanisms are still poorly understood [10].
On the other hand, it is important to realise that the establishment of epigenetic patterns is fundamentally a dynamic process, where biochemical tags are constantly deposited, removed and degraded on histones, which can themselves be displaced during transcription or replaced after replication [6,[30][31][32][33]. For this reason, simple models where epigenetic marks are stably deposited along chromatin [11, 13-15, 34, 35] are only crude approximations of a much more complex and dynamic scenario. Crucially, these "static" models fail to address key questions such as, how epigenetic patterns are first established along chromosomes, and how these change, for instance, with cellular ageing [8] or during disruptive events in the cell cycle [31]. In addition, an understanding of cell-to-cell variability in genome organisation [36,37] and of efficiency of the induced-pluripotency pathway [38] may be achieved only through models which can take into account the plasticity of the epigenetic landscape.
In recent years, wide-spread technological advances in the field of molecular biology allowed the biophysical community to identify some of the key players in the dynamics of epigenetic patterning [5,8,31,39]. At the heart of this process are proteins which can "read" and "write" biochemical tags along chromatin [6,31,40]. Importantly, some of these proteins are found in the same complex [31] or are known to directly recruit one another [41]. For instance, the heterochromatin binding protein HP1 (a reader) possesses a chromodomain recognising tri-methylation of Lysine 9 on histone 3 (H3K9me3) [42][43][44]; at the same time, it can recruit the methyltransferase complex SUV39h1 [31] (a writer), therefore engaging a positive feedback loop which self-perpetuates this repressive mark [41,45,46].
Computational models of this positive feedback loop in 1D [47,48] have shown that accounting for long-range contacts is necessary to allow spreading of repressive marks. On the other hand, 3D models coupling epigenetic and chromatin dynamics [49] have shown that spreading of silencing marks on a mobile 3D fibre can be viewed as an "all-or-none" transition, where a local fluctuation triggers an epigenetic wave of repressive marks which stabilises a compact globular state [49].
Existing experimental and computational studies suggest that the feedback between 3D chromatin structure and epigenetic dynamics along the chromosome may be a key potential mechanism underlying the establishment of epigenetic patterns which can display memory of their landscape. However, direct experimental observations of dynamical epigenetic and chromatin conformational changes in the nucleus are largely missing, in view of the difficulty to achieve enough spatio-temporal resolution in vivo. In light of this, here we propose, and investigate in silico, a novel avenue to directly observe the coupling between epigenetic dynamics and 3D chromatin folding.
The system we consider can be recreated in vitro via single-molecule experiments: it is inspired both by experiments studying the stretching of reconstituted chromatin [50][51][52], and by investigations of protein-DNA and protein-chromatin interactions in vitro [53,54]. Specifically, we envisage a set-up where reconstituted chromatin is attached at both ends to macroscopic beads so as to exert stretching forces on the strand via optical or magnetic tweezers [50,52]. Further, we imagine that the chromatin fibre is embedded in a solution where reader and writer proteins are added and activated. In practice, one may achieve this by including "readers" such as heterochromatin HP1 proteins [55] (which are known to be multivalent, hence act as bridges which fold chromatin), and "writers" such as SUV39h [31]. These proteins can respectively read and write the repressive H3K9me3 mark [41,55]. An alternative option is to use Polycomb PRC complexes and Enhancer-of-zeste (E(z)) proteins which are respectively able to bind and deposit the H3K27me3 mark on histones [31,39,[56][57][58].
Our results show that by tuning the strength of the external stretching force f and the binding affinity of reader proteins , the system can display either a stretched, epigenetically disordered, state (SD) or a compact, epigenetically ordered, one (CO). The two regimes are separated by an abrupt transition line which we quantitatively locate in the (f, ) parameter space. We further show that when knots are tied along the chromatin, these can readily be identified by looking at the local order of epigenetic patterning. We argue that the topology of the knots localise and enhance 3D chromatin interactions, thereby "protecting" the epigenetically ordered region, which remains localised within the knotted arc. We dub these novel and remarkable states "epigenetic knotted solitons".
We finally envisage that our findings will inform the design of novel single-molecule experiments in vitro and illuminate the interplay between epigenetic dynamics and chromatin topology in vivo.

Chromatin Model and Experimental Set-up
We model the chromatin as a "recolourable" bead-spring polymer (for further details see Methods). Each bead has a size σ -for definiteness, we take σ = 3 kb or 30 nm [59][60][61], but our results do not depend on this precise choice -and it bears a "colour" representing a specific biochemical tag (see Fig. 1). For instance, one may think of blue beads as representing chromatin regions with an excess of heterochromatin mark, H3K9me3, and red beads as ones with excess of Polycomb mark, H3K27me3. We further include grey beads which represent unmarked regions of chromatin. For simplicity, in this work we consider only these three colours, in agreement with the generic experimental observation that few epigenetic "states" are sufficient to well capture the overall epigenetic landscape in several organisms [28,34,[62][63][64]. [92] Beads interact through a Weeks-Chandler-Andersen potential, which strongly suppresses any overlap between beads. This interaction is further modified to account for attractive interactions (with affinity ) between beads bearing the same colour. This effectively mimics the action of reading proteins, or "bridges", which can bring together distant chromatin segments bearing the same epigenetic mark. Unmarked (grey) or differently marked beads are therefore considered to be solely sterically interacting. The truncated Lennard-Jones (LJ) potential can be written as a function of distance x between, and the colour q of, the interacting beads a and b as follows (1) for x ≤ x qaq b c , whereas U ab LJ (x) = 0 otherwise (N is a normalisation constant, see Materials and Methods). In Eq. (1) x qaq b c is the cut-off distance between beads a and b with colors q a and q b respectively. This potential is a computationally efficient way to model repulsive and attractive interactions between beads. Specifically, it can display an attractive (negative) region by setting the cut-off larger than the minimum of the potential or can yield pure repulsion otherwise. In practice, we use the knowledge of the colours q a and q b to define a colour-dependent cut-off as The shift (last two terms in eq. (1)) ensure that there is no discontinuity in the potential. The binding affinity ab = when q a = q b = 0 and is a free parameter of our model, otherwise ab = k B T L when q a = q b or one of the two is unmarked.
The "recolouring" process is modelled via a Monte Carlo procedure which occurs at an inverse rate of k −1 R = τ R = 10 3 τ Br 10 s (see Materials and Methods). Every time a bead (a) is selected for a recolouring attempt, its colour Chromatin is modelled as a coarse-grained "recolourable" bead-spring polymer; each bead represents few nucleosomes and it is coloured either red, blue or grey (unmarked bead). We implicitly model reader and writer proteins interacting with chromatin: the former allow polymer folding by bridging between segments bearing the same colour, the latter deposit new marks on nucleosomes and allow beads to change colour. (B) As red/blue marks are deposited on beads, respective reader proteins can bind to them; in turn, they recruit writer proteins, which deposit new red/blue marks. This mutual recruitment triggers a positive feedback loop. (C) Reader and writer proteins for red and blue marks compete over the same chromatin segments thereby generating dynamic epigenetic marking. (D) Experimental set-up: macroscopic beads (green) are attached at both chromatin ends and allow the strand to be stretched by external tweezers with force f [50,52]. Starting from a configuration in equilibrium with the external force and with no proteins in solution, the system is flooded with reader/writer proteins thereby triggering the competition between red and blue marks (in the model, this is effectively done by allowing attractive interactions between chromatin segments). The system may eventually evolve towards an equilibrium state where chromatin is compact and only one epigenetic mark dominates ( > c(f ), shown in the figure), or towards one where chromatin is stretched and no coherent mark is established ( < c(f ), not shown).
is randomly changed into one of the two remaining colours and the new energy is computed. Because the only colourdependent potential employed in the simulations is the one written in Eq. (1), the difference in energy is given by The recolouring attempt is then accepted with a Metropolis probability where T R is a generic temperature that can, in principle, be tuned according to the efficiency of the writing process and is therefore independent on the solution temperature T L , which regulates the stochastic dynamics of the beads in 3D [49] (see Materials and Methods). The recolouring process and the motion of the chromatin in 3D are engaged in a positive feedback loop: beads that are near each other in space are more likely to become equally coloured and beads bearing the same colour are more likely to stick together. Thus, both species (red and blue) compete over the chromatin strand in a similar way that "up" and "down" spins compete over a mobile string, or a "magnetic polymer" [66].
Here, we study the behaviour of the system upon tuning the stretching force f and the attraction strength , since this protocol may be realised in vitro by using reconstituted chromatin [50,52] and proteins such as HP1 and SUV39 [45,46], as previously mentioned. For simplicity, we limit to the case T L = T R which ensures that the epigenetic read-write mechanism and the chromatin folding are governed by transition rules between different microstates that obey detailed balance and that can be described in terms of an effective free energy. Considering out-of-equilibrium conditions [49] leads to transitions between states with similar epigenetic patterns, i.e. from swollen ordered (or disordered) to compact ordered (or disordered), which can be understood as the homopolymer (or heteropolymer) limit of our system. Because these transitions do not shed light onto the interplay between epigenetics and chromatin conformations, we here decide not to pursue them.
In order to quantify the behaviour of the system we perform Brownian Dynamics simulations for a typical runtime of 5 10 6 τ Br -corresponding to 1000 τ R , or Monte Carlo sweeps -of a chain L = 1000 3 Mb beads long (if not specified otherwise). Although chromatin has not yet been reconstituted to such a length, the trends we uncover are generic and also hold for smaller values of L.

Epigenetic and Conformational Transitions of Stretched Chromatin
We initialise the systems by assuming that no reader or writer protein is present at t < 0 and by letting the chromatin equilibrate while subject to the stretching force f and in a bath at temperature T L . In practice, we do this by lowering the cut-off for the Weeks-Chandler-Andersen potential written in Eq. (1) to 2 1/6 σ for any pair of beads so that the interactions between beads are purely repulsive.
At time t = 0 the system is flooded with reader and writer proteins which can bridge beads with the same epigenetic mark (i.e., colour) and to attempt recolouring of the beads (i.e., writing) at a rate τ −1 R . To visualise the epigenetic evolution of the system, in Fig. 2(A-B) we report typical "kymographs", which show the colour (q) of each bead along the polymer at a given time t. Kymographs readily capture the epigenetic "ordering" of the chromatin strand as a function of time for a specific choice of and f . These two latter quantities are the main free parameters of the system: the stretching force f can be controlled in vitro through either optical or magnetic tweezers, while the binding affinity between readers and epigenetic marks may be varied by considering mutants of bridging proteins such as HP1 [45,55].
In order to further quantify the behaviour of the system in parameter space, we systematically vary f and and quantify the equilibrium states by measuring the values of the absolute "epigenetic magnetisation" and of the average chromatin extension via its radius of gyration where r i and r CM are the positions of segment i and of the centre of mass of the chain, respectively. In order to construct an equilibrium phase diagram of the system, we estimate the average R g ≡ R 2 g 1/2 and m , first by time-averaging over a trajectory at steady state and then averaging the results across several trajectories. By measuring these two observables for different choices of and f we report the phase diagram shown Fig. 2(C). One can readily notice that the typical equilibrium configurations can be grouped into two distinct phases: one is compact and epigenetically ordered (CO), i.e. with whereas the other is stretched and epigenetically disordered (SD), i.e. with R g ∼ Lσ, m 0 .
These two regimes are separated by a transition line and one can readily appreciate from the profiles of R g and m (Figs. 2(D-E)) that this transition is abrupt. In the limit of stretching forces f → 0, we retrieve the firstorder transition observed in Ref. [49] for un-stretched chromatin. When larger stretching forces are applied, the system therefore displays a force-dependent critical line c (f ) (or f c ( )) which retains the first-order-like features observed at f = 0. Within the CO region, there is a parameter range where a multi-domain epigenetic structure emerges dynamically. This can be readily seen from the profile of the epigenetic magnetisation m in Fig. 2(E): near the transition line (from the compact ordered side) this quantity displays a sharp peak which then decreases for smaller values of the force. In the region labelled as "multi" in Fig. 2(E), the magnetisation is lower than unity because multiple macroscopic ordered regions populated by different epigenetic marks compete with one another, thereby lowering the overall magnetisation. These domains must be metastable (though long-lived) as in steady state a single domain is preferred since it removes domain walls, which have a free energy cost. The existence of multi-domain patterns arises because, at stretching forces far from critical f c ( ), there may be multiple nucleation points along the fibre for the spreading of epigenetic marks. These nucleation points are generated by a local increase in chromatin density, for instance through the transient formation of loops and 3D interactions, which then trigger local spreading of different epigenetic marks (see Fig. 3). As the stretching force gets closer to the critical value f c ( ), the nucleation probability declines, hence only one epigenetic mark is able to take over the whole chromatin strand, thereby enhancing its overall epigenetic magnetisation.
It is useful to compare the epigenetic transition between the CO and the SD state in Fig. 2 with the equilibrium transition between a compact and a stretched state which can be observed in a stretched homopolymer with selfattractive interactions (i.e., in a poor solvent). While the transition is first order in the homopolymer case as well [69,70], here the epigenetic degrees of freedom increase the entropy of the disordered phase, thereby the value of the critical force, for a given , is lower. Another key difference is that multi-pearl structures are only observed transiently close to the transition for the homopolymer case [69], whereas in the epigenetic case multi-domain  (41) knot [67], whose escape from the fibre is avoided by the macroscopic beads at the terminal ends. Superimposed on the kymograph, we also show the boundaries of the knot computed from the 3D polymer configurations through the knot identification algorithm described in Ref. [68] (black lines). Intriguingly, from this plot one can immediately realise that the knotted portion of the chain can be identified with the epigenetically ordered region, i.e. the chromatin segment where the beads are homogeneously coloured. (B) Snapshot of the chromatin fibre at the time-step highlighted in the kymograph. (C) Analysis of the knotted soliton for different system replicas. In the main panel we show the time-averaged magnetisation m(x) at position x relative to the centre of the knot (see Materials and Methods); this displays a localised increase which we identify as an "epigenetic soliton". In the inset we show the probability of observing a knotted arc of length lK through the knot identification algorithm of Ref. [68]. As one can notice, the length of the epigenetically ordered region and that of the knotted arc are in near-perfect agreement. This strongly suggests that the localised knotted segment leads to the local coherency in the epigenetic marks. states arise far from the transition, and they are long-lived. Both these differences should be experimentally detectable, as the homopolymer case may be recreated by reconstituting chromatin fibres with histone octamers with controlled biochemical tags or in the presence of linker histones (H1), while avoiding the action of writers.

Long-Lived Epigenetic Domains and Boundary Diffusion
The competition of epigenetic marks over a chromatin segment is a phenomenon relevant to many biological systems. For example, genes that are positioned near telomeres may switch between transcriptionally active and silenced states from one generation to the next [71]. This phenomenon, generally referred to as "telomeric position effect", is particularly relevant in yeast [71] but it has also been observed in human cells [72]. It occurs because heterochromatin marks largely populating the telomeric gene-poor regions of the genome invade gene-rich stretches of chromatin, effectively switching off gene transcription.
In our proposed single molecule set-up, the emergence of multiple competing epigenetic domains in reconstituted chromatin under (small) tension allows a quantitative study of stability and diffusive dynamics of boundaries between neighbouring epigenetic domains, processes which may be relevant to understand the mechanisms underlying competition and spreading of epigenetic marks in vivo.
In practice, an epigenetic domain in the stretched chromatin fibre is defined as a macroscopic region (e.g., consisting of w ≥ 50 beads) over which more than θ = 90% of the beads are homogeneously coloured. A boundary between domains is found if the beads on either side of the boundary have opposite colours (i.e., red and blue) and the signed magnetisation difference ∆m of the left and right domain is such that |∆m| > 2θ (i.e. the domains have opposite signed magnetisation). time, one can measure their mean square displacement

By tracking the position of the domain boundaries over
where x(t) is the position of a given boundary at time t and T the total measurement time. Plots of δ 2 x(t) for different replicas and two choices of interaction energies are reported in Fig. 3(C-D). Data points are grouped according to the applied stretching force f , showing that the dependence on f is much weaker compared to that on . We extract the diffusion coefficient D as the long time behaviour of δ 2 x(t), i.e.
which is reported in Fig. 3(F) and shows an exponential decay of D as a function of . In our stretched chromatin assay, the diffusive dynamics of epigenetic boundaries is therefore mainly controlled by the strength of the attractive interactions between beads bearing the same epigenetic mark. In vivo, another possible important factors affecting diffusivity of epigenetic boundaries may be the presence of insulators and architectural proteins such as CTCF and cohesins [1,53,54]. These features may be included in future studies focused on understanding the actions of these architectural elements.

Epigenetic Knotted Solitons
The abrupt, first order, phase transition between the compact-ordered and the stretched-disordered state shown in Figure 2 is here due to the coupling between global 3D structure and 1D epigenetic dynamics [49], which gives rise to a positive feedback loop where compaction aids spreading, which leads to further crumpling at small enough f , or sufficiently large .
It is interesting to ask whether one may design an experimentally realisable situation whereby this positive feedback is only realised locally, rather than globally as done in Ref. [49]. To do so, we consider a knotted chromatin fibre: the idea behind this construct is that a polymeric knot tightly localises upon stretching [74,75], and that the knotted region is more compact, thereby providing a natural nucleation point for epigenetic ordering. In the single molecule experiment set-up we consider (Fig. 2), a physical knot may be tied within chromatin by using, for instance, micro-manipulating techniques used to generate knots on DNA [76][77][78]. Alternatively, a knotted chromatin fibre may be self-assembled by first knotting the naked DNA strand, and later allowing nucleosome formation by adding histone octamers to the solution [79].
In Figure 4 (A-B) we show the kymograph for a chromatin strand along which a figure-of-eight (4 1 ) knot is tied. As one can see, even though we choose values of and f so that the system is in the stretched-disordered region of the phase diagram (see Fig. 2), the kymograph clearly displays a localised epigenetically ordered region. To confirm whether this region corresponds to the knotted arc, tightened upon stretching and self-attraction, we pinpoint and monitor the time evolution of the knotted region by resorting to a top-down search of the smallest portion of chromatin that yields, upon suitable closure, a ring with the same topology of the whole chain [68]. Remarkably, the knotted portion found with this well-established algorithm perfectly matches the epigenetic domain in the kymograph (black lines in Fig. 4 identify the boundary of the knotted arc, see also Methods).
We call this remarkable configuration an "epigenetic knotted soliton", because, as a soliton, it travels along the fibre by keeping near-constant shape, hence displaying particle-like behaviour. Similarly to solitons recently observed in meta-materials [80] and minimal surfaces [81], it is the non-trivial topology of the localised knot which keeps the structure together. The epigenetic state of the A B Figure 6. Epigenetics and 3D topology. This figure intuitively shows the reason behind the formation of a stable epigenetic domain within the knotted chromatin segment: the local 3D topology of a knot "protects" the underlying epigenetic information by enhancing local 3D interactions (B). On the contrary, an unknotted stretched segment lacks long-range interactions required to sustain epigenetic coherency (A). knotted soliton can "switch" from red to blue, or vice versa -switching events occur, in Figure 4, at about 2 10 5 τ Br and 3.5 10 5 τ Br ). Switching between the two stable states is observed because the soliton has a finite size, hence the effective free energy barrier separating the red and blue states can be occasionally overtaken by fluctuations.
The physical properties of the soliton are tunable: for instance, by increasing the stretching force the knot tightens and becomes more compact, so that the soliton shrinks in size but becomes more stable due to the increase of 3D interactions which locally "protect" the epigenetically ordered state. By increasing values of , such that the system remains overall disordered ( < c (f )) the switching rate decreases (i.e., the epigenetic ordering is more robust, see Fig. 4(D)) -the diffusion coefficient, however, does not depend on appreciably (unlike for the case of epigenetic domains, see Fig. 3(E-F)). It would be of interest to quantify how the soliton diffusion coefficient depends on knot type, and whether, as in knots in swollen polymers, twist knots and more complex knots are less mobile [67,77,82,83].
For values of > c (f ), we observe the same SD-CO phase transition reported in Figure 2, where the localised knot now acts as a nucleation point for the spreading of its epigenetic mark.
Finally, we investigate how knotted solitons interact on a chromatin strand. We consider a system where two knots, a pentafoil 5 1 and a trefoil 3 1 , are initially tied along the chromatin. As the system evolves each knot forms its independent epigenetic soliton. As shown by the kymographs in Figure 5(A-D), both solitons diffuse along the chain [77]; they also interact and may merge to form a single domain, and later on split again (Fig. 5(A)), or even cross each other [73] (Fig. 5(D)). Following every collision or interaction, each soliton retains its overall structure and reestablishes its locally coherent epigenetic mark.

DISCUSSION AND CONCLUSIONS
In several areas of chromosome biology, researchers assume the existence of a tight coupling between 3D chromatin structure and epigenetic dynamics, as this coupling provides an appealing mechanism for the de novo establishement and maintainance of epigenetic patterns [47,49,55,56,58]. However, as both epigenetic spreading and chromatin dynamics may occur on comparable and relatively fast timescales in vivo (minutes to hours [31,84]), it is extremely difficult to design an experiment to demonstrate this coupling dynamically inside the cell. Here, instead, we proposed, and investigated in silico, a set-up for an experiment to test and quantify the coupling directly in vitro. The set-up involves a reconstituted chromatin fibre [50,52], which is stretched (e.g., by laser tweezers) and interacts with an ensemble of reader and writer proteins -such as heterochromatin HP1 and SUV39h1 complexes [45,55] or polycomb [56,58] enzymes able to bind to and deposit their respective marks.
We showed that by varying, for instance, the stretching force, one can trigger a phase transition between a compact-ordered phase, where a single epigenetic mark invades the whole chromatin fibre, and a stretcheddisordered phase, where no single mark can take over the system. It is important that this transition may be observed by varying only the stretching force, as this is simpler to fine tune continuously in vitro with respect to temperature (which can inactivate readers or writers altogether) and effective self-attraction between chromatin segments (which may be achieved through the use of mutant reader proteins).
By using single-molecule imaging [53,54,85] or superresolution techniques [86] one may envisage to employ antibodies recognising H3K9me3 and H3K27me3 epigenetic states to separately stain the marks and thereby discriminate between an epigenetically ordered and an epigenetically disordered state. In principle, one may even be able to observe the transition from one to the other by increasing the stretching force in the single molecule assay. Another possibility to assess the degree of epigenetic order is to use single-cell chromatin immunoprecipitation [87] (ChIP) on the reconstituted chromatin, to quantify histone modifications along the fibres. Note that, for this approach to be viable, one would have to use genomic DNA to reconstitute chromatin (rather than repeating, such as 601, sequences), so that locations along the DNA can be mapped uniquely. Therefore, whilst an experimental investigation of the compact-ordered to stretched-disordered transition is certainly challenging, and will require stateof-the-art experimental techniques, it is in principle feasible, and would constitute the first direct measurement of the so far elusive coupling between chromatin structure and epigenetic dynamics. A more experimentally accessible insight (although less informative regarding the epige-netic state of the system) may also be achieved via opticaltweezers by measuring force-extension curves [50,52,88]. These can discriminate between a (partially) collapsed and a stretched coil [88], while cannot determine whether the chromatin is epigenetically ordered or disordered. On the other hand, the typical force-extension curves obtained with this set-up may display sensitivity on the concentration or efficiency of epigenetic readers and writers thereby offering indirect quantification of the recolouring process.
By tying a knot along the reconstituted chromatin fibre we observed that the system can harbour a novel structure, an "epigenetically knotted soliton" (see Figs. 4-6). This is a tight knot which is locally ordered epigenetically and diffuses freely, within an epigenetically disordered background. The size of the soliton is tunable by varying the stretching force applied to the chromatin fibre, and different solitons interact in a variety of ways: they may bounce off one another or cross through each other when they collide. We also expect similar topological solitons to be universally found in knotted magnetic polymer, a new kind of topological soft matter which has not yet been realised experimentally [66].
The findings we have reported here may also be employed to detect knots in chromatin strands. Starting from an ordered and crumpled state (such as heterochromatin [1]), one may imagine to apply and continuously increase an external stretching force f and monitor the evolution of the system; whereas ordered unknotted regions will undergo an abrupt phase transition and become epigenetically disordered (see phase diagram in Fig. 2), knotted chromatin will instead preserve a localised epigenetically ordered region that can be identified with the knotted segment (Fig. 4).
Besides all this, the epigenetic solitons provide another observable consequence of the dynamic coupling between 3D chromatin structure and epigenetics which can be tested in future single-molecule experiments.

METHODS
The datasets generated during and/or analysed during the current study are available from the corresponding author on reasonable request.

Chromatin Model
Chromatin is modelled as a bead-spring polymer chain [89] where each bead has nominal size σ. Attractive and repulsive interactions are controlled by the potential described in the main text (Eq. (1)) and plotted in Fig. S1. This is a truncated-and-shifted Lennard-Jones potential which is broadly used in Molecular Dynamics [90] to model short-ranged interactions. We include a colourdependence in order to model like-colour attraction and different-colour repulsion as described in eq. (1) and shown in Fig. S1. Harmonic springs between beads a and b are imposed as U ab H (x) = (δ b,a+1 + δ b,a−1 )k B T L κ(x − x 0 )/2 with x 0 = 1.1σ and κ = 100k B T L to ensure the connectivity of the backbone. The chain stiffness is regulated by a Kratky-Porod potential between triplets of beads forming an angle θ = (t a · t b ) / (|t a | |t b |) where t a is the vector joining beads a and a + 1 as U ab KP (θ) = (δ b,a+1 + δ b,a−1 )k B T L l p (1 + cos θ) /σ with l p = 3σ. The total potential U a (x) experienced by each bead is given by the sum over all the possible interacting pairs and triplets, i.e.

U a(x)
The dynamics of each bead therefore obeys the Langevin equation where γ = 1 (in dimensionless LJ units) is a friction coefficient and ξ a (t) is stochastic delta-correlated noise which obeys the fluctuation dissipation relationship ξ α,a (t)ξ β,b (t ) = 2γk B T L δ(t − t )δ αβ δ ab where the Latin indexes run over particles and Greek indexes over Cartesian components. The mass of the beads are taken to be unity in the dimensionless Lennard-Jones units [90]. The beads employed at the ends of the chromatin are five times larger than the beads forming the polymer and have the same mass. In order to simulate the stretching of the chromatin, we directly apply a force f on these beads. Finally, Eq. (11) is integrated with a velocity-Verlet scheme within the LAMMPS [91] engine. Simulations units can be mapped to real ones by considering that the polymer beads can coarse-grain groups of nucleosomes. The nominal size of the beads can is here considered to be 30 nm, as capturing the thickness of reconstituted chromatin fibre [50,[59][60][61]. The typical timescales over which each bead diffuses its own size is therefore where we used the Einstein relation D self = k B T L /3πησ and we considered η = 150cP as the effective viscosity. The recolouring process occurs at inverse rate τ R = 10 3 τ Br . In each attempt, a bead is selected at random and a colour change is proposed; the move is then accepted if satisfying the Metropolis criterion as described in the text. In other words, on average, all beads will have been attempted to change their colour every τ R = 10 3 τ Br steps.
For simplicity, we model reader and writer proteins implicitly. This is effectively done by setting attractive pairwise interactions between beads bearing the same colour and by assuming a uniform selection probability during the recolouring attempts. In other words, the physical proximity of reader and writers is not required for bridging segments and changing the beads colour. This model is therefore a good approximation for the case in which reader and writer proteins saturate the solution. On the other hand, lowering the concentration of reader and writer proteins would effectively slow down the kinetics towards the equilibrium state: this situation is more complicated and we defer it to future studies.

Knotted Soliton
To test the stability of the knotted soliton we consider the time-averaged degree of epigenetic coherence m(x) at position x relative to the centre of knotted arc, c K . This can be defined as where ∆τ is the typical time between epigenetic "switching", where the magnetisation of the whole knot changes sign (such as the one occurring in Fig. 4 at about 2 10 5 τ Br and 3.5 10 5 τ Br ) and m(l, t) is the magnetisation at position l and time t. In Fig. 4(C) we report m(x) for three system replicas and compare it with the distribution of knot sizes obtained from the 3D polymer configurations (inset). As one can notice, the length of the ordered domain agrees with the length of the knotted arc. Further, the shape of the magnetisation profile m(x) displays a broad peak, which is characteristic of soliton-type solutions.