Cross-linker mediated compaction and local morphologies in a model chromosome
CCross-linker mediated compaction and local morphologies in a model chromosome
Amit Kumar
1, 2, ∗ and Debasish Chaudhuri
1, 2, † Institute of Physics, Sachivalaya Marg, Bhubaneswar 751005, India Homi Bhaba National Institute, Anushaktigar, Mumbai 400094, India (Dated: May 22, 2019)Chromatin and associated proteins constitute the highly folded structure of chromosomes. Weconsider a self-avoiding polymer model of the chromatin, segments of which may get cross-linkedvia protein binders that repel each other. The binders cluster together via the polymer mediatedattraction, in turn, folding the polymer. Using molecular dynamics simulations, and a mean fielddescription, we explicitly demonstrate the continuous nature of the folding transition, characterizedby unimodal distributions of the polymer size across the transition. At the transition point thechromatin size and cross-linker clusters display large fluctuations, and a maximum in their negativecross-correlation, apart from a critical slowing down. Along the transition, we distinguish the localchain morphologies in terms of topological loops, inter-loop gaps, and zippering. The topologies aredominated by simply connected loops at the criticality, and by zippering in the folded phase.
I. INTRODUCTION
Chromosomes consist of DNA and associated proteins.Long DNA chains with lengths that vary from millime-ters (bacteria) to meters (mammals) are compacted andorganized within micron sized bacterial cells or cell nucleiof eukaryotes. The DNA must be compacted by severalorders of magnitude, and still allow information process-ing in terms of gene expression in both pro- and eukary-otes [1, 2], and also replication in bacteria [3]. DNA asso-ciated proteins play a crucial role in such processes [1–4].In the smallest scale, DNA double helix wraps aroundhistone octamers forming a bead on string chromatinstructure with a connected set of nucleosomes [5]. Inbacteria, histone like nucleoid structuring (H-NS) pro-tein dimers bind to DNA [6, 7]. It may be noted that thepositive charges on most DNA binding proteins providenon-specific affinity to negatively charged DNA [1, 2].Higher order structure formation involves bringing to-gether of contour-wise distant parts of the chromatin intospatial contact to form loops [8–13]. This is observed inall domains of life, in bacteria [14, 15], archea [16] andeukaryotes [12, 17–19]. Such loops may be maintained byproteins cross-linking spatially proximal chromatin seg-ments [20–22, 22–31], or by extrusion [32–35]. A numberof proteins are identified that stabilize these loops intoseparate topological domains [36–38]. In eukaryotes co-hesin and CTCF are identified as chromatin loop regu-lators [12, 18, 19, 32]. In bacteria, loops are stabilizedin part by non-specific cross-linkers such as H-NS, Lrpand SMC proteins [13, 39, 40]. First direct evidence ofthe chromosomal loops were found in electron microscopyexperiments [41–44]. Complementary experiments us-ing chromosome conformation capture techniques pro-vide contact maps that exhibit spatial contacts betweendifferent genes along the DNA contour [45–47]. The bio- ∗ Electronic address: [email protected] † Electronic address: [email protected] logical function of chromosomes are often related to theirlocal morphology [3, 36], as spatially proximal genes, irre-spective of their location along the DNA contour, can beregulated together [48–53]. Given their structural com-plexity, morphologies of long folded chains are often an-alyzed in terms of their generic topological features [54–57].During interphase chromosomes display several univer-sal properties, e.g., scale free nature in subchain exten-sion [58], and average contact probabilities [45] similarto homopolymers but with exponents that differ fromsimple chains [59, 60]. At small separations the chromo-somal architecture is determined by its bending rigidity,while at long range they show behavior typical of frac-tal globules [45, 61–66]. The contact maps reveal topo-logically associated domains at smaller genomic separa-tions ( (cid:46) a r X i v : . [ c ond - m a t . s o f t ] M a y ecent mean field approach that incorporates fluctuationsof co-solvent density showed that the nature of the coil-globule transition with the increase in co-solvent densitydepends on the the kind of polymer- co-solvent interac-tion. When the interaction strength between the polymerand co-solvent is purely repulsive the predicted transitionis continuous, whereas it turns out to be a first ordertransition if the interaction is purely attractive [88].In this paper, we model the chromatin as a self-avoiding chain and explicitly consider its non-specific at-traction with diffusing binder proteins that repel eachother. We perform molecular dynamics simulations inthe presence of a Langevin heat bath to fully character-ize the binder mediated folding transition, and associatedlocal topologies of the chromosome. The binders cross-link different segments of the chromatin, bringing themtogether. As a result more binders accumulate, formingclusters. Such clusters, in turn, fold the polymer. Us-ing numerical simulations and a mean field descriptionwe show that the folding is a continuous transition, me-diated by a linear instability towards formation of largeclusters of cross-linkers. While the linear stability pre-diction shows reasonable agreement with simulations forgrowth of cluster size before the transition, the mean fieldprediction for chromosome size display better agreementafter the transition, as fluctuations get suppressed in theglobule. The polymer size distribution shows a singlemaximum across the transition, signifying that there isno metastable phase on the other side of the transition.The criticality is characterized by large and slow fluctua-tions – the fluctuation amplitude and time-scale increasewith system size. At criticality, the cross-correlation be-tween chromosome size and cross-linker density shows anegative maximum.We further analyze topologies of the chromatin loopsthat are formed due to binder cross-linking, identifyingthe simply connected and higher order loops. The av-erage number of simply connected loops show a maxi-mum at the critical point, while the relative probabili-ties of loops of different orders change qualitatively withincreasing cross-linker density. The first order loop-sizesshow power law distributions with exponents that changemonotonically across the transition. The gaps betweensuch loops, in contrast, follow exponential distributions.The mean gap size hits a minimum at the critical point.Apart from forming loops, the binders may also zippercontiguous segments belonging to different parts of thechromatin. As we show, the mean number of zipperingdisplays a sigmoidal behavior along the coil-globule tran-sition, with saturation at large binder densities.In Sec. II, we present the model and details of numer-ical simulations. In Sec.III we discuss simulation resultsidentifying the coil-globule phase transition in the modelchromatin chain, and clustering of polymer-bound cross-linkers. The transition and clustering are interpretedin terms of a mean field description and linear stabil-ity analysis. Finally, in Sec. IV we characterize the localmorphology of chromosomes in terms of contact proba- s(a) (b) (c)(d)(e) FIG. 1: (color online) Representative snapshots of the modelchromosome with N = 256 bead chain and binders at thetransition point φ c = 1 . × − , where large conformationalfluctuations are observed. ( a ) A relatively compact conforma-tion. The chromatin is shown by blue monomers connected bybonds. The cross-linkers attached to the chromatin are shownas red beads, while the freely diffusing binders are shown asgreen beads. ( b ) One relatively open conformation. ( c ) Amagnified portion of ( b ) shows a contact formation denoted bythe aqua-green bar. A monomer pair, contour-wise separatedby s , have come within r c forming the contact. ( d ) A magni-fied portion of ( a ) shows loop formation by a polymer boundcross-linker (red bead). The red bars indicate the bonds thatit forms. The line with arrow-heads identifies a simply con-nected loop (for further details see Sec. IV A). ( e ) Clustersof polymer bound cross-linkers in ( a ). For better visibilityof cross-linkers, the chromatin is represented by a faded line.The thick dashed circle identifies one cluster. bility, loop topologies of various orders, and zippering.We conclude presenting a connection of our results toexperimentally verifiable predictions in Sec. V. II. MODEL
We use a self-avoiding flexible chain model of thechromatin. The bead size is assumed to be largerthan the Kuhn length. The chain connectivity is main-tained by finitely extensible nonlinear elastic (FENE)bonds between consecutive beads, U FENE ( r i +1 ,i ) = − k R ln[1 − ( r i +1 ,i /R ) ] where k and R fix the bond,and r ij = | r i − r j | denotes separation between i -thand j -th bead. The self avoidance is implemented viathe Weeks-Chandler-Anderson potential, U WCA ( r ij ) =4 (cid:15) [( σ/r ij ) − ( σ/r ij ) + 0 .
25] between beads separatedby distance r ij < / σ , else U WCA ( r ij ) = 0 [89]. Thus2 . . . . a ) . . .
28 0 0 . . . ( i ) φ r = 0 . φ r = 2 . ii ) h R g i φ r ∆ R g R g φ r R g t × − . . . . b ) . . . . . . ( i ) φ r = 0 . φ r = 3 . ii ) h C s i φ r ∆ C s C s φ r C g t × − . . . . c ) − − −
50 0 . . . . d ) τ c × − φ r C R g , n a φ r FIG. 2: (color online) The coil-globule transition as a function of ambient cross-linker density φ c , expressed in terms of φ r = φ c × . ( a ) The decrease in mean radius of gyration of polymer (cid:104) R g (cid:105) with φ c , the data are shown by (cid:5) and errorbars, captures a coil-globule transition. The blue dashed line is a guide to eye. The mean field prediction (cid:104) R g (cid:105) = constant= (cid:104) R g (cid:105) ( φ c = 0) at φ < φ ∗ c = 1 . × − . At higher densities, simulation results for (cid:104) R g (cid:105) show reasonable agreement withEq.4 with fitting parameter u/v = 0 .
1. The two curves are shown by green solid lines. At the transition point φ ∗ c , relativefluctuation of polymer size ∆ R g /R g shows a maximum (inset ( i )). The equilibration of R g with time t at two densities φ c = 0 . × − (pink), 2 . × − (green) are shown in inset ( ii ). ( b ) The mean size of the polymer bound clusters of cross-linkers (cid:104) C s (cid:105) increases with φ c . The data are shown by (cid:5) and error bars. The blue dashed line is a guide to eye. At φ c < φ ∗ c ,the data show reasonable agreement with (cid:104) C s (cid:105) = A [(1 − φ c /φ ∗ c ] − / , using A = 1 .
9. The relative cluster size fluctuation∆ C s /C s shows a sharp maximum at the transition point φ ∗ c (inset ( i )). Inset ( ii ) shows how the instantaneous mean clustersize C g equilibrates with time t at two cross-linker densities φ c = 0 . × − (pink) , π × − (green). ( c ) Correlation times τ c = τ R g ( (cid:5) ) , τ n a ( ◦ ) are obtained from auto-correlation functions of polymer radius of gyration R g , and the total number ofchromatin-bound cross-linkers n a . They reach their maximum values at the transition point φ ∗ c . The solid (red) line is a shiftedplot of the scaling form (1 − φ c /φ ∗ c ) − / added to a constant background, with the shift aimed at better visibility. ( d ) Thenegative values of the cross-correlation coefficient C R g ,n a between fluctuations in R g and n a show anti-correlation, with theamplitude maximizing at the transition point φ ∗ c . U = U FENE + U WCA defines the polymer [90]. The repul-sion between cross-linkers are modeled through the same U W CA interaction. The energy and length scales are setby (cid:15) and σ respectively. The FENE potential is set by k = 30 . (cid:15)/σ , R = 1 . σ .The interaction between cross-linkers and monomers ismodeled through a truncated and shifted Lennard-Jonespotential, U shift ( r ) = βU LJ ( r ) − U LJ ( r c ) for r < r c and U shift ( r ) = 0 otherwise, where U LJ ( r ) = 4 (cid:15) m [( σ/r ) − ( σ/r ) ], with (cid:15) m = 3 . (cid:15) , and r c = 1 . σ . The choice of (cid:15) m is stronger than the typical hydrogen bonds (1 . k B T )and provides better stability [91], e.g., as for transcrip-tion factors [28], however, allows equilibration throughattachment- detachment kinematics over the simulationtime scales. The bond between a cross-linker and amonomer is formed if they come within the range ofattraction r c . A single cross-linker may bind to multi-ple monomers, capturing the presence of multiple DNAbinding domains in a number of regulatory proteins [28].The molecular dynamics simulations are performed us-ing the standard velocity-Verlet algorithm [92] using timestep δt = 0 . τ , where τ = σ (cid:112) m/(cid:15) is the character-istic time scale. The mass of the particles are chosento be m = 1. The temperature of the system is keptconstant at T = 1 . (cid:15)/k B by using a Langevin thermo-stat [90] characterized by an isotropic friction constant γ = 1 /τ , as implemented by the ESPResSo moleculardynamics package [93]. Similar methods have been suc-cessfully used earlier in simulation of polymers in variouscontexts [94]. Note that the diffusion of a single beadover its size σ takes a time γσ /k B T , which is the sameas the characteristic time τ .Unless stated otherwise, in this paper, we consider a N = 256 bead chain. Its typical size in absence of bindersis given by the radius of gyration R g = (13 . ± . σ .The largest fluctuations in its end to end separation arerestricted within 80 σ (data not shown). To avoid anypossible boundary effect, we perform simulations in acubic volume of significantly bigger size with sides of L = 200 σ , and implement periodic boundary condition.We vary the total number of cross-linkers from N c = 0to 6000 that changes the dimensionless cross-linker den-sity from φ c = πσ N c /L = 0 to π × − . The ap-proach to equilibrium is followed over 10 τ , longer thanthe longest time taken for equilibration near the transi-tion point. The analyses are performed over further runsof 10 − τ . A couple of representative equilibriumconfigurations are shown using VMD [95] in Fig.1 illus-trating polymer contacts, loop formation, and clusteringof cross-linkers. The system size dependence is studiedusing a restricted set of simulations, as simulating longerchains requires longer equilibration, larger simulation box3nd larger number of cross-linkers, increasing the simu-lation time significantly. III. RESULTSA. The coil-globule transition
The passive binders diffuse in three dimensions andattach to polymer segments following the Boltzmannweight. They are multi-valent, typically cross-linkingmultiple polymer segments. The probability of numberof chromatin segments that a binder can cross-link si-multaneously shows a maximum that increases from 4 to6 as the average binder concentration is increased (seeAppendix-A). This range overlaps with the typical mul-tiplicity of binding factors like CTCF and transcriptionfactories [28].As different polymer segments start attaching to across-linker the local density of monomers increases, gen-erating a positive feedback recruiting more cross-linkersand as a result localizing more monomers. Such a poten-tially runaway process gets stabilized, within our model,due to the inter-binder repulsion that ensures the binder-clusters are spatially extended. These clusters are iden-tified using the clustering algorithm in Ref. [96], and thecluster-size is given by the total number of binders in acluster. Concomitant with such clustering, the polymergets folded undergoing a coil-globule transition.Fig. 2( a ) shows the transition in terms of the decreas-ing radius of gyration (cid:104) R g (cid:105) of the model chromatin withincrease in the average cross-linker density φ c in theenvironment. The solid (green) line shows the meanfield prediction that we present in Sec. III B. The tran-sition point, φ ∗ c = 1 . × − , is characterized by amaximum in relative fluctuations of the polymer size∆ R g /R g = (cid:113) (cid:104) R g (cid:105) − (cid:104) R g (cid:105) / (cid:104) R g (cid:105) , shown in the inset ( i )of Fig. 2( a ). The equilibrations of R g at two representa-tive binder concentrations φ c are illustrated in the inset( ii ). As we show in Fig. 10 in Appendix-C, the rela-tive fluctuations ∆ R g /R g near phase transition increaseswith polymer size N , suggesting divergence in the ther-modynamic limit, a characteristic of continuous phasetransitions. In Fig.1, the large fluctuations at the phasetransition point are further illustrated with the help oftwo representative conformations: a relatively compactconformation in Fig.1( a ), and and a more open confor-mation in Fig.1( b ).The coil-globule transition occurs concomitantly withthe formation of polymer-bound cross-linker clusters. Ata given instant, several disjoined clusters may form alongthe model chromatin (see Fig.1( e ) ). The cluster size (cid:104) C s (cid:105) is the average number of binders constituting the clusters.It grows significantly as φ c approaches phase transitionfrom below (Fig. 2( b )). The linear stability estimate ofcluster size, as discussed in Sec. III B 2, is representedby the (pink) solid line in Fig. 2( b ). The relative fluc- tuations in cluster size ∆ C s /C s = (cid:112) (cid:104) C s (cid:105) − (cid:104) C s (cid:105) / (cid:104) C s (cid:105) show a sharp maximum at the phase transition point φ ∗ c (inset ( i )). Equilibration of the mean cluster size C g , in-stantaneous average over all clusters, at two φ c values areshown in the inset( ii ).The dynamical fluctuations at equilibrium are char-acterized by the auto-correlation functions C R g ( t ) = (cid:104) δR g ( t ) δR g (0) (cid:105) / (cid:104) δR g (cid:105) of chromatin size R g , and C n a ( t ) = (cid:104) δn a ( t ) δn a (0) (cid:105) / (cid:104) δn a (cid:105) of the total number ofchromatin-bound cross-linkers, where δR g ( t ) and δn a ( t )denote instantaneous deviations of the two quantitiesfrom their respective mean values (see Appendix-B). Forthe finite sized chain, the corresponding correlation times τ c = τ R g , τ n a show sharp increase at φ ∗ c (Fig. 2( c )), rem-iniscent of the critical slowing down [97]. As is shownin Appendix-C, the correlation time τ R g grows with thechain-length as τ R g ∼ N / suggesting divergence in thethermodynamic limit.The fluctuations in n a and R g are anti-correlated, asa larger number of attached cross-linkers reduces thechromatin size. Thus the cross-correlation coefficient C R g ,n a = (1 /τ p ) (cid:82) τ p dt (cid:104) δR g ( t ) δn a ( t ) (cid:105) <
0. Remarkably,the amount of anti-correlation maximizes at the criticalpoint φ ∗ c signifying a large reduction in polymer size as-sociated with a small increase of attached cross-linkers,and vice versa (Fig. 2( d )). A living cell may utilize thisphysical property for easy conformational reorganization,useful for providing access to DNA-tracking enzymes inan otherwise folded chrmosome.The probability distribution of the radius of gyra-tion P ( R g ) shows clear unimodal shape across the tran-sition (Fig. 3( a ) ). This clearly displays absence ofmetastable phase on the other side of the transition, char-acteristic of the continuous transition. Note that thisobservation is in contrast to the mean field prediction ofRef. [88], while is in agreement with the numerical simu-lations in Ref [21].Remarkably, the associated distribution of cross-linkercluster sizes P ( C s ), shows clear bimodality in much ofthe φ c range scanned across the transition (Fig.3( b ) ),capturing coexistence of clusters of small and large sizes.However, such clusters have similar densities (data notshown) and do not suggest coexistence of two phases. Infact, as we show in Sec.III B 2, the assumption of constantbinder density within the clusters provides a reasonabledescription of the growth of mean cluster size throughEq.(5) (see Fig.2( b ) ). B. Mean field description
In view of the above phenomenology, we present amean field model based on two coupled fields, the cross-linker density φ c ( r ), and the deviation of monomer den-sity due to cross-linkers ρ ( r ) = ρ m ( r ) − ρ b , where ρ b = σ N/ ( R g ) with R g denotes the radius of gyration ofthe open chain in absence of cross-linkers. A fraction ofcross-linkers are in polymer bound state, φ ( r ), and the4 − − − a ) P ( R g ) R g φ r = 0 φ r = 0 . φ r = 1 . φ r = 1 . φ r = 1 . φ r = 2 . φ r = 2 . φ r = 3 .
14 10 − − − − b ) P ( C s ) C s φ r = 0 . φ r = 1 . φ r = 1 . φ r = 2 . φ r = 2 . φ r = 3 . FIG. 3: (color online) Probability distributions of the polymer radius of gyration, P ( R g ), are plotted at different cross-linkerdensities denoted by φ r = φ c × . ( b ) Corresponding probability distributions of cluster size of cross-linkers, P ( C s ), showcoexistence of small and large clusters at densities φ c ≥ . × − . rest constitutes the detached fraction. We adopt a freeenergy density [109] βf = 12 u (cid:18) − φφ ∗ (cid:19) ρ + v ρ + κ ∇ ρ ) + 12 wφ . (1)The direct repulsion between polymer segments and be-tween cross-linkers are captured by free energy costs uρ / wφ / ρ φ ρ with strength − u/ φ ∗ . The quartic energy cost vρ / κ in thegradient term adds free energy cost to the formation ofsharp interfaces in local monomer-density. The evolutionof coupled fields are represented by [97], ∂ρ∂t = M ρ ∇ (cid:20) u (cid:18) − φφ ∗ (cid:19) ρ + vρ − κ ∇ ρ (cid:21) ∂φ∂t = M φ ∇ (cid:20) − u φ ∗ ρ + wφ (cid:21) − r ( φ − φ ) , (2)where, the second term in the right hand side of the sec-ond equation accounts for the turnover between the at-tached and detached fractions of the cross-linkers. Here r = ( r a + r d ), φ = Ω φ c with Ω = r a / ( r a + r d ). Theattachment detachment rates r a,d are determined by theinteraction and detailed balance condition. The coeffi-cients M ρ,φ denote mobilities of the two conserved fields ρ and φ . A similar approach was used earlier in Ref. [22].In the uniform equilibrium state φ = φ , and ρ = ρ .Using φ = Ω φ c and φ ∗ = Ω φ ∗ c , if φ c < φ ∗ c the solution ρ = 0 , else ρ = uv ( φ c − φ ∗ c ) φ ∗ c . (3)
1. Chromosome size
The mean monomer density ρ m = σ N/ (cid:104) R g (cid:105) = ρ + ρ b . As φ c ≥ φ ∗ c , using Eq.(3) one obtains (cid:104) R g (cid:105) = R g (cid:34) N / (cid:18) uv φ c − φ ∗ c φ ∗ c (cid:19) / (cid:35) − / . (4)This shows reasonable agreement with simulation re-sults with fitting parameter u/v = 0 . a )), asfluctuations are suppressed in the globule phase [86].In the limit of φ c (cid:29) φ ∗ c , (cid:104) R g (cid:105) ≈ N / σ [( u/v )( φ c − φ ∗ c ) /φ ∗ c ] − / , i.e., an equilibrium globule with (cid:104) R g (cid:105) ∼ N / σ gets further compacted with cross-linker densityas [( u/v )( φ c − φ ∗ c ) /φ ∗ c ] − / . The solution ρ = 0 at φ c < φ ∗ c corresponds to an open chain following Floryscaling R g ≈ σN / . Allowing for a bilinear coupling be-tween the monomer and the binder density fields leadsto ρ ∼ φ suggesting a non-linear decrease in (cid:104) R g (cid:105) = R g [1 + N / ( z Ω /u ) φ c ] − / (see Appendix-D).
2. Cluster size
An estimate of the increase in the cluster size of thepolymer-bound cross-linkers can be obtained by perform-ing linear stability analysis of Eq.(2) around a uniformstate of ρ = ¯ ρ and φ = ¯ φ . This analysis is presentedin detail in Appendix-E. It shows that the uniform stategets unstable towards formation of clusters as an effectivecoupling strength χ = u ¯ ρ/ Ω φ ∗ c crosses a threshold value.The mean spatial extension of such clusters is given by (cid:96) = 2 π M φ κwru (cid:16) − φ c φ ∗ c (cid:17) + 3 vu ¯ ρ / . (cid:104) C s (cid:105) ∼ (cid:96) leading to (cid:104) C s (cid:105) = A (cid:20)(cid:18) − φ c φ ∗ c (cid:19) + 3 vu ¯ ρ (cid:21) − / . (5)Replacing ¯ ρ = ( z Ω /u ) φ c , the dependence (cid:104) C s (cid:105) = A [(1 − φ c /φ ∗ c + B φ c )] − / reasonably captures the growth inmean cluster size with A = 1 . B = (3 vz Ω /u ) such that B φ c (cid:28)
1, as the coil-globuletransition is approached from below (Fig.2( b )).
3. Time scale
The diverging time-scale observed in simulations canbe understood using the following scaling argument basedon Eq.(2). For this purpose, we use the length scaleassociated with the unstable mode (cid:96) . Eq.(2) suggestsa relaxation time τ r ≈ ( (cid:96) /M ρ u )(1 − φ c /φ ∗ c ) − . Using γφ c (cid:28) τ r ≈ π M ρ u (cid:16) M φ κwru (cid:17) / (cid:20) − φ c φ ∗ c (cid:21) − / , (6)suggesting a divergence of correlation times as (1 − φ c /φ ∗ c ) − / near the critical point. For finite sized chains,while the time scales do not diverge, they show signifi-cant increase near criticality (Fig.2( c )). Added with aconstant background, Eq.6 gives a reasonable descrip-tion of the simulation results. As is shown in Fig.10( c ) ofAppendix-C, the correlation time at criticality increaseswith chain length with an approximate power law ∼ N / indicating divergence. IV. LOCAL MORPHOLOGY
The binder mediated chromosomal compaction is as-sociated with local morphological changes. The cross-linking due to binders may cause loop formation. Inchromosomes, formation of such loops are expected tobe highly complex, involving polydispersity of loop-sizes.The cross-linkers may also form zipper between contigu-ous polymeric segments. These, in turn, would enhancecontact formation, and as a result modify subchain ex-tensions. In this section, we discuss the change in all ofthese three aspects along the phase transition describedabove.
A. Loops
We describe the possible loop-topologies with the helpof Fig.4. A simply connected, or, first order loop isformed by a cross-linker binding two segments of thepolymer in such a way that if one moves along the chainfrom one such segment to the other, no other cross-linker lsdso = 1 :o = 2 :o = 3 :
FIG. 4: (color online) Schematics of loop topologies of order o :Polymer segments are indicated by blue beads and polymer-bound cross-linkers are shown by red open circles. ( a ) Simplyconnected loops of order o = 1. Two first order loops of size l s are separated by a gap of size d s . ( b ) Three examples of o = 2 loops. In the first two cases, the second order loop hasone o = 1 loop embedded inside. The third case shows twoembedded o = 1 loops. ( c ) Three examples of o = 3 loops.In the first two cases, the third order loop has a o = 1 and a o = 2 loop embedded. The third case shows two first orderloops and a second order loop embedded inside the o = 3loop. is encountered on the way. With removal of the cross-linker- bond stabilizing such a loop, the first order loopitself disappears (Fig.4: o = 1). In the figure, l s and d s denote loop-size and gap-size between such loops, respec-tively. In numerical evaluation of mean d s , all interme-diate higher order loops are disregarded.A higher order loop denoted by order o = n , embedsall possible lower order loops o = 1 , . . . , ( n −
1) withinit. In Fig.4: o = 2, three examples of second order loopsare shown. In the first two examples removing one cross-linker reduces the second order loop to a first order loop.In the third example of o = 2 loop, three bonds of a singlecross-linker maintains the loop, and with its removal thewhole loop structure disappears. In Fig.4: o = 3 we showthree examples of third order loops. Note that the firstorder and higher order loops identified here are related tothe serial and parallel topologies described in Ref. [57].As it has been shown before, consideration of chromo-somal loops is crucial in understanding of its emergentbehavior [98–101]. In this paper, we restrict ourselvesto the relative importance of different orders of loops inlocal chromosomal morphology.In Fig.5( a ) the mean number of loops (cid:104) n o (cid:105) of order o = 1 , , φ c . All through, (cid:104) n (cid:105) remains larger than (cid:104) n o =2 , (cid:105) cor-responding to higher order loops that show a sigmoidaldependence on φ c . Interestingly, (cid:104) n (cid:105) maximizes at thephase transition point φ ∗ c . Thus at the critical point the6 . . . . a ) h n o i φ r h n ih n ih n i − − − − − b ) Π o oφ r = 0 . φ r = 1 . φ r = 1 . φ r = 2 . φ r = 2 . φ r = 3 . FIG. 5: (color online) ( a ) Mean number of o -th order loops (cid:104) n o (cid:105) as a function of density of cross-linkers φ r = φ c × .Here (cid:104) n , , (cid:105) denote the mean number of first, second andthird order loops. ( b ) Probability of o -th order loop Π o isplotted on semi-log scale, for various cross-linker densitiesdenoted in the labels. At φ r = 3.14, the probability ofhigher order loops decays with an approximate Gaussian formexp( − o / g ) where the standard deviation g = 9 .
83 (solidbrown line). For φ r = 1.04, the probability of higher or-der loops decays exponentially as exp( − o/ ¯ o ) with ¯ o = 1 . local morphology of the model chromosome is dominatedby the first order loops.Fig.5( b ) shows the probability Π o of a loop to be of o -th order. At small cross-linker densities φ c < φ ∗ c , theprobability of higher order loops fall exponentially asΠ o = exp( − o/ ¯ o ). This behavior changes qualitativelyafter the coil-globule transition ( φ r = 1 .
57) to a Gaus-sian profile exp( − o / g ), as is shown in Fig.5( b ) .Given that loop sizes could be measured from elec-tron microscopy [44], we further analyze the statisticsof loop-sizes and inter-loop gaps corresponding to thefirst order loops in Fig. 6. With increasing cross-linkerdensity φ c , the mean size of first order loops (cid:104) l s (cid:105) de-creases (Fig. 6( a )), as their number increases (Fig.5( a ))reducing the mean gap size (cid:104) d s (cid:105) (Fig.6( b )). However, in-creased φ c stabilizes the loops, shown by decreased fluc-tuation of loop-sizes δl s = (cid:112) (cid:104) l s (cid:105) − (cid:104) l s (cid:105) . The mean gapsize (cid:104) d s (cid:105) and its fluctuation δd s = (cid:112) (cid:104) d s (cid:105) − (cid:104) d s (cid:105) reachtheir minimum at the transition point φ ∗ c = 1 . × − .The increase in the inter-loop separation (cid:104) d s (cid:105) beyond thispoint is due to the increase in probability of higher order
024 0 1 2 3( a ) 10203040 0 1 2 3( b )10 − − − − − c ) 10 − − − − − d ) h l s i , δ l s φ r h l s i δl s h d s i , δ d s φ r h d s i δd s P ( l s ) l s P ( d s ) d s FIG. 6: (color online) ( a ) Decrease of the mean size of firstorder loops (cid:104) l s (cid:105) and its fluctuations δl s with φ r = φ c × .( b ) Non-monotonic variation of mean separation between firstorder loops (cid:104) d s (cid:105) and its fluctuations δd s with φ r . Prob-ability distributions of the size of first order loops P ( l s )and gaps between them P ( d s ) are plotted in ( c ) and ( d )at φ r = 0 .
26 ( (cid:50) ) , .
57 ( ◦ ) , π ( (cid:52) ). At the transition point φ r = 1 . P ( l s ) ∼ l − . s , and P ( d s ) ∼ exp( − d s /λ ) with λ = 13 . σ , shown by the solid (brown) lines in ( c ) and ( d )respectively. loops in the local morphology of the model chromatin.Fig.6( c ) and ( d ) show the probability distributions offirst order loop sizes P ( l s ), and separation between con-secutive first order loops P ( d s ), respectively. For all φ c values, P ( l s ) ∼ l − µs , with µ increasing with φ c in asigmoidal fashion, giving µ = 3 . φ ∗ c = 1 . × − . The power law distribution of P ( l s )shows that their is no characteristic loop size, and loopsof all possible lengths are present. On the other hand, thegap size distributions follow an approximate exponentialform P ( d s ) ≈ (1 / (cid:104) d s (cid:105) ) exp( − d s / (cid:104) d s (cid:105) ). B. Zippering
The binders can also zipper different segments of thepolymer. The inset of Fig.7 shows one such zipper main-tained by cross-linkers. The zipper fraction of a con-formation is given by Z p = (1 /N ) (cid:80) ξ,i N ξi , where N ξi are the number of monomers involved in forming ξ -thzipper, and N is the total number of monomers in thechain. Fig.7 shows variation of ensemble averaged zip-pered fraction with the cross-linker density. The zipperfraction increases non-linearly to saturate in the equilib-rium globule phase to a value that remains within 60%of the completely zippered filament (cid:104) Z p (cid:105) = 1. Near thecritical point of the coil-globule transition (cid:104) Z p (cid:105) ≈ . FIG. 7: (color online) The zippered fraction of model chro-mosome (cid:104) Z p (cid:105) as a function of the cross-linker density φ r = φ c × . The inset shows two contiguous segments contain-ing N and N monomers (blue beads) forming a zipper viabinders (red beads). The corresponding zipper fraction is Z p = ( N + N ) /N . C. Contacts and subchain extensions
Conformational relaxation of a polymer brings contourwise distant parts of the chromatin in contact with eachother, even in absence of binders. An example of sucha contact in our model was shown in Fig.1( c ). The pro-cesses of cross-linking and zippering increase the prob-ability of contact formation that shows an asymptoticbehavior Π c ( s ) ∼ s − α between segments separated bya contour length s . With increasing cross-linker den-sity φ c , the exponent α decreases non-linearly to vanish,capturing the asymptotic plateauing of Π c ( s ) at large φ (Appendix-F). At criticality, α ≈ .
1, similar to thefractal globule and human chromosomes [45, 63, 64].This behavior is in agreement with an earlier latticemodel [28]. The mean subchain extension shows asymp-totic power law (cid:104) r ( s ) (cid:105) ∼ s ν . The exponent ν reducesfrom 3 / φ c = 0 to the fractal-globule like ν ≈ / ν = 0, a behavior typical of equilibrium glob-ules (Appendix-G). Finally, the structure of the detailedcontact-map changes with φ (Appendix-H). At small φ c ,only contour-wise neighbors participate in contact forma-tion. However, at the transition point φ c = φ ∗ c , contactsbegin to percolate to monomers separated by long con-tour lengths. V. DISCUSSION
In summary, using an off-lattice model of self avoid-ing polymer and diffusing protein binders cross-linkingdifferent segments of the chromatin fibre, we have pre-sented an extensive characterization of the continuous chromatin folding transition, and analyzed the associatedchanges in chromatin morphology in terms of formationof loops, zippering and contacts. The criticality is charac-terized by unimodal distributions, divergent fluctuationsand critical slowing down. The negative maximum in thecross-correlation between the number of attached bindersand chromosome size, at criticality, might be utilized byliving cells for easy switching between folded and openconformations, providing easy access to DNA-trackingenzymes. This is suggestive of a possibility that chro-mosomes might be poised at criticality [102], vindicatedfurther by the similarity of the calculated contact prob-ability at the critical point with the average behaviorof human chromosomes. Although the local chromatinmorphology does show highly complex loop structures,at criticality, it is dominated by simply connected loops.Each coarse-grained chromatin bead in our model canbe considered as 10 −
12 closely packed nucleosomes con-taining around 2 − . σ ≈ −
40 nm [63, 103]. The dimensionlesscritical volume fraction φ c is equivalent to a concentra-tion [ φ ∗ c / (4 πσ / ∼
60 nmol/l −
470 nmol/l. The mean size of the first or-der loops observed at criticality translates to 4 − (cid:104) l s (cid:105) : (cid:104) d s (cid:105) ≈ η ,the dissipation constant γ = 3 πησ . As it has been ob-served, the nucleoplasm viscosity η felt by objects withinthe cell nucleus depends on their size [104, 105]. Usingthe measured viscosity ∼
10 Pa-s felt by solutes having ∼
10 nm size [105] for the σ = 20 nm beads, the char-acteristic time which is the same as the time requiredto diffuse over the length-scale σ can be determined byusing the relation τ = γσ /k B T = 0 . τ R g denoting chromosomal relax-ation over ∼ . − . ≈
22 minutes atthe critical point.Here we should reemphasize that our study representsan average description of chromosomes using a coarsegrained homopolymer model. This approach did not aimto distinguish interaction between specific protein typesand gene sequences. While some of our predictions ap-pear to compare well with experiments, others involvingcross-linker clusters, relaxation time, and loop morphol-ogy are amenable to experimental verifications.
Acknowledgments
The simulations were performed on SAMKHYA,the high performance computing facility at IOP,Bhubaneswar. DC acknowledges SERB, India, for finan-cial support through grant number EMR/2016/001454,ICTS-TIFR, Bangalore for an associateship, and thanksBela M. Mulder for useful discussions. We thank Ab-8ishek Chaudhuri for a critical reading of the manuscript.
Appendix A: Valency of cross-linkers . . .
27 0 4 8 12 16 Π v n v φ r = 0 . φ r = 1 . φ r = 3 . FIG. 8: (color online) Π v represents probability that a cross-linker is simultaneously attached to n v number of monomers,at cross-linker densities φ r = φ c × denoted in the figure.The maximum shifts from n v = 4 to 6 across the coil-globuletransition. At the transition point, φ r = 1 .
57, the maximumprobability corresponds to the valency n v = 5. The cross-linkers used in the simulations are poten-tially multivalent. Here the question we ask is how manymonomers of the chain, a cross-linker attaches to simulta-neously? From the simulations, we identify the polymerbound cross-linkers, count the number of monomers thatlie within the range of attraction r c = 1 . σ identifyingthe instantaneous valency of a cross-linker, and computethe histogram over all the cross-linkers and time. Thisleads to the probability Π v of valency v , normalized to (cid:80) v Π v = 1. The maximum of the probability indicatesthe typical valency of cross-linkers at an ambient density φ c (see Fig.8). As φ c increases, the peak shifts towardslarger values. It means with increase of φ c , a single cross-linker on an average binds to more number of monomersof the chain. In the fully compact state, at the largest φ c , the typical valency we find is 6. Appendix B: Correlation function
In Fig.9, we present normalized auto-correlation func-tions of the polymer radius of gyration C R g ( t ) = (cid:104) δR g ( t ) δR g (0) (cid:105) / (cid:104) δR g (cid:105) , and the total number of boundcross-linkers C n a ( t ) = (cid:104) δn a ( t ) δn a (0) (cid:105) / (cid:104) δn a (cid:105) at variouscross-linker densities φ c across the coil-globule transi-tion. The correlations show approximate exponential de-cay exp( − t/τ c ) with correlation time τ c denoted by τ R g for the polymer radius of gyration, and τ n a for the totalnumber of polymer bound cross-linkers. The fitted val-ues show a maximum at the phase transition point φ ∗ c as − . . . .
751 0 10000 20000 30000 40000( a ) C R g ( t ) t φ r = 0 . φ r = 1 . φ r = 3 . − . . . .
751 0 10000 20000 30000 40000( b ) C n a ( t ) t φ r = 0 . φ r = 1 . φ r = 3 . FIG. 9: (color online) The auto-correlation functions( a ) C R g ( t ) of polymer radius of gyration R g , and ( b ) C n a ( t ) ofthe number of cross-linkers attached to the chain n a , at threecross-linker densities φ r = φ c × ; t is expressed in unit of τ . Fitting them to exponential forms exp( − t/τ c ) gives cor-relation time τ c = τ R g , τ n a for R g and n a respectively. Twosuch fittings are shown in each plot by solid lines. The fittedcorrelation times are τ R g = 7371 τ , τ n a = 5402 τ at φ r = 1 . τ R g = 307 τ , τ n a = 720 τ at φ r = π . is shown in Fig. 2( c ). Appendix C: System size dependence at thecoil-globule transition
We performed simulations with various chain lengths N to check the system size dependence on the coil-globuletransition. We observe a continuous change of (cid:104) R g (cid:105) withcross-linker density φ c = φ r × − (Fig. 10( a ) ) with thetransition becoming sharper at larger N . The relativefluctuations ∆ R g /R g at the transition point increaseswith N (Fig.10( b ) ). The correlation time τ R g character-izing fluctuations in R g at the transition point divergeswith the increase in polymer size N . The simulationssuggest a dynamical scaling τ R g ∼ N ζ with ζ ≈ / Appendix D: Before transition
The most general mean field theory would also con-tain terms bilinear in ρ and φ , which we neglected inthe discussion of phase transition. In presence of such acoupling, keeping terms only up to quadratic order, βf = 12 u ρ + 12 wφ − zρφ + κ ∇ ρ ) . . . . .
81 0 1 . . a ) 00 . . . . b )10
50 300 10 ( c ) h R g i / R g φ r ∆ R g / h R g i φ r τ R g N FIG. 10: (color online) Variation of ( a ) the scaled radius ofgyration (cid:104) R g (cid:105) /R g , and ( b ) its relative fluctuations ∆ R g / (cid:104) R g (cid:105) with cross-linker density φ r = φ c × for chain lengths N =128 ( (cid:13) ) ,
256 ( (cid:53) ) ,
512 ( (cid:5) ) , (cid:52) ). ( c ) The correlation timeat the transition point τ R g increases with chain length N .The dash-dotted line denotes τ R g ∼ N / . This does not describe any phase transition, however,suggests a uniform mean field solution ρ = zφ /u = z Ω φ c /u . Thus, before the transition, the mean radiusof gyration is expected to decrease with φ c as (cid:104) R g (cid:105) = R g [1 + N / ρ ] − / = R g [1 + N / ( z Ω /u ) φ c ] − / . Appendix E: Linear stability analysis
The formation of cross-linker clusters mediates foldingof the model chromosome. To characterize the dynamics,we use linear stability analysis for small deviations froma homogeneous state ρ = ¯ ρ + δρ ( r ), φ = ¯ φ + δφ ( r ). Thedynamics in Eq.(2) for these small deviations become ∂ t δρ = D ρ ∇ δρ − M ρ κ ∇ δρ − M ρ χ ∇ δφ∂ t δφ = D φ ∇ δφ − M φ χ ∇ δρ − r δφ, where D ρ = M ρ u (cid:20)(cid:18) − ¯ φφ ∗ (cid:19) + 3 vu ¯ ρ (cid:21) , and D φ = M φ w are the effective diffusion constants ofthe two components, and χ = u ¯ ρ/φ ∗ is the strength ofcross-coupling. In the above equations ∂ t denotes thepartial derivative with respect to time t . Expressing timein units of inverse turnover rate, τ u = 1 /r , and lengthsin units of x u = (cid:112) M φ w/r , one finds ∂ τ δρ = D ∇ ξ δρ − K∇ ξ δρ − C∇ ξ δφ∂ τ δφ = ∇ ξ δφ − C (cid:48) ∇ ξ δρ − δφ, (E1)with control parameters of the dynamics D = D ρ /M φ w , K = M ρ M φ κrw , C = M ρ M φ χw , and C (cid:48) = χw . The dimensionless time and length scales are denoted by τ = t/τ u , and ξ = x/x u , respectively.Fourier transform of this equation gives evolution ofmodes as matrix equations ∂ τ ( δρ q , δφ q ) = M ( δρ q , δφ q ),where, M = (cid:18) − q ( D + K q ) C q C (cid:48) q − ( q + 1) (cid:19) . The eigenvalues of M are given by λ ( q ) = 12 (cid:110) Tr M ± (cid:112) (Tr M ) − M (cid:111) As the trace of this matrixTr . M = − q ( D + K q ) − ( q + 1) < , the only way of having instability (one of the eigenvaluesbecomes positive) is if the determinantdet M = q ( q + 1)( D + K q ) − CC (cid:48) q < . This last criterion leads to CC (cid:48) > F ( q ) where F ( q ) =(1 + q )( D + K q ). This will be satisfied for any q ifeven the minimum of F ( q ) obeys this inequality. Onecan easily show that F ( q ) is minimized at q = (cid:112) D / K ,and F ( q ) = ( √D + √K ) . Thus the instability criterionbecomes √CC (cid:48) > ( (cid:112) D + √K ) . Note from Eq.(E1) that C and C (cid:48) denote coupling coef-ficients between the evolution of the two fields ρ and φ .The inequality suggests a minimal coupling strength χ is required to generate instability towards formation ofcross-linker clusters, χ > (cid:114) κrM φ + (cid:115) uw (cid:20)(cid:18) − ¯ φφ ∗ (cid:19) + 3 vu ¯ ρ (cid:21) . Once this condition is satisfied, instability in the formof clustering of cross-linkers, mediated by the attractiveinteraction with monomers, arise. The fastest growingmode q = ( D / K ) / predicts the most unstable lengthscale (cid:96) /x u = 2 π/q = 2 π ( K / D ) / , which gives themean extension of the clusters (cid:96) = 2 π M φ κwru (cid:16) − ¯ φφ ∗ (cid:17) + 3 vu ¯ ρ / . Appendix F: Contact probability
To analyze contact formation from simulations one re-quires a finite cutoff length such that if two monomersfall within such a separation they are defined to be incontact. Here we use r c = 1 . σ for this purpose. We10 − − − − − − a ) Π c ( s ) sφ r = 0 φ r = 1 . φ r = 3 . s − . s − . . . . . . . . b ) α φ r FIG. 11: (color online) ( a ) Contact probabilities Π c ( s ) atdifferent cross-linker densities φ r = φ c × . They followasymptotic power law profiles Π c ( s ) ∼ s − α at all φ c , with α being a function of φ c . ( b ) The decrease of exponent α withincreasing φ r is related to the coil-globule transition. have checked that our main results do not depend on theprecise choice of this length scale.The contour wise separation between two monomers, s , defines the genomic distance between chromatin seg-ments. The contact probability Π c ( s ) is a measure ofsuch two segments to be in contact. In absence of binders,we get Π c ( s ) ∼ s − α with α ≈ .
1, as expected for self-avoiding chains [59]. Even in presence of cross-linkers,the asymptotic power law persists with φ c -dependent α (Fig.11.( a )). At the transition point φ ∗ c = 1 . × − ,the simulation results are consistent with α ≈ .
1, anumber that agrees well with the prediction of the frac-tal globule model [63, 64]. It is interesting to note that α ≈ . φ c values, after the comple-tion of the coil-globule transition, contact probabilitiesat large s plateaus to a constant, indicating α = 0. As afunction of φ c , the asymptotic exponent α reveals a con-tinuous decrease (see Fig.11( b ) ), capturing the change inpolymeric organization in the course of the coil-globuletransition. Appendix G: Extension of subchains
Here we consider the scaling behavior of subchain ex-tensions, measured in terms of the mean squared end to h r ( s ) i sφ r = 0 . φ r = 1 . φ r = 3 . s / s / FIG. 12: (color online) The scaling behavior of subchainextension (cid:104) r ( s ) (cid:105) , at three cross-linker concentrations φ r = φ c × , before, at and after the coil-globule transition. Atlow densities, it approximately follows Flory scaling (cid:104) r ( s ) (cid:105) ∼ s / (red solid line). At the transition point φ c = 1 . × − ,the asymptotic behavior agrees with the fractal globule esti-mate (cid:104) r ( s ) (cid:105) ∼ s / (blue dashed line). At the highest con-centrations we find asymptotic plateauing, a characteristic ofequilibrium globules. end distance (cid:104) r ( s ) (cid:105) in subchains of contour length s .We observe three different scaling behaviors across thecoil-globule transition.A sub-chain inside a compact equilibrium globule is ex-pected to behave like a random walk due to strong screen-ing of interaction by large monomeric density. Thus (cid:104) r ( s ) (cid:105) ∼ s , before the globule boundary is encoun-tered. Multiple reflections from the globule boundary, as s > (cid:104) r ( s ) (cid:105) ∼ N / , fills the space inside the globule uni-formly, so that it becomes equally likely to find the otherend of the subchain anywhere inside the globule, satu-rating (cid:104) r ( s ) (cid:105) to a constant. Thus in equilibrium glob-ules (cid:104) r ( s ) (cid:105) ∼ s up to s < N / , and saturates beyondthat length scale [59, 63]. The random loop model, withfixed probability of attraction between monomers, showsall the features of equilibrium globule in final configura-tions [106, 107]. On the other hand, the fractal globule isspace filling at all scales, such that (cid:104) r ( s ) (cid:105) ∼ s / [63, 64].At small φ c (= 0 . × − ), we find a behavior typicalof open chains, (cid:104) r ( s ) (cid:105) ∼ s / , that follows Flory scaling.In the fully folded compact phase at high φ c (= π × − ), (cid:104) r ( s ) (cid:105) shows plateauing at large s as in compact equilib-rium globules, and random loop models [63, 85, 106, 107].Such plateauing was earlier related to folding of chromo-some into territories [108]. In the compact phase, themolecular cross-linkers may not only pull different seg-ments close to each other, by doing so, they may displacewell separated parts further away from each other [29], re-flected in the eventual increase of (cid:104) r ( s ) (cid:105) as s approachesthe full length N , e.g., at highest φ c . At the critical point, φ ∗ c (= 1 . × − ), simulation results for subchain ex-11 IG. 13: (color online) From left to right, contact maps at φ c = 1 . × − , . × − , and 1 . × − are plotted. Thecolor code captures contact frequency and is shown in log scale. tensions is consistent with (cid:104) r ( s ) (cid:105) ∼ s / as in fractalglobules [64]. This is close to the threshold-exponentpredicted in [28] (cid:104) r ( s ) (cid:105) ∼ s ν with ν = 0 .
39. Thus withincreasing cross-linker density, the model chromatin mor-phology changes from an open chain to compact equilib-rium globule, via an intermediate fractal globule behaviorobserved at the critical point. The sustenance of fractalglobule like non-equilibrium behavior at the critical pointcan be understood in terms of the super-slow relaxation.
Appendix H: Contact maps
In Fig.13 we present ensemble averaged contact mapsover equilibrium configurations at different cross-linker densities. Such maps represent probability measures oftwo chromatin segments to be in spatial proximity. Atthe coil-globule transition φ c = φ ∗ c , the contact showsemergence of local pattern, indicating enhanced proba-bility of contour-wise well separated segments to comeinto spatial proximity. In the compact phase at φ c =1 . × − , the chromosomal contacts spread over thewhole chromatin chain. [1] B. Alberts, D. Bray, K. Hopkin, A. Johnson, J. Lewis,M. Raff, K. Roberts, and P. Walter, Essential Cell Bi-ology , 4th ed. (Garland Science, New York, USA, 2015).[2] A. Wolffe,
Chromatin: Structure and function , 3rd ed.(Academic press, London, 2012).[3] M. Thanbichler, P. H. Viollier, and L. Shapiro, Curr.Opin. Genet. Dev. , 153 (2005).[4] W. Wang, G.-W. Li, C. Chen, X. S. Xie, and X. Zhuang,Science , 1445 (2011).[5] K. Luger, A. W. M¨ader, R. K. Richmond, D. F. Sargent,and T. J. Richmond, Nature , 251 (1997).[6] R. T. Dame, M. C. Noom, and G. J. Wuite, Nature ,387 (2006).[7] P. A. Wiggins, R. T. Dame, M. C. Noom, and G. J.Wuite, Biophys. J. , 1997 (2009).[8] A. Worcel and E. Burgi, J. Mol. Biol. 71: , 127 (1972).[9] V. Holmes and N. Cozzarelli, Proc. Natl. Acad. Sci. ,1322 (2000).[10] S. B. Zimmerman, Shape and compaction of Escherichiacoli nucleoids, 2006.[11] R. A. Van Der Valk, J. Vreede, F. Cr´emazy, and R. T.Dame, J. Mol. Microbiol. Biotechnol. , 344 (2014).[12] S. S. Rao, M. H. Huntley, N. C. Durand, E. K. Stamen-ova, I. D. Bochkov, J. T. Robinson, A. L. Sanborn, I.Machol, A. D. Omer, E. S. Lander, and E. L. Aiden,Cell , 1665 (2014).[13] R. T. Dame and M. Tark-Dame, Curr. Opin. Cell Biol. , 60 (2016).[14] M. S. Luijsterburg, M. C. Noom, G. J. Wuite, and R. T.Dame, J. Struct. Biol. , 262 (2006).[15] R. T. Dame, Mol. Microbiol. , 858 (2005).[16] E. Peeters, R. P. C. Driessen, F. Werner, and R. T.Dame, Nat. Rev. Microbiol. , 333 (2015).[17] W. A. Bickmore and B. van Steensel, Cell , 1270(2013).[18] K. Nasmyth and C. H. Haering, Annu. Rev. Genet. ,525 (2009).[19] J. E. Phillips and V. G. Corces, Cell , 1194 (2009).[20] C. Annunziatella, A. M. Chiariello, A. Esposito, S.Bianco, L. Fiorillo, and M. Nicodemi, Methods ,81 (2018).[21] G. Le Treut, F. K´ep`es, and H. Orland, Biophys. J. ,51 (2016).[22] C. A. Brackley, B. Liebchen, D. Michieletto, F. Mouvet,P. R. Cook, and D. Marenduzzo, Biophys. J. , 1085(2017).[23] C. A. Brackley, S. Taylor, A. Papantonis, P. R. Cook,and D. Marenduzzo, Proc. Natl. Acad. Sci. U. S. A. ,E3605 (2013).[24] J. Johnson, C. A. Brackley, P. R. Cook, and D. Maren-duzzo, J. Phys. Condens. Matter , 064119 (2015).[25] C. A. Brackley, J. Johnson, S. Kelly, P. R. Cook, andD. Marenduzzo, Nucleic Acids Res. , 3503 (2016).[26] D. Marenduzzo, I. Faro-Trindade, and P. R. Cook, rends Genet. , 126 (2007).[27] A. B. Kolomeisky, Biophys. J. , 459 (2015).[28] M. Barbieri, M. Chotalia, J. Fraser, L.-M. Lavitas, J.Dostie, A. Pombo, and M. Nicodemi, Proc. Natl. Acad.Sci. U. S. A. , 16173 (2012).[29] M. Nicodemi and A. Prisco, Biophys. J. , 2168 (2009).[30] M. Nicodemi and A. Pombo, Curr. Opin. Cell Biol. ,90 (2014).[31] G. Fudenberg and L. A. Mirny, Curr. Opin. Genet. Dev. , 115 (2012).[32] G. Fudenberg, M. Imakaev, C. Lu, A. Goloborodko, N.Abdennur, and L. A. Mirny, Cell Rep. , 2038 (2016).[33] A. L. Sanborn, S. S. P. Rao, S.-C. Huang, N. C. Durand,M. H. Huntley, A. I. Jewett, I. D. Bochkov, D. Chinnap-pan, A. Cutkosky, J. Li, K. P. Geeting, A. Gnirke, A.Melnikov, D. McKenna, E. K. Stamenova, E. S. Lan-der, and E. L. Aiden, Proc. Natl. Acad. Sci. , E6456(2015).[34] E. Alipour and J. F. Marko, Nucleic Acids Res. ,11202 (2012).[35] C. A. Brackley, J. Johnson, D. Michieletto, A. N. Mo-rozov, M. Nicodemi, P. R. Cook, and D. Marenduzzo,Phys. Rev. Lett. , 138101 (2017).[36] D. J. Brocken, M. Tark-Dame, and R. T. Dame, Curr.Opin. Syst. Biol. , 137 (2018).[37] A. Goloborodko, M. V. Imakaev, J. F. Marko, and L.Mirny, Elife , 1 (2016).[38] D. Song and J. J. Loparo, Trends Genet. , 164 (2015).[39] R. T. Dame, Mol. Microbiol. , 858 (2005).[40] R. T. Dame, in Bact. Chromatin , edited by R. T. Dame(Springer Netherlands, Dordrecht, 2018).[41] H. Delius and A. Worcel, Cold Spring Harb. Symp.Quant. Biol. , 53 (1974).[42] R. Kavenoff and B. Bowen, Chromosoma , 89 (1976).[43] N. Trun and J. Marko, ASM News , 276 (1998).[44] L. Postow, C. D. Hardy, J. Arsuaga, and N. R. Coz-zarelli, Genes Dev. , 1766 (2004).[45] E. Lieberman-Aiden, N. L. van Berkum, L. Williams,M. Imakaev, T. Ragoczy, A. Telling, I. Amit, B. R. La-joie, P. J. Sabo, M. O. Dorschner, R. Sandstrom, B.Bernstein, M. A. Bender, M. Groudine, A. Gnirke, J.Stamatoyannopoulos, L. A. Mirny, E. S. Lander, and J.Dekker, Science (80-. ). , 289 (2009).[46] T. B. K. Le, M. V. Imakaev, L. A. Mirny, and M. T.Laub, Science (80-. ). , 731 (2013).[47] M. Marbouty, A. Cournac, J.-F. Flot, H. Marie-Nelly,J. Mozziconacci, and R. Koszul, Elife , e03318 (2014).[48] T. Cremer and C. Cremer, Nat. Rev. Genet. , 292(2001).[49] J. Dekker, Science (80-. ). , 1793 (2008).[50] S. C. Dillon and C. J. Dorman, Nat. Rev. Microbiol. ,185 (2010).[51] S. Meyer, S. Reverchon, W. Nasser, and G. Muskhel-ishvili, Curr. Genet. , 555 (2018).[52] M. Fritsche, S. Li, D. W. Heermann, and P. A. Wiggins,Nucleic Acids Res. , 972 (2012).[53] P. M. Llopis, A. F. Jackson, O. Sliusarenko, I. Surovt-sev, J. Heinritz, T. Emonet, and C. Jacobs-Wagner, Na-ture , 77 (2010).[54] K. Mizuguchi and N. Go, Curr. Opin. Struct. Biol. ,377 (1995).[55] M. H. Bailor, X. Sun, and H. M. Al-Hashimi, Science(80-. ). , 202 (2010).[56] G. Cavalli and T. Misteli, Nat. Struct. Mol. Biol. , 290 (2013).[57] A. Mashaghi, R. J. van Wijk, and S. J. Tans, Structure , 1227 (2014).[58] C. M¨unkel, R. Eils, S. Dietzel, D. Zink, C. Mehring, G.Wedemann, T. Cremer, and J. Langowski, J. Mol. Biol. , 1053 (1999).[59] P.-G. de Gennes, Scaling concepts in polymer physics (Cornell University Press,, Ithaca, NY, 1979).[60] M. Rubinstein and R. H. Colbi,
Polymer Physics (Ox-ford University Press,, Oxford, NY, 2003).[61] D. Jost, A. Rosa, C. Vaillant, and R. Everaers,
NuclearArchitecture and Dynamics (Elsevier, London, 2018),Vol. 6, pp. 149–169.[62] A. Rosa and R. Everaers, PLoS Comput. Biol. ,e1000153 (2008).[63] L. A. Mirny, Chromosom. Res. , 37 (2011).[64] A. Y. Grosberg, S. Nechaev, and E. Shakhnovich, J.Phys. , 2095 (1988).[65] A. Lesage, V. Dahirel, J.-M. Victor, and M. Barbi,bioRxiv (2018).[66] J. Victor, J. Imbert, and D. Lhuillier, J. Chem. Phys. , 5372 (1994).[67] J. R. Dixon, S. Selvaraj, F. Yue, A. Kim, Y. Li, Y. Shen,M. Hu, J. S. Liu, and B. Ren, Nature , 376 (2012).[68] N. Ganai, S. Sengupta, and G. I. Menon, Nucleic AcidsRes. , 4145 (2014).[69] D. Jost, P. Carrivain, G. Cavalli, and C. Vaillant, Nu-cleic Acids Res. , 9553 (2014).[70] J. D. Olarte-Plata, N. Haddad, C. Vaillant, and D. Jost,Phys. Biol. , 26001 (2016).[71] M. Tark-Dame, H. Jerabek, E. M. M. Manders, D. W.Heermann, and R. van Driel, PLoS Comput. Biol. ,e1003877 (2014).[72] B. Zhang and P. G. Wolynes, Proc. Natl. Acad. Sci. , 6062 (2015).[73] A. M. Chiariello, C. Annunziatella, S. Bianco, A. Es-posito, and M. Nicodemi, Sci. Rep. , 29775 (2016).[74] T. Agarwal, G. P. Manjunath, F. Habib, P. L. Vad-davalli, and A. Chatterji, Journal of Physics: Con-densed Matter , 034003 (2017).[75] T. Agarwal, G. P. Manjunath, F. Habib, and A. Chat-terji, EPL (Europhysics Letters) , 18004 (2018).[76] T. Agarwal, G. P. Manjunath, F. Habib, and A. Chat-terji, arXiv:1808.09396 (2018).[77] T. Agarwal, G. P. Manjunath, F. Habib, and A. Chat-terji, arXiv:1810.05869 (2018).[78] A. Agrawal, N. Ganai, S. Sengupta, and G. I. Menon,bioRxiv: 10.1101/315812 (2018).[79] N. Ramakrishnan, K. Gowrishankar, L. Kuttippurathu,P. B. S. Kumar, and M. Rao, arXiv:1510.04157 (2015).[80] M. Tark-Dame, M. S. Luijsterburg, D. W. Heermann,and R. van Driel, in Genome Organization and Func-tion in the Nucleus , edited by K. Rippe (Wiley, Win-heim, Germany, 2012), Chap. Understanding genomefunction: Quantitative modelling of chromatin foldingand chromatin-associated processes, pp. 535– 555.[81] V. F. Scolari and M. Cosentino Lagomarsino, Soft Mat-ter , 1677 (2015).[82] V. F. Scolari, G. Mercy, R. Koszul, A. Lesne, and J.Mozziconacci, Phys. Rev. Lett. , 057801 (2018).[83] O. B. Ptitsyn and Y. E. Eisner, Biofizika , 3 (1965).[84] Y. E. Eisner, Polym. Sci. U.S.S.R. , 409 (1969).[85] P. G. de Gennes, J. Phys. (Paris) Lett. , 55 (1975).[86] I. M. Lifshitz, A. Y. Grosberg, and A. R. Khokhlov, ev. Mod. Phys. , 683 (1978).[87] I. M. Lifshitz, A. Y. Grosberg, and A. R. Khokhlov, Zh.Eksp. Teor. Fiz. , 1634 (1976).[88] Y. A. Budkov, A. L. Kolesnikov, N. Georgi, and M. G.Kiselev, J. Chem. Phys. , 014902 (2014).[89] J. D. Weeks, D. Chandler, and H. C. Andersen, J.Chem. Phys. , 5237 (1971).[90] G. Grest and K. Kremer, Phys. Rev. A , 3628 (1986).[91] J. D. Watson, T. A. Baker, S. P. Bell, A. Gann, M.Levine, and R. Losick, Molecular Biology of the Gene ,7th ed. (Benjamin-Cummings, Menlo Park, CA., 2003).[92] D. Frenkel, B. Smit,
Understanding molecular simula-tion: from algorithms to applications , (Academic press,NY, 2002).[93] H. Limbach, A. Arnold, B. Mann, and C. Holm, Com-put. Phys. Commun. , 704 (2006).[94] B. Duenweg, M. Stevens, and K. Kremer, in MonteCarlo and Molecular Dynamics Simulations in PolymerScience, edited by K. Binder (Oxford University Press,New York, 1995).[95] W. Humphrey, A. Dalke, and K. Schulten, Journal ofMolecular Graphics , 33 (1996).[96] M. P. Allen and D. J. Tildesley, Computer simulationof liquids (Oxford University Press, Oxford, 1989).[97] P. M. Chaikin and T. C. Lubensky,
Principles of Con-densed Matter Physics (Cambridge University Press,Cambridge, 2012).[98] D. Chaudhuri and B. M. Mulder, Physical Review Let-ters , 268305 (2012). [99] D. Chaudhuri and B. M. Mulder, in
Bacterial Chro-matin , edited by R. T. Dame (Springer Protocols, Hu-mana Press, New York, 2018), Chap. Molecular Dynam-ics Simulation of a Feather-Boa Model of a BacterialChromosome, p. 403.[100] P. Swain, B. M. Mulder, and D. Chaudhuri, Soft Matter , 2677 (2019).[101] F. Wu, A. Japaridze, X. Zheng, J. W. Kerssemakers,and C. Dekker, bioRxiv 10.1101/348052 (2018).[102] T. Mora and W. Bialek, J. Stat. Phys. , 268 (2011).[103] A. Routh, S. Sandin, and D. Rhodes, Proc. Natl. Acad.Sci. , 8872 (2008).[104] G. L. Lukacs, P. Haggie, O. Seksek, D. Lechardeur, N.Freedman, and A. S. Verkman, J. Biol. Chem. , 1625(2000).[105] Y. Tseng, J. S. H. Lee, T. P. Kole, I. Jiang, and D.Wirtz, J. Cell Sci. , 2159 (2004).[106] M. Bohn, D. W. Heermann, and R. van Driel, Phys.Rev. E , 051805 (2007).[107] J. Mateos-Langerak, M. Bohn, W. de Leeuw, O. Giro-mus, E. M. M. Manders, P. J. Verschure, M. H. G. In-demans, H. J. Gierman, D. W. Heermann, R. van Driel,and S. Goetze, Proc. Natl. Acad. Sci. , 3812 (2009).[108] T. Cremer and C. Cremer, Nat. Rev. Genet. , 292(2001).[109] In the analysis of phase transition, near the transitionpoint, a possible bilinear coupling − zρφ is neglected infavor of algebraic simplicity.is neglected infavor of algebraic simplicity.