OpenNucleome for high resolution nuclear structural and dynamical modeling

Zhuohan Lao; Kartik Kamat; Zhongling Jiang; Bin Zhang

doi:10.7554/eLife.93223.2

eLife assessment

This important work significantly advances the field of computational modelling of genome organisation through the development of OpenNucleome. The evidence supporting the tool's effectiveness is compelling, as the authors compare their predictions with experimental data. It is anticipated that OpenNucleome will attract significant interest from the biophysics and genomics communities.

https://doi.org/10.7554/eLife.93223.2.sa3

Significance of findings

important: Findings that have theoretical or practical implications beyond a single subfield

landmark
fundamental
important
valuable
useful

Strength of evidence

compelling: Evidence that features methods, data and analyses more rigorous than the current state-of-the-art

exceptional
compelling
convincing
solid
incomplete
inadequate

During the peer-review process the editor and reviewers write an eLife assessment that summarises the significance of the findings reported in the article (on a scale ranging from landmark to useful) and the strength of the evidence (on a scale ranging from exceptional to inadequate). Learn more about eLife assessments

Abstract

The intricate structural organization of the human nucleus is fundamental to cellular function and gene regulation. Recent advancements in experimental techniques, including high-throughput sequencing and microscopy, have provided valuable insights into nuclear organization. Computational modeling has played significant roles in interpreting experimental observations by reconstructing high-resolution structural ensembles and uncovering organization principles. However, the absence of standardized modeling tools poses challenges for furthering nuclear investigations. We present OpenNucleome—an open-source software designed for conducting GPU-accelerated molecular dynamics simulations of the human nucleus. OpenNucleome offers particle-based representations of chromosomes at a resolution of 100 KB, encompassing nuclear lamina, nucleoli, and speckles. This software furnishes highly accurate structural models of nuclear architecture, affording the means for dynamic simulations of condensate formation, fusion, and exploration of non-equilibrium effects. We applied OpenNucleome to uncover the mechanisms driving the emergence of “fixed points” within the nucleus—signifying genomic loci robustly anchored in proximity to specific nuclear bodies for functional purposes. This anchoring remains resilient even amidst significant fluctuations in chromosome radial positions and nuclear shapes within individual cells. Our findings lend support to a nuclear zoning model that elucidates genome functionality. We anticipate OpenNucleome to serve as a valuable tool for nuclear investigations, streamlining mechanistic explorations and enhancing the interpretation of experimental observations.

Introduction

The highly complex structural organization of the human nucleus plays a crucial role in the functioning and regulation of our cells.^1–10 The complexity arises from the diverse range of nuclear landmarks, such as nucleoli,¹¹ nuclear speckles,^8,12 and the nuclear lamina,¹³ each serving distinct functions. These landmarks provide specialized environments for various nuclear processes, allowing for efficient coordination and regulation of gene expression. Moreover, the spatial arrangement of chromosomes within the nucleus, intertwined with the nuclear landmarks, is critical for proper gene regulation and communication between different genome regions. Disruptions or abnormalities in the nuclear organization can have profound consequences on cellular function and can contribute to the development of diseases, including cancer and genetic disorders.^14,15

Recent advancements in experimental techniques have significantly enhanced our understanding of nuclear organization.^3,16–20 The advent of high-throughput sequencing-based methods, such as genome-wide chromosome-conformation capture (Hi-C), has unveiled crucial structural elements of the genome,^21,22 including chromatin loops,²³ topologically associating domains,^24,25 and compartments.²² Additionally, sequencing-based techniques such as DamID,²⁶ Chip-Seq,²⁷ and TSA-Seq²⁸ have revealed valuable information regarding interactions between chromosomes and nuclear landmarks. However, it is worth noting that these sequencing methods often offer averaged contacts, which can mask the heterogeneity present across populations, although single-cell techniques are also emerging.^29–31 Moreover, translating contact data into spatial positions can be challenging, adding complexity to interpreting experimental findings.

To complement these sequencing approaches, microscopic imaging techniques directly probe the spatial positions within individual nuclei.^3,13,32–34 Recent advancements in DNA FISH (Fluorescence In Situ Hybridization) have enabled high-throughput imaging of thousands of loci simultaneously.^35,36 These imaging studies have not only confirmed the structural features observed through sequencing techniques but have also provided valuable insights into the heterogeneity present at the single-cell level.

The abundance of available experimental data in the field of nuclear organization provides a fertile ground for structural modeling.^37–63 To make sense of this wealth of information, various computational approaches have been introduced, with polymer simulation approaches being extensively utilized. These simulation techniques aid in reconstructing structural ensembles that closely replicate experimental data, offering valuable insights into the mechanisms underlying chromosome folding. In recent studies, these approaches have also been employed to investigate the interplay between the genome and the nuclear lamina, ^64–67 as well as nucleoli,⁶⁸ shedding light on their dynamic relationships.

Despite the progress made in computational modeling, the absence of well-documented software with easy-to-follow tutorials pose a challenge. Many research groups develop their own independent software, which complicates cross-validation and hinders the establishment of best practices for genome modeling. ^40,69,70 Moreover, comprehensive models of the entire nucleus, especially at high resolution, remain scarce. Addressing these limitations and fostering collaboration in the scientific community can be achieved through the development of open-source tools. By promoting transparency and accessibility, such tools have the potential to greatly facilitate nuclear modeling and contribute to a more unified and collaborative research environment.

We present OpenNucleome, an open-source software designed for conducting molecular dynamics (MD) simulations of the human nucleus. This software streamlines the process of setting up whole nucleus simulations through just a few lines of Python scripting. Open-Nucleome can unveil intricate, high-resolution structural and dynamic chromosome arrangements at a 100 KB resolution. It empowers researchers to track the kinetics of condensate formation and fusion while also exploring the influence of chemical modifications on condensate stability. Furthermore, it facilitates the examination of nuclear envelope deformation’s impact on genome organization. The software’s modular architecture enhances its adaptability and extensibility. Leveraging the power of OpenMM,⁷¹ a GPU-accelerated MD engine, OpenNucleome ensures efficient simulations.

Our work demonstrates the fidelity of the simulated nuclear organizations by faithfully reproducing Hi-C, Lamin B DamID, TSA-Seq, and DNA-MERFISH data. The dynamic insights extracted from this model are pivotal in advancing our understanding of nuclear organization mechanisms. Our findings reveal that inherent heterogeneity in chromosome contacts naturally emerges within single cells. Interestingly, robust contacts between chromosomes and nuclear bodies can also be established due to a coupled self-assembly mechanism. Notably, the resilience of contacts involving nuclear bodies supports a nuclear zoning model for genome function. In the realm of nuclear investigations, we anticipate OpenNucleome to serve as an invaluable tool, seamlessly complementing experimental techniques.

Results

Non-equilibrium nucleus model at 100 KB resolution

We present an open-source implementation of a computational framework that facilitates the structural and dynamical characterization of the human nucleus. This framework builds upon a previous investigation⁶⁵ but incorporates several significant modifications. Firstly, we enhance the model resolution by a factor of ten, enabling the precise determination of the spatial positioning of each chromatin segment measuring 100 KB in length. Secondly, we present a kinetic scheme for speckles that accounts for the phosphorylation of protein molecules. This inclusion captures the influence of chemical reactions on the stability and dynamics of nuclear bodies. Thirdly, we incorporate explicit nuclear envelope dynamics to explore the impact of large-scale deformations on genome organization. Finally, our implementation into OpenMM offers the advantages of Python Scripting and GPU acceleration, facilitating easy extension and customization. These features will facilitate the broad applicability and adoption of the proposed model.

The nucleus model provides particle-based representations for chromosomes, nucleoli, speckles, and the nuclear envelope. As shown in Fig. 1A and B, each of the 46 chromosomes is represented as a beads-on-a-string polymer, where each bead represents a 100-KB-long genomic segment. Based on Hi-C data, we further assign each bead as compartment A, B, or C to signify euchromatin, heterochromatin, or pericentromeric regions. The lamina was modeled as a spherical enclosure with 10 µm diameter, using discrete particles arranged to represent a mesh grid with covalent bonds linking together nearest neighbors.⁷² We modeled nucleoli and speckles as liquid droplets that emerge through the spontaneous phase separation of coarse-grained particles, representing protein and RNA molecule aggregates.^8,11 These particles exhibited attractive interactions within the same type to promote condensation. More details about the various components of the system can be found in the Supporting Information Section: Components of the whole nucleus model.

Computer model of the human nucleus for structural and dynamical characterizations.
(A) 3D rendering of the nucleus model with particle-based representations for the 46 chromosomes shown as ribbons, the nuclear lamina (grey), nucleoli (cyan), and speckles (yellow). As shown on the right, chromosomes are modeled as beads-on-a-string polymers at a 100 KB resolution, with the beads further categorized into compartment A (red), compartment B (light blue), or centromeric regions (green). (B) Speckle particles undergo chemical modifications concurrent to their spatial dynamics, and the dephosphorylated (dP) particles contribute to droplet formation. (C) Illustration of the ideal and compartment potential that promotes chromosome compaction and microphase separation. Specific interactions between chromosomes and nuclear landmarks are shown on the right.

The energy function of the nucleus model includes three components that account for the self-assembly of chromosomes, the assembly of nuclear bodies, and the coupling between chromosomes and nuclear landmarks. Therefore, the model approximates nuclear organization as a coupled self-assembly process. The chromosome energy function (see Eq. S7 in Supporting Information Section: Hi-C inspired interactions for the diploid human genome) includes terms that account for the polymer connectivity and excluded volume effect, an ideal potential, compartment-specific interactions, and specific interchromosomal interactions. As shown in Fig. 1C, the ideal potential is only applied for beads from the same chromosome to approximate the effect of loop extrusion by Cohesin molecules ^73,74 for chromosome compaction and territory formation.^75,76 Compartment-specific interactions, on the other hand, promote microphase separation and compartmentalization of euchromatin and heterochromatin. Finally, interchromosomal interactions account for sequence-specific effects that compartment-dependent potentials cannot capture.

Interactions among coarse-grained particles that form nuclear bodies were designed to promote and stabilize the formation of liquid droplets, as has been revealed by many experiments.^77–79 We adopted the Lennard-Jones potential for nucleolar particles to mimic the weak, multivalent interactions that arise from protein and RNA molecules that make up the nucleoli. As a first attempt to approximate their complex dynamics, we considered two types of particles that form speckles: phosphorylated (P) and de-phosphorylated (dP). The two types can interconvert via chemical reactions^80–82 and dP particles share attractive interactions modeled with the Lennard-Jones potential.

As shown in Fig. 1C, to recognize specific interactions between chromosomes and nuclear landmarks, we introduced contact potentials between them. These potentials are inspired by the experimental techniques that probe the corresponding contacts. The Supporting Information Section: Chromosome-nuclear landmark interactions and Nuclear landmark-nuclear landmark interactions contain more details about all the nuclear landmark related energy functions.

Optimization of model parameters with experimental data

The nucleus model was designed to be interpretable such that energy terms represent physical processes. Furthermore, the expressions of the interaction potentials were also designed such that their parameters can be determined from experimental data via the maximum entropy optimization algorithm.^9,83,84 Below, we briefly outline the procedure used for parameter optimization and further details can be found in the Supporting Information Section: Optimization of the whole nucleus model parameters.

As illustrated in Fig. 2A, starting from a given set of parameters, we first perform MD simulations to produce a collection of 3D structures for the diploid genome and various nuclear bodies. These structures are then transformed into a contact map or contact probabilities between chromatin beads and nuclear landmarks by averaging over homologous chromosomes. Constraints corresponding to different energy terms could be obtained from the simulated results and compared with those estimated from Hi-C, SON TSA-Seq, and Lamin B DamID profiles. Finally, the model parameters were updated based on the difference between simulated and experimental constraints using the adaptive moment estimation (Adam) optimization algorithm.⁸⁵ The three steps can be repeated with updated parameters to improve the simulation-experiment agreement further.

Overview of the iterative algorithm for parameterizing the nucleus model with experimental data.
Starting from an initial set of parameters, we perform MD simulations to produce an ensemble of nuclear structures. These structures can be transformed into contacts between chromosomes or between chromosomes and nuclear landmarks for direct comparison with experimental data. Differences between simulated and experiment contacts are used to update parameters for additional rounds of optimization if needed.

No quantitative experimental data exists for interactions among nuclear body particles to serve as constraints. We varied the strength of the interaction potential to produce 2-3 nucleoli and ∼ 30 speckle clusters during the simulations (Fig. S1) while ensuring the fluidity of the resulting droplets.

Molecular dynamics simulations with GPU acceleration

We implemented the nucleus model into the MD engine OpenMM.⁷¹ OpenMM offers an excellent interface with Python scripting, significantly improving the readability and customizability of the model. The code was designed into functional modules, with different components, such as chromosomes and nuclear landmarks, written as separate classes. This design further facilitates the introduction of additional nuclear components, if desired, with minimal changes to existing code. We provide examples of simulation set up, trajectory analysis, parameter optimization, and introducing new features in the GitHub repository.

Fig. 3A illustrates the workflow for setting up and executing whole nucleus simulations. A configuration file that provides the position of individual particles in the PDB file format is needed to initialize the simulations. This file also contains topological information regarding whether a particle represents chromosomes or nuclear landmarks and the identity of specific chromosomes. The input file can be generated with provided Python scripts by randomly distributing the positions of chromosomes, speckles, and nucleoli, though optimized configurations are also included in the GitHub repository. By default, the lamina particles will be uniformly placed on a sphere of 10 µm in diameter. Upon parsing the configuration file, interactions among various components can be set up with optimized parameters. This step will produce an object that can be used for MD simulations. As shown in Fig. 3B, the workflow only requires a few lines of code. The package also includes analysis scripts to compute contact maps, monitor conformational dynamics, and track nuclear bodies.

OpenNucleome facilitates GPU-accelerated simulations of the human nucleus.
(A) Illustration of workflow for setting up, performing, and analyzing MD simulations. (B) Python scripts setting up whole nucleus simulations. (C) Performance of MD simulations on different number of CPU cores and a single GPU.

A significant benefit of OpenMM⁷¹ is its native support of GPU acceleration. As shown in Fig. 3C, the simulation speed with one Nvidia Volta V100 GPU is 150 times faster than that of the 4 Intel Xeon Platinum 8260 CPU cores. Notably, this performance enhancement cannot be achieved by simply increasing the CPU core numbers. For example, the simulation speed with 32 CPU cores is less than twice that of 4 CPU cores, potentially due to the system’s heterogeneous distribution of particles.

Simulations reproduce and predict diverse experimental data

We extensively validated the parameterized nucleus model to examine its biological relevance. MD simulations initialized from 50 different initial configurations were performed to build an ensemble of structures. As mentioned in the following section, a diverse set of initial configurations is essential for reproducing interchromosomal contacts probed in Hi-C. From the simulated structures, we computed various quantities for direct comparison with experimental measurements. Given that the majority of experimental data were analyzed for the haploid genome, we adopted a similar approach by averaging over paternal and maternal chromosomes to facilitate direct comparison. More details on data analysis can be found in the Supporting Information Section: Details of simulation data analysis.

We compared the simulated contact probabilities among chromosomes with Hi-C data. As shown in Fig. 4A and Fig. S2, the simulated and experimental contact maps are highly correlated. The squares along the diagonal support the formation of chromosome territories that promote intrachromosomal contacts, and the apparent checkboard patterns follow the compartmentalization of various chromatin types. We further examined the decay of intra-chromosomal contacts as a function of the sequence separation, which is known to deviate from that of an equilibrium globule.²² As shown in Fig. 4B, the simulated results overlap well with the Hi-C data (orange curve). In addition, the simulated average contact probabilities between various compartment types match values estimated from Hi-C data. Moreover, the simulated and experimental average contact probabilities between pairs of chromosomes agree well, and the Pearson correlation coefficient between the two datasets reaches 0.89.

Simulated structures reproduce contact frequencies between chromosomes and between chromosomes and nuclear landmarks.
(A) Comparison between simulated (top right) and experimental (bottom left) whole-genome contact probability maps with Pearson correlation coefficient r = 0.89. Zoom-ins of various regions are provided in Fig. S2. (B) Comparison between simulated and experimental average contact frequencies, including average contacts between genomic loci from the same chromosomes at a given separation (top), average contacts between genomic loci classified into different compartment types (middle), and average contacts between various chromosome pairs (bottom). (C) Comparison between simulated and experimental Lamin-B DamID (top) and SON TSA-Seq signals (bottom), with Pearson correlation coefficients of haploid chromosomes shown on the right.

We further examined the contacts between chromosomes and nuclear landmarks. As illustrated in Fig. 4C, the simulated Lamin-B DamID signals for chromosome 7 match well with the experimental results, capturing the complex contact pattern that weaves chromatin towards and away from the nuclear envelope. Similarly, SON TSA-Seq data that quantify the contact between chromosomes and speckles are well captured by simulated structures. The anti-correlation between DamID and TSA-Seq is clearly visible. The observed agreement between simulation and experimental results is not limited to any particular chromosome. Good agreements are achieved for all chromosomes.

The simulations also provide 3D representations of the nucleus that can be compared with DNA-MERFISH data.³⁵ We found that the simulated radius of gyration of individual chromosomes matches well with experimental values (Fig. 5A). The simulated and experimental average normalized chromosome radial positions also correlate strongly, as shown in Fig. 5B. We note that while the sequencing results presented in Fig. 4 were used for model parameterization, the MERFISH data were not. Therefore, the simulation results here are de novo predictions, and their agreement with experimental data strongly supports the coupled assembly mechanism used for designing the energy function.

Structural and dynamical predictions of the nucleus model match results from microscopy imaging.
(A) Comparison between the simulated and experimental radius of gyration, R_g, for haploid chromosomes. The Pearson correlation coefficient between the two, r, is shown in the legend. (B) Comparison between the simulated and experimental normalized radial positions for haploid chromosomes, with their Pearson correlation coefficient shown in the legend. Detailed definition of the normalized radial positions is provided in the Supporting Information *Section: Computing simulated normalized chromosome radial positions*. (C) Mean-squared displacements, MSDs, as a function of time are shown for selected telomeres. (D) The probability distribution of the anomalous exponent, α, obtained from fitting the MSDs curves for all telomeres with the expression, ⟨r²(Δt)⟩ = *D_α*Δ*t^α*.

A significant advantage of MD simulation-based models is the dynamical information they naturally produce. We measured the dynamics of telomeres by tracking the mean square displacements (MSDs), ⟨r²(Δt)⟩, as a function of time. In Fig. 5C, we plot representative MSD trajectories over a one-hour timescale. In line with previous research, ^86–88 telomeres display anomalous subdiffusive motion. When fitted with the equation ⟨r²(Δt)⟩ = D_αΔt^α, these trajectories yield a spectrum of α values, with a peak around 0.59. The exponent and the diffusion coefficient D_α = (27±11)×10⁻⁴µm²·s^−α both match well with the experimental values,^89,90 upon setting the nucleoplasmic viscosity as 1Pa · s (see Supporting Information Section: Mapping the reduced time unit to real time for more details).

The good agreement in the dynamics of individual loci further inspired us to examine the diffusion of whole chromosomes. In particular, we plotted the normalized chromosome radial positions as a function of time in Fig. 6A. Remarkably, we found that chromosomes appear arrested and no significant changes in their positions are observed over timescales comparable to the cell cycle (see also Fig. S12). Therefore, our simulations predict that large-scale movements of chromosomes are unlikely during the G1 phase.

Heterogeneity and conserved features of nuclear organizations.
(A) Normalized chromosome radial positions as a function of simulation time. (B) Contacts between chromosome 1 and 2 from two independent simulation trajectories show significant variations. (C) Genome-wide *in silico* Lamin B DamID (top) and SON TSA-Seq (bottom) profiles computed from two independent trajectories. Pearson correlation coefficients, r, are provided on each plot. (D) Pairwise Person correlation coefficients between interchromosomal contact matrices (left), genome-wide Lamin B DamID profiles (middle), and genome-wide SON TSA-Seq profiles (right) determined from independent trajectories. The averages excluding the diagonals of the three datasets are 0.06, 0.53, and 0.72.

Heterogeneity and robustness of the simulated conformational ensemble

The lack of relaxation of chromosome radial positions suggests the importance of starting configurations used to initialize the simulations. Statistical averages of the resulting ensemble of nuclear structures depend crucially on these starting configurations. Using an optimization procedure, we selected them from 1000 configurations to maximize the agreement with experimental lamin-B DamID and interchromosomal contact probabilities. The Supporting Information Section: Initial configurations for simulations provides more details on preparing the 1000 initial configurations.

We selected a total of 50 starting configurations to initiate independent simulations. Smaller sets of starting configurations are not sufficient to reproduce the interchromosomal contact probabilities, as shown in Fig. S3B. Notably, different sets of 50 configurations selected from independent trials show significant overlap (Fig. S3D), supporting the robustness of the selection protocol in detecting conserved features of genome organization.

While the ensemble as a whole is relatively robust, individual configurations with the ensemble exhibit significant differences. For example, the Lamin B DamID profiles produced from different trajectories are only weakly correlated (Fig. 6C), with an average correlation coefficient of 0.53. These weak correlations result from significant differences in the normalized radial positions of chromosomes, as can be seen in representative configurations from two simulation trajectories (Fig. 6B). The fluctuations of normalized radial positions cause changes in contacts between chromosomes as well, resulting in little correlation between interchromosomal contact matrices (Fig. 6D).

We examined genome organizations reported by Su et al. ³⁵ and found a similar variation of interchromosomal contact probabilities across individual cells (Fig. S3A, D). Notably, the simulated configurations capture the fluctuations of interchromosomal contacts observed in DNA-MERFISH data, further supporting the biological relevance of the reported in silico structures.

Despite the differences in interchromosomal contacts across trajectories, high conservation of connections between chromosomes and speckles can be observed in individual simulations. For example, the average correlation coefficient between in silico SON TSA-Seq profiles produced from different trajectories is 0.72, much higher than the corresponding value for Lamin B DamID profiles. Conservation of contacts between chromosomes and nuclear bodies (zones) across individual cells has indeed been reported in a previous study that simultaneously images chromatin and various subnuclear structures.³⁶

Nuclear deformation preserves chromosome-nuclear body contacts

Numerous studies have highlighted the remarkable influence of nuclear shape on the positioning of chromosomes and the regulation of gene expression.^53,91 The nucleus, once regarded as a mere compartment for DNA storage, is increasingly recognized as a dynamic and intricately structured organelle. To better understand the interplay between nuclear shape and genome organization as a fundamental mechanism that shapes the transcriptional landscape, we performed additional simulations in which the nuclear lamina was altered from a sphere into more ellipsoidal shapes by applying a force along the z-axis (Fig. 7A). More details about these simulations can be found in the Section: Nuclear envelope deformation simulations of the Supporting Information.

Nuclear deformations influence genome organization while preserving chromatin-speckle contacts.
(A) Illustration of force induced nuclear envelope deformation. The nuclear lamina is modeled as a particle mesh where neighboring lamina particles are covalently bonded together. (B) Example nucleus conformations at different strengths of applied force. (C) Pearson correlation coefficients between results from simulations of deformed nuclei and those from a spherical nucleus for interchromosomal contacts (left), DamID profiles (middle), and TSA-Seq (right). The values at zero force were computed from two independent simulations starting from the same initial configurations.

As illustrated in Fig. 7B, the presence of external forces resulted in significant alterations in nuclear shape. We conducted two independent simulations with different force strengths, leading to varying degrees of deformation in the nuclear lamina. This deformation, in turn, caused a reorganization of chromosomes, affecting their normalized radial positions and pairwise contacts (see Fig. S4 and Fig. 7C). We observed that more deformed nuclei exhibited lower correlation coefficients for interchromosomal contacts compared to results obtained from simulations in a spherical nucleus. Similarly, the DamID profiles exhibited significant variations upon nucleus deformation, whereas TSA-Seq signals were much less affected and remained highly correlated with the results from the spherical nucleus simulations.

Therefore, it appears that speckles, and potentially other nuclear condensates, can dynamically reorganize in response to changes in chromosome conformations to maintain contacts with genomic loci. This robustness in nuclear body contacts may be essential for ensuring the robust functioning of the genome in a population of cells with significant variability in nuclear shape.

Conclusions and Discussion

We introduced a computational model, OpenNucleome, to facilitate simulations for the human nucleus at high structural and temporal resolution. We conducted extensive crossvalidation with experimental data to support the biological relevance of simulated 3D structures. Implementing the model into the MD package, OpenMM enables GPU acceleration for long-timescale simulations. Tutorials in the format of Python Scripts with extensive documentation are provided to facilitate the adoption of the model by the community.

Our software enhances the capabilities of existing genome simulation tools. ^40,69,70 Specifically, OpenNucleome aligns with the design principles of Open-MiChroM, ⁷⁰ prioritizing open-source accessibility while expanding simulation capabilities to the entire nucleus. Similar to software from the Alber lab,⁶⁹ OpenNucleome offers high-resolution genome organization that faithfully reproduces a diverse range of experimental data. Furthermore, beyond static structures, OpenNucleome facilitates dynamic simulations with explicit representations of various nuclear condensates, akin to the model developed by Fujishiro and Sasai ⁴⁰.

A significant advantage of OpenNucleome lies in its predictive power for dynamical information. For example, the model succeeded in reproducing the subdiffusive behavior of telomeres. We further showed that the dynamics of individual chromosomes are slow and their radial positions do not relax over the time course of a cell cycle. This is consistent with previous theoretical estimations on chromosome dynamics ⁹² and recent observations of solid behavior of chromatin in vivo.⁹³ Live cell experiments that directly track the positions of multiple chromosomes could further validate/falsify this prediction. We anticipate the model will greatly facilitate the investigation of the dynamics of genomic loci and nuclear bodies and the interpreting of live cell imaging results.

Slow chromosome dynamics and a lack of conformational relaxation naturally result in the heterogeneity of chromosome radial positions across individual cells. This heterogeneity raises doubts about the notion that chromosome radial positions provide robust and reliable mechanisms for gene regulation.^2,94–96 Instead, our results support the nuclear zoning model for gene regulation,³⁶ where specific loci function as “fixed points” anchored to certain nuclear bodies in all cells. This anchoring mechanism robustly creates the desired molecular environment surrounding these genomic segments. Unlike chromosome radial positions, contacts between genomic loci and speckles can be robustly established in individual cells, as shown in our simulations. It was achieved through a nucleation process that attracts speckle particles towards specific loci due to specific interactions. Nucleation occurs much more rapidly than chromosome rearrangement due to the smaller size of speckle particles. The coupled self-assembly mechanism for chromosomes and nuclear bodies can similarly facilitate the formation of other nuclear zones for different kinds of fixed points.

Despite the heterogeneity in chromosome positions and interchromosomal contacts, the ensemble of nuclear structures as a whole is not random and exhibits conserved features. For example, on average, certain chromosomes remain closer to the nuclear envelope than others (see Fig. 5B). Similarly, the average contact frequency between certain chromosome pairs is higher than others, though this trend can be frequently violated in individual cells. How such conserved features arise as cells exit from the mitotic phase remains unclear and would be interesting for further explorations.

Methods

Molecular dynamics simulation details

We used the software package OpenMM⁷¹ to perform MD simulations in reduced units at constant temperature (T = 1.0). Unless otherwise specified, we froze the lamina particles and only propagated the dynamics of chromatin, nucleoli, and speckles.

Two integration schemes were used with a time step of dt = 0.005 to efficiently generate structural ensembles and produce realistic dynamical information, respectively. For simulations used in parameter optimization and building structural ensembles, we employed the Langevin integrator with a damping coefficient of γ⁻¹ = 10.0. In the case of MSD calculations shown in Fig. 6, we utilized Brownian dynamics with a damping coefficient of γ⁻¹ = 0.01. The higher damping coefficient provides a better approximation to the viscous nucleus environment, while the smaller value in the Langevin integrator facilitates conformational sampling with faster diffusion rates.

We employed the semi-grand Monte Carlo technique⁹⁷ to simulate chemical transitions between two types of speckle particles. At every 4000 simulation steps, we attempt a total of N_Sp chemical reactions that converts one type of speckle particles to the other type with a probability of 0.2. N_Sp corresponds to the total number of speckle particles, and the switching probability was chosen to be comparable to the experimental phosphorylation rate. More details on the speckle dynamics are provided in the Supporting Information Section: Speckles as phase-separated droplets undergoing chemical modifications.

When deforming the nuclear envelope, we unfroze the lamina particles and evolved them dynamically as the rest of the nucleus. Bonded interactions among lamina particles held the nuclear envelope together as a particle mesh. A harmonic force along the z-axis was introduced to compress the particle mesh. More details are provided in the Supporting Information Section: Nuclear envelope deformation simulations.

For simulations used to optimize parameters, a total of 50 independent 3-million-step-long trajectories were performed. Configurations were recorded at every 2000 simulation steps for analysis. The first 500,000 steps of each trajectory were discarded as equilibration. For production simulations, we performed 50 independent 10-million-step long trajectories starting from different initial configurations. Nuclear structures were again recorded at every 2000 steps to determine statistical averages presented in the paper. An additional 8 simulations of 30 million steps in length were performed to compute telomere MSDs.

We mapped the reduced units to real units with the conversion of length scale σ = 385 nm and the time scale in Brownian dynamics simulations τ = 0.65 s. These conversions were determined by as detailed in the Supporting Information Section: Unit Conversion.

Experimental data processing and analysis

We obtained the in situ Hi-C data, SON TSA-seq data, and Lamin-B DamID data of HFF cell lines from the 4DN data portal. The intra and interchromosomal interactions were calculated at 100 KB resolution with VC SQRT normalization applied to the interaction matrices. Hi-C data extraction and normalization were performed using Juicer tools.⁹⁸ We followed the same processing and normalization method described in Ref. 99 to analyze TSA-seq data. Two biological replicates of Lamin-B DamID data were merged and the normalized counts over Dam-only control were used for analysis. The SON TSA-Seq and Lamin-B DamID data were processed at the 25 KB resolution and the average values at the 100 KB resolution were used in Fig. 4 for model validation.

Acknowledgements

This work was supported by the National Institutes of Health (Grant R35GM133580).

Competing interests

Authors declare that they have no competing interests.

Data and materials availability

Hi-C data (https://data.4dnucleome.org, accession number: 4DNFIB59T7NN). SON TSA-seq data (https://data.4dnucleome.org, accession number: pulldown data 4DNEX6U8TS3Y, control data 4DNEXI7XUWFK). Lam-inB DamID data (https://data.4dnucleome.org, accession number 4DNESXZ4FW4T). The software is available at: available at https://github.com/ZhangGroup-MITChemistry/OpenNucleome.

References

(1)
1. Dekker J.
2. Belmont A. S.
3. Guttman M.
4. Leshyk V. O.
5. Lis J. T.
6. Lomvardas S.
7. Mirny L. A.
8. O’shea C. C.
9. Park P. J.
10. Ren B.
11. et al.
2017The 4D nucleome projectNature 549:219–226
(2)
1. Hübner M. R.
2. Eckersley-Maslin M. A.
3. Spector D. L.
2013Chromatin organization and transcriptional regulationCurrent opinion in genetics & development 23:89–95
(3)
1. Bickmore W. A
2013The spatial organization of the human genomeAnnual review of genomics and human genetics 14:67–84
(4)
1. Gorkin D. U.
2. Leung D.
3. Ren B
2014The 3D genome in transcriptional regulation and pluripotencyCell stem cell 14:762–775
(5)
1. Dekker J.
2. Mirny L
2016The 3D genome as moderator of chromosomal communicationCell 164:1110–1121
(6)
1. Furlong E. E.
2. Levine M.
2018Developmental enhancers and chromosome topologyScience 361:1341–1345
(7)
1. Finn E. H.
2. Misteli T
2019Molecular basis and biological function of variability in spatial genome organizationScience 365
(8)
1. Chen Y.
2. Belmont A. S
2019Genome organization around nuclear specklesCurrent opinion in genetics & development 55:91–99
(9)
1. Lin X.
2. Qi Y.
3. Latham A. P.
4. Zhang B
2021Multiscale modeling of genome organization with maximum entropy optimizationJ. Chem. Phys 155
(10)
1. Liu S.
2. Athreya A.
3. Lao Z.
4. Zhang B
2024From Nucleosomes to Compartments: Physicochemical Interactions Underlying Chromatin OrganizationAnnual Review of Biophysics 53
(11)
1. Lafontaine D. L.
2. Riback J. A.
3. Bascetin R.
4. Brangwynne C. P
2021The nucleolus as a multiphase liquid condensateNature reviews Molecular cell biology 22:165–182
(12)
1. Lamond A. I.
2. Spector D. L
2003Nuclear speckles: a model for nuclear organellesNature reviews Molecular cell biology 4:605–612
(13)
1. Van Steensel B.
2. Belmont A. S
2017Lamina-associated domains: links with chromosome architecture, heterochromatin, and gene repressionCell 169:780–791
(14)
1. Seruga B.
2. Zhang H.
3. Bernstein L. J.
4. Tannock I. F
2008Cytokines and their relationship to the symptoms and outcome of cancerNature Reviews Cancer 8:887–899
(15)
1. Schuster-Böckler B.
2. Lehner B
2012Chromatin organization is a major influence on regional mutation rates in human cancer cellsnature 488:504–507
(16)
1. Schmitt A. D.
2. Hu M.
3. Ren B
2016Genome-wide mapping and analysis of chromosome architectureNature reviews Molecular cell biology 17:743–755
(17)
1. McCord R. P.
2. Kaplan N.
3. Giorgetti L
2020Chromosome conformation capture and beyond: toward an integrative view of chromosome structure and functionMolecular cell 77:688–708
(18)
1. Parmar J. J.
2. Woringer M.
3. Zimmer C
2019How the genome folds: the biophysics of four-dimensional chromatin organizationAnnual review of biophysics 48:231–253
(19)
1. Jerkovic I.
2. Cavalli G
2021Understanding 3D genome organization by multidisciplinary methodsNature Reviews Molecular Cell Biology 22:511–528
(20)
1. Chen B.
2. Guan J.
3. Huang B.
2016Imaging specific genomic DNA in living cellsAnnu Rev Biophys
(21)
1. Dekker J.
2. Rippe K.
3. Dekker M.
4. Kleckner N
2002Capturing chromosome conformationscience 295:1306–1311
(22)
1. Lieberman-Aiden E.
2. Van Berkum N. L.
3. Williams L.
4. Imakaev M.
5. Ragoczy T.
6. Telling A.
7. Amit I.
8. Lajoie B. R.
9. Sabo P. J.
10. Dorschner M. O.
11. et al.
2009Comprehensive mapping of long-range interactions reveals folding principles of the human genomescience 326:289–293
(23)
1. Rao S. S.
2. Huntley M. H.
3. Durand N. C.
4. Stamenova E. K.
5. Bochkov I. D.
6. Robinson J. T.
7. Sanborn A. L.
8. Machol I.
9. Omer A. D.
10. Lander E. S.
11. et al.
2014A 3D map of the human genome at kilobase resolution reveals principles of chromatin loopingCell 159:1665–1680
(24)
1. Dixon J. R.
2. Gorkin D. U.
3. Ren B
2016Chromatin domains: the unit of chromosome organizationMolecular cell 62:668–680
(25)
1. Dekker J.
2. Heard E
2015Structural and functional diversity of topologically associating domainsFEBS letters 589:2877–2884
(26)
1. Greil F.
2. Moorman C.
3. van Steensel B
2006[16] DamID: mapping of in vivo protein– genome interactions using tethered DNA adenine methyltransferaseMethods in enzymology 410:342–359
(27)
1. Park P. J
2009ChIP–seq: advantages and challenges of a maturing technologyNature reviews genetics 10:669–680
(28)
1. Chen Y.
2. Zhang Y.
3. Wang Y.
4. Zhang L.
5. Brinkman E. K.
6. Adam S. A.
7. Goldman R.
8. Van Steensel B.
9. Ma J.
10. Belmont A. S
2018Mapping 3D genome organization relative to nuclear compartments using TSA-Seq as a cytological rulerJournal of Cell Biology 217:4025–4048
(29)
1. Wen S.
2. Ma D.
3. Zhao M.
4. Xie L.
5. Wu Q.
6. Gou L.
7. Zhu C.
8. Fan Y.
9. Wang H.
10. Yan J
2020Spatiotemporal single-cell analysis of gene expression in the mouse suprachiasmatic nucleusNature neuroscience 23:456–467
(30)
1. Ramani V.
2. Deng X.
3. Qiu R.
4. Gunderson K. L.
5. Steemers F. J.
6. Disteche C. M.
7. Noble W. S.
8. Duan Z.
9. Shendure J
2017Massively multiplex single-cell Hi-CNature methods 14:263–266
(31)
1. Nagano T.
2. Lubling Y.
3. Stevens T. J.
4. Schoenfelder S.
5. Yaffe E.
6. Dean W.
7. Laue E. D.
8. Tanay A.
9. Fraser P
2013Single-cell Hi-C reveals cell-to-cell variability in chromosome structureNature 502:59–64
(32)
1. Chen K. H.
2. Boettiger A. N.
3. Moffitt J. R.
4. Wang S.
5. Zhuang X
2015Spatially resolved, highly multiplexed RNA profiling in single cellsScience 348
(33)
1. Boettiger A. N.
2. Bintu B.
3. Moffitt J. R.
4. Wang S.
5. Beliveau B. J.
6. Fudenberg G.
7. Imakaev M.
8. Mirny L. A.
9. Wu C.-t.
10. Zhuang X.
2016Super-resolution imaging reveals distinct chromatin folding for different epigenetic statesNature 529:418–422
(34)
1. Shachar S.
2. Voss T. C.
3. Pegoraro G.
4. Sciascia N.
5. Misteli T
2015Identification of gene positioning factors using high-throughput imaging mappingCell 162:911–923
(35)
1. Su J.-H.
2. Zheng P.
3. Kinrot S. S.
4. Bintu B.
5. Zhuang X
2020Genome-scale imaging of the 3D organization and transcriptional activity of chromatinCell 182:1641–1659
(36)
1. Takei Y.
2. Yun J.
3. Zheng S.
4. Ollikainen N.
5. Pierson N.
6. White J.
7. Shah S.
8. Thomassie J.
9. Suo S.
10. Eng C.-H. L.
11. et al.
2021Integrated spatial genomics reveals global architecture of single nucleiNature 590:344–350
(37)
1. Qi Y.
2. Reyes A.
3. Johnstone S. E.
4. Aryee M. J.
5. Bernstein B. E.
6. Zhang B
2020Data-driven polymer model for mechanistic exploration of diploid genome organizationBiophysical Journal 119:1905–1916
(38)
1. Qi Y.
2. Zhang B
2019Predicting three-dimensional genome organization with chromatin statesPLoS computational biology 15
(39)
1. Boninsegna L.
2. Yildirim A.
3. Polles G.
4. Zhan Y.
5. Quinodoz S. A.
6. Finn E. H.
7. Guttman M.
8. Zhou X. J.
9. Alber F
2022Integrative genome modeling platform reveals essentiality of rare contact events in 3D genome organizationsNature Methods 19:938–949
(40)
1. Fujishiro S.
2. Sasai M
2022Generation of dynamic three-dimensional genome structure through phase separation of chromatinProceedings of the National Academy of Sciences 119
(41)
1. Shi G.
2. Thirumalai D
2021From Hi-C contact map to three-dimensional organization of interphase human chromosomesPhysical Review X 11
(42)
1. Dekker J.
2. Marti-Renom M. A.
3. Mirny L. A
2013Exploring the three-dimensional organization of genomes: interpreting chromatin interaction dataNature Reviews Genetics 14:390–403
(43)
1. Jost D.
2. Carrivain P.
3. Cavalli G.
4. Vaillant C
2014Modeling epigenome folding: formation and dynamics of topologically associated chromatin domainsNucleic acids research 42:9553–9561
(44)
1. Giorgetti L.
2. Galupa R.
3. Nora E. P.
4. Piolot T.
5. Lam F.
6. Dekker J.
7. Tiana G.
8. Heard E
2014Predictive polymer modeling reveals coupled fluctuations in chromosome conformation and transcriptionCell 157:950–963
(45)
1. Di Pierro M.
2. Cheng R. R.
3. Lieberman Aiden E.
4. Wolynes P. G.
5. Onuchic J. N
2017De novo prediction of human chromosome structures: Epigenetic marking patterns encode genome architectureProceedings of the National Academy of Sciences 114:12126–12131
(46)
1. Buckle A.
2. Brackley C. A.
3. Boyle S.
4. Marenduzzo D.
5. Gilbert N
2018Polymer simulations of heteromorphic chromatin predict the 3D folding of complex genomic lociMolecular cell 72:786–797
(47)
1. Nuebler J.
2. Fudenberg G.
3. Imakaev M.
4. Abdennur N.
5. Mirny L
2018Chromatin organization by an interplay of loop extrusion and compartmental segregationBiophysical Journal 114
(48)
1. Bianco S.
2. Lupiáñez D. G.
3. Chiariello A. M.
4. Annunziatella C.
5. Kraft K.
6. Schöpflin R.
7. Wittler L.
8. Andrey G.
9. Vingron M.
10. Pombo A.
11. et al.
2018Polymer physics predicts the effects of structural variants on chromatin architectureNature genetics 50:662–667
(49)
1. Shi G.
2. Liu L.
3. Hyeon C.
4. Thirumalai D
2018Interphase human chromosome exhibits out of equilibrium glassy dynamicsNature communications 9
(50)
1. MacPherson Q.
2. Beltran B.
3. Spakowitz A. J
2018Bottom–up modeling of chromatin segregation due to epigenetic modificationsProceedings of the National Academy of Sciences 115:12739–12744
(51)
1. Shin S.
2. Shi G.
3. Thirumalai D
2023From Effective Interactions Extracted Using Hi-C Data to Chromosome Structures in Conventional and Inverted NucleiPRX Life 1
(52)
1. Amiad-Pavlov D.
2. Lorber D.
3. Bajpai G.
4. Reuveny A.
5. Roncato F.
6. Alon R.
7. Safran S.
8. Volk T
2021Live imaging of chromatin distribution reveals novel principles of nuclear architecture and chromatin compartmentalizationSci. Adv 7
(53)
1. Brahmachari S.
2. Contessoto V. G.
3. Di Pierro M.
4. Onuchic J. N
2022Shaping the genome via lengthwise compaction, phase separation, and lamina adhesionNucleic Acids Res 50:1–14
(54)
1. Jiang Z.
2. Qi Y.
3. Kamat K.
4. Zhang B
2022Phase Separation and Correlated Motions in Motorized GenomeJ. Phys. Chem. B 126:5619–5628
(55)
1. Ganai N.
2. Sengupta S.
3. Menon G. I
2014Chromosome positioning from activity-based segregationNucleic Acids Res 42:4145–4159
(56)
1. Liu S.
2. Zhang L.
3. Quan H.
4. Tian H.
5. Meng L.
6. Yang L.
7. Feng H.
8. Gao Y. Q
2018From 1D sequence to 3D chromatin dynamics and cellular functions: A phase separation perspectiveNucleic Acids Res 46:9367–9383
(57)
1. Laghmach R.
2. Di Pierro M.
3. Potoyan D. A
2020Mesoscale Liquid Model of Chromatin Recapitulates Nuclear Order of EukaryotesBiophys. J 118:2130–2140
(58)
1. Chu X.
2. Wang J
2021Deciphering the molecular mechanism of the cancer formation by chromosome structural dynamicsPLOS Comput. Biol 17
(59)
1. Lappala A.
2. Wang C.-Y.
3. Kriz A.
4. Michalk H.
5. Tan K.
6. Lee J. T.
7. Sanbon-matsu K. Y.
2021Four-dimensional chromosome reconstruction elucidates the spatiotemporal reorganization of the mammalian X chromosomeProc. Natl. Acad. Sci 118
(60)
1. Chu X.
2. Wang J
2022Quantifying Chromosome Structural Reorganizations during Differentiation, Reprogramming, and TransdifferentiationPhys. Rev. Lett 129
(61)
1. Goychuk A.
2. Kannan D.
3. Chakraborty A. K.
4. Kardar M
2023Polymer folding through active processes recreates features of genome organizationProceedings of the National Academy of Sciences 120
(62)
1. Sun Q.
2. Perez-Rathke A.
3. Czajkowsky D. M.
4. Shao Z.
5. Liang J
2021High-resolution single-cell 3D-models of chromatin ensembles during Drosophila embryogenesisNat. Commun 12:1–12
(63)
1. Kadam S.
2. Kumari K.
3. Manivannan V.
4. Dutta S.
5. Mitra M. K.
6. Padinhateeri R
2023Predicting scale-dependent chromatin polymer properties from systematic coarse-grainingNat. Commun 14
(64)
1. Bajpai G.
2. Pavlov D. A.
3. Lorber D.
4. Volk T.
5. Safran S
2021Mesoscale phase separation of chromatin in the nucleusElife 10
(65)
1. Kamat K.
2. Lao Z.
3. Qi Y.
4. Wang Y.
5. Ma J.
6. Zhang B
2023Compartmentalization with nuclear landmarks yields random, yet precise, genome organizationBiophysical Journal 122:1376–1389
(66)
1. Laghmach R.
2. Di Pierro M.
3. Potoyan D. A
2021The interplay of chromatin phase separation and lamina interactions in nuclear organizationBiophysical Journal 120:5005–5017
(67)
1. Stephens A. D.
2. Banigan E. J.
3. Marko J. F
2018Separate roles for chromatin and lamins in nuclear mechanicsNucleus 9:119–124
(68)
1. Qi Y.
2. Zhang B
2021Chromatin network retards nucleoli coalescenceNature Communications 12
(69)
1. Yildirim A.
2. Hua N.
3. Boninsegna L.
4. Zhan Y.
5. Polles G.
6. Gong K.
7. Hao S.
8. Li W.
9. Zhou X. J.
10. Alber F
2023Evaluating the role of the nuclear microenvironment in gene function by population-based modelingNature Structural & Molecular Biology :1–14
(70)
1. Junior A. B. O.
2. Contessoto V. G.
3. Mello M. F.
4. Onuchic J. N
2021A scalable computational approach for simulating complexes of multiple chromosomesJournal of molecular biology 433
(71)
1. Eastman P.
2. Swails J.
3. Chodera J. D.
4. McGibbon R. T.
5. Zhao Y.
6. Beauchamp K. A.
7. Wang L.-P.
8. Simmonett A. C.
9. Harrigan M. P.
10. Stern C. D.
11. et al.
2017OpenMM 7: Rapid development of high performance algorithms for molecular dynamicsPLoS computational biology 13
(72)
1. Strom A. R.
2. Biggs R. J.
3. Banigan E. J.
4. Wang X.
5. Chiu K.
6. Herman C.
7. Collado J.
8. Yue F.
9. Ritland Politz J. C.
10. Tait L. J.
11. et al.
2021HP1α is a chromatin crosslinker that controls nuclear and mitotic chromosome mechanicsElife 10
(73)
1. Sanborn A. L.
2. Rao S. S.
3. Huang S.-C.
4. Durand N. C.
5. Huntley M. H.
6. Jewett A. I.
7. Bochkov I. D.
8. Chinnappan D.
9. Cutkosky A.
10. Li J.
11. et al.
2015Chromatin extrusion explains key features of loop and domain formation in wild-type and engineered genomesProceedings of the National Academy of Sciences 112:E6456–E6465
(74)
1. Fudenberg G.
2. Imakaev M.
3. Lu C.
4. Goloborodko A.
5. Abdennur N.
6. Mirny L. A
2016Formation of chromosomal domains by loop extrusionCell reports 15:2038–2049
(75)
1. Di Pierro M.
2. Zhang B.
3. Aiden E. L.
4. Wolynes P. G.
5. Onuchic J. N
2016Transferable model for chromosome architectureProceedings of the National Academy of Sciences 113:12168–12173
(76)
1. Zhang B.
2. Wolynes P. G.
2017Genomic energy landscapesBiophysical journal 112:427–433
(77)
1. Handwerger K. E.
2. Cordero J. A.
3. Gall J. G
2005Cajal bodies, nucleoli, and speckles in the Xenopus oocyte nucleus have a low-density, sponge-like structureMolecular biology of the cell 16:202–211
(78)
1. Caragine C. M.
2. Haley S. C.
3. Zidovska A
2018Surface fluctuations and coalescence of nucleolar droplets in the human cell nucleusPhysical review letters 121
(79)
1. Caragine C. M.
2. Haley S. C.
3. Zidovska A
2019Nucleolar dynamics and interactions with nucleoplasm in living cellsElife 8
(80)
1. Brackley C. A.
2. Liebchen B.
3. Michieletto D.
4. Mouvet F.
5. Cook P. R.
6. Marenduzzo D
2017Ephemeral protein binding to DNA shapes stable nuclear bodies and chromatin domainsBiophysical journal 112:1085–1093
(81)
1. Söding J.
2. Zwicker D.
3. Sohrabi-Jahromi S.
4. Boehning M.
5. Kirschbaum J
2020Mechanisms for active regulation of biomolecular condensatesTrends in Cell Biology 30:4–14
(82)
1. Carrero G.
2. Hendzel M.
3. De Vries G
2006Modelling the compartmentalization of splicing factorsJournal of theoretical biology 239:298–312
(83)
1. Xie W. J.
2. Zhang B
2019Learning the Formation Mechanism of Domain-Level Chromatin States with Epigenomics DataBiophys. J 116:2047–2056
(84)
1. Schuette G.
2. Ding X.
3. Zhang B
2023Efficient Hi-C inversion facilitates chromatin folding mechanism discovery and structure predictionBiophys. J 122:3425–3438
(85)
1. Kingma D. P.
2. Ba J.
2014Adam: A method for stochastic optimizationarXiv
(86)
1. Di Pierro M.
2. Potoyan D. A.
3. Wolynes P. G.
4. Onuchic J. N
2018Anomalous diffusion, spatial coherence, and viscoelasticity from the energy landscape of human chromosomesProceedings of the National Academy of Sciences 115:7753–7758
(87)
1. Bronstein I.
2. Israel Y.
3. Kepten E.
4. Mai S.
5. Shav-Tal Y.
6. Barkai E.
7. Garini Y
2009Transient anomalous diffusion of telomeres in the nucleus of mammalian cellsPhysical review letters 103
(88)
1. Lee D. S.
2. Wingreen N. S.
3. Brangwynne C. P
2021Chromatin mechanics dictates subdiffusion and coarsening dynamics of embedded condensatesBiophysical Journal 120
(89)
1. Bronshtein I.
2. Kepten E.
3. Kanter I.
4. Berezin S.
5. Lindner M.
6. Redwood A. B.
7. Mai S.
8. Gonzalo S.
9. Foisner R.
10. Shav-Tal Y.
11. et al.
2015Loss of lamin A function increases chromatin dynamics in the nuclear interiorNature communications 6
(90)
1. Jack A.
2. Kim Y.
3. Strom A. R.
4. Lee D. S.
5. Williams B.
6. Schaub J. M.
7. Kellogg E. H.
8. Finkelstein I. J.
9. Ferro L. S.
10. Yildiz A.
11. et al.
2022Compartmentalization of telomeres through DNA-scaffolded phase separationDevelopmental cell 57:277–290
(91)
1. Contessoto V. G.
2. Dudchenko O.
3. Aiden E. L.
4. Wolynes P. G.
5. Onuchic J. N.
6. Di Pierro M
2023Interphase chromosomes of the Aedes aegypti mosquito are liquid crystalline and can sense mechanical cuesNature Communications 14
(92)
1. Rosa A.
2. Everaers R
2008Structure and dynamics of interphase chromosomesPLoS Comput. Biol 4
(93)
1. Strickfaden H.
2. Tolsma T. O.
3. Sharma A.
4. Underhill D. A.
5. Hansen J. C.
6. Hendzel M. J
2020Condensed chromatin behaves like a solid on the mesoscale in vitro and in living cellsCell 183:1772–1784
(94)
1. Maeshima K.
2. Hihara S.
3. Eltsov M
2010Chromatin structure: does the 30-nm fibre exist in vivo?Current opinion in cell biology 22:291–297
(95)
1. Fraser P.
2. Bickmore W
2007Nuclear organization of the genome and the potential for gene regulationNature 447:413–417
(96)
1. Takizawa T.
2. Meaburn K. J.
3. Misteli T
2008The meaning of gene positioningCell 135:9–13
(97)
1. Sadigh B.
2. Erhart P.
3. Stukowski A.
4. Caro A.
5. Martinez E.
6. Zepeda-Ruiz L
2012Scalable parallel Monte Carlo algorithm for atomistic simulations of precipitation in alloysPhysical Review B 85
(98)
1. Durand N. C.
2. Shamim M. S.
3. Machol I.
4. Rao S. S.
5. Huntley M. H.
6. Lander E. S.
7. Aiden E. L
2016Juicer provides a one-click system for analyzing loop-resolution Hi-C experimentsCell systems 3:95–98
(99)
1. Zhang L.
2. Zhang Y.
3. Chen Y.
4. Gholamalamdari O.
5. Wang Y.
6. Ma J.
7. Belmont A. S
2021TSA-seq reveals a largely conserved genome organization relative to nuclear speckles with small position changes tightly correlated with gene expression changesGenome research 31:251–264

Article and author information

Author information

Zhuohan Lao
Department of Chemistry, Massachusetts Institute of Technology, Cambridge, MA, USA
ORCID iD: 0000-0001-5404-2183
Kartik Kamat
Department of Chemistry, Massachusetts Institute of Technology, Cambridge, MA, USA
ORCID iD: 0000-0002-0909-8839
Zhongling Jiang
Department of Chemistry, Massachusetts Institute of Technology, Cambridge, MA, USA
Bin Zhang
Department of Chemistry, Massachusetts Institute of Technology, Cambridge, MA, USA
ORCID iD: 0000-0002-3685-7503
- E-mail: binz@mit.edu

Version history

Sent for peer review: October 11, 2023
Preprint posted: October 18, 2023
Reviewed Preprint version 1: February 8, 2024
Reviewed Preprint version 2: July 1, 2024
Version of Record published: August 15, 2024

Copyright

This article is distributed under the terms of the Creative Commons Attribution License, which permits unrestricted use and redistribution provided that the original author and source are credited.

Revised: This Reviewed Preprint has been revised by the authors in response to the previous round of peer review; the eLife assessment and the public reviews have been updated where necessary by the editors and peer reviewers.

Reviewing Editor
Rosana Collepardo
University of Cambridge, Cambridge, United Kingdom
Senior Editor
Yamini Dalal
National Cancer Institute, Bethesda, United States of America

Reviewer #1 (Public Review):

Summary:

In this paper the authors develop a comprehensive program to investigate the organization of chromosome structures at 100 kb resolution. It is extremely well executed. The authors have thought through all aspects of the problem. The resulting software will be most useful to the community. Interestingly they capture many experimental observations accurately. I have very little complaints.

Strengths:

A lot of details are provided. The success of the method is well illustrated. Software is easily available,

Weaknesses:

The number of parameters in the energy function is very large. Any justification? Could they simply be the functions?

What would the modification be if the resolution is increased?

They should state that the extracted physical values are scale dependent. Example, viscosity.

https://doi.org/10.7554/eLife.93223.2.sa2

Reviewer #2 (Public Review):

Summary:

In this work, Lao et al. develop an open-source software (OpenNucleome) for GPU-accelerated molecular dynamics simulation of the human nucleus accounting for chromatin, nucleoli, nuclear speckles, etc. Using this, the authors investigate the steady-state organization and dynamics of many of the nuclear components.

Strengths:

This is a comprehensive open-source tool to study several aspects of the nucleus, including chromatin organization, interactions with lamins and organization, and interactions with nuclear speckles and nucleoli. The model is built carefully, accounting for several important factors and optimizing the parameters iteratively to achieve experimentally known results. Authors have simulated the entire genome at 100kb resolution (which is a very good resolution to simulate and study the entire diploid genome) and predict several static quantities such as the radius of gyration and radial positions of all chromosomes, and time-dependent quantities like the mean-square displacement of important genomic regions.

Weaknesses:

One weakness of the model is that it has several parameters. Some of them are constrained by the experiments. However, the role of every parameter is not clear in the manuscript.

https://doi.org/10.7554/eLife.93223.2.sa1

Reviewer #3 (Public Review):

Summary:

The authors present OpenNucleome, a computational tool for simulating the structure and dynamics of the human nucleus. The software models nuclear components, including chromosomes and nuclear bodies, and incorporates GPU acceleration for potential performance gains. The authors aim to advance the understanding of nuclear organization by providing a tool that aligns with experimental data and is accessible to the genome architecture research community.

Strengths:

OpenNucleome provides a model of the nucleus, contributing to the advancement of computational biology.
Utilizing GPU acceleration with OpenMM may offer potential performance improvements.

Weaknesses:

It could still take advantage of clearer explanations regarding the generation and usage of input and output files and compatibility with other tools.

https://doi.org/10.7554/eLife.93223.2.sa0

Author response:

The following is the authors’ response to the original reviews.

Reviewer 1:

Comment 0: In this paper, the authors develop a comprehensive program to investigate the organization of chromosome structures at 100 kb resolution. It is extremely well executed. The authors have thought through all aspects of the problem. The resulting software will be most useful to the community. Interestingly they capture many experimental observations accurately.

I have very few complaints.

We appreciate the reviewer’s strong assessment of the paper’s significance, novelty, and broad interest, and we thank them for the detailed suggestions and comments.

Comment 1: The number of parameters in the energy function is very large. Is there any justification for this? Could they simplify the functions?

We extend our gratitude to the reviewer for their insightful remarks. The parameters within our model can be categorized into two groups: those governing chromosome-chromosome interactions and those governing chromosome-nuclear landmark interactions.

In terms of chromosome-chromosome interactions, the parameter count is relatively modest compared to the vast amount of Hi-C data available. For instance, while the whole-genome Hi-C matrix at the 100KB resolution encompasses approximately 303212 contacts, our model comprises merely six parameters for interactions among different compartments, along with 1000 parameters for the ideal potential. As outlined in the supporting information, the ideal potential is contingent upon sequence separation, with 1000 chosen to encompass bead separations of up to 100MB. While it is theoretically plausible to reduce the number of parameters by assuming interactions cease beyond a certain sequence separation, determining this scale a priori presents a challenge.

During the parameterization process, we observed that interchromosomal contacts predicted solely based on compartmental interactions inadequately mirrored Hi-C data. Consequently, we introduced 231 additional parameters to more accurately capture interactions between distinct pairs of autosomes. These interactions may stem from factors such as non-coding RNA or proteins not explicable by simple, non-specific compartmental interactions.

Regarding parameters concerning chromosome-nuclear landmark interactions, we have 30321 parameters for speckles and 30321 for the nuclear lamina. To streamline the model, we opted to assign a unique parameter to each chromatin bead. However, it is conceivable that many chromatin beads share a similar mechanism for interacting with nuclear lamina or speckles, potentially allowing for a common parameter assignment. Nonetheless, implementing such simplification necessitates a deeper mechanistic understanding of chromosome-nuclear landmark interactions, an aspect currently lacking.

As our comprehension of nuclear organization progresses, the interpretability of parameter counts may improve, facilitating their reduction.

Comment 2: What would the modification be if the resolution is increased?

To increase the resolution of chromatin, we can in principle keep the same energy function as defined in Eq. S6. In this case, we only need to carry out further parameter optimization.

However, transitioning to higher resolutions may unveil additional features not readily apparent at 100kb. Notably, chromatin loops with an average size of 200kb or smaller have been identified in high-resolution Hi-C data [1]. To effectively capture these loops, new terms in the energy function must be incorporated. For instance, Qi and Zhang [2] employed additional contact potentials between CTCF sites to account for loop formation. Alternatively, an explicit loop-extrusion process could be introduced to model loop formation more accurately.

Comment 3: They should state that the extracted physical values are scale-dependent. For example, viscosity.

We thank the reviewer for the comment and would like to clarify that our model does not predict the viscosity. The nucleoplasmic viscosity was set as 1Pa · s to produce a diffusion coefficient that reproduces experimental value. The exact value for the nucleoplasmic viscosity is still rather controversial, and our selected value falls in the range of reported experimental values from 10−1Pa·s to 102Pa · s.

We have modified the main text to clarify the calculation of the diffusion coefficient.

“The exponent and the diffusion coefficient Dα = (27±11)×10−4μm2 · s−α both match well with the experimental values [cite], upon setting the nucleoplasmic viscosity as 1Pa · s (see Supporting Information Section: Mapping the reduced time unit to real time for more details).”

Reviewer 2:

Comment 0: In this work, Lao et al. develop an open-source software (OpenNucleome) for GPU-accelerated molecular dynamics simulation of the human nucleus accounting for chromatin, nucleoli, nuclear speckles, etc. Using this, the authors investigate the steady-state organization and dynamics of many of the nuclear components.

We thank the reviewer for summary of our work.

Comment 1: The authors could introduce a table having every parameter and the optimal parameter value used. This would greatly help the reader.

We would like to point out that model parameters are indeed provided in Table S1, S2, S3, S4, and Fig. S7. In these tables, we further provided details on how the parameters were determined.

Given the large number of parameters for the ideal potential (1000), we opted to plot it rather than listing out all the numbers. We added three new figures to plot the interaction parameters between chromosomes, between chromosomes and speckles, and between chromosomes and the nuclear lamina. Numerical values can be found online in the GitHub repository (parameters).

Comment 2: How many total beads are simulated? Do all beads have the same size?

The total number of the coarse-grained beads is 70542, including 60642 chromatin beads, 300 nucleolus beads, 1600 speckle beads, and 8000 nuclear lamina beads. The radius of the chromatin, nucleolus, and speckle beads is 0.25, while that of the lamina bead is 0.5. More information of the size and number of the beads are discussed in the Section: Components of the whole nucleus model.

Comment 3: In Equation S17, what is the 3rd and 4th powers mean? What necessitates it?

The potential defined in Equation S17 follows the definition of class2 bond in the LAMMPS package (LAMMPS docs). Compared to a typical harmonic potential, the presence of higher order terms produces sharper increase in the energy at large distances (Author response image 1). This essentially reduces the flucatuation of bond length in simulations.

Author response image 1.

Comparison between the Class2 potential (defined in Eq. S17) and the Harmonic potential (K(r − r0)2, with K = 20 and r0 = 0.5).

Comment 4: What do the X-axis and Y-axis numbers in Figure 5A and 5B mean? What are their units?

We apologize for the lack of clarify in our original figure. In Fig. 5A, the X and Y axis depicts the simulated and experimental radius of gyration (Rg) for individual chromosomes, as indicated in the title of the figure. Similarly, in Fig. 5B, the X and Y axis depicts the simulated and experimental radial position of individual chromosomes.

We have converted the chromosome Rg values into reduced units and labeled the corresponding axes in the updated figure (Fig. 5). The normalized radial position is unitless and its detailed definition is included in the supporting information Section: Computing simulated normalized chromosome radial positions. We updated the figure caption to provide an explicit reference to the SI text.

Reviewer 3:

Comment 0: In this work, the authors present the development of OpenNucleome, a software for simulating the structure and dynamics of the human nucleus. It provides a detailed model of nuclear components such as chromosomes and nuclear bodies, and uses GPU acceleration for better performance based on the OpenMM package. The work also shows the model’s accuracy in comparisons with experimental data and highlights the utility in the understanding of nuclear organization. While I consider this work a good tool for the genome architecture scientific community, I have some comments and questions that could further clarify the usage of this tool and help potential users. I also have a few questions that would help to clarify the technique and results and some suggestions for references.

We appreciate the reviewer’s strong assessment of the paper’s significance, novelty, and broad interest, and we thank them for the detailed suggestions and comments.

Comment 1: Could the authors elaborate on what they consider to be ’well-established and easily adoptable modeling tools’?

By well established, we meant that models that have been extensively validated and verified, and are highly regarded by the community.

By easily adoptable, we meant that tools that are well documented and can be relatively easily learned by new groups without help from the developers.

We have revised the text to clarify our meaning.

“Despite the progress made in computational modeling, the absence of well-documented software with easy-to-follow tutorials pose a challenge.”

Comment 2: Recognizing the value of a diverse range of tools in the community, the Open-MiChroM tool is also an open-source platform built on top of OpenMM. The documentation shows various modeling approaches and many tutorials that contain different approaches besides the MiChroM energy function. How does OpenNucleome compare in terms of facilitating crossvalidation and user accessibility? The two tools seem to be complementary, which is a gain to the field. I recommend adding one or two sentences in the matter. Also, while navigating the OpenNucleome GitHub, I have not found the tutorials mentioned in the text. I also consider a barrier in the process of generating necessary input files. I would suggest expanding the tutorials and documentation to help potential users.

We thank the reviewer for the excellent comments. We agree that while many of the tutorials were included in the original package, they were not as clearly documented. We have revised them extensively to to now present:

• A tutorial for optimizing chromosome chromosome interactions.

• A tutorial for optimizing chromosome nuclear landmark interactions.

• A tutorial for building initial configurations.

• A tutorial for relaxing the initial configurations.

• A tutorial for selecting the initial configurations.

• A tutorial for setting up performing Langevin dynamics simulations.

• A tutorial for setting up performing Brownian dynamics simulations.

• A tutorial for setting up performing simulations with deformed nucleus.

• A tutorial for analyzing simulation trajectories.

• A tutorial for introducing new features to the model.

These tutorials and our well-documented and open source code (https://zhanggroup-mitchemistry.github.io/OpenNucleome) should significantly promote user accessibility. Our inclusion of python scripts for analyzing simulation trajectorials shall allow users to compute various quantities for evaluating and comparing model quality.

We added a new paragraph in the Section: Conclusions and Dicussion of the main text to compare OpenNucleosome with existing software for genome modeling.

“Our software enhances the capabilities of existing genome simulation tools [cite]. Specifically, OpenNucleome aligns with the design principles of Open-MiChroM [cite], prioritizing open-source accessibility while expanding simulation capabilities to the entire nucleus. Similar to software from the Alber lab [cite], OpenNucleome offers highresolution genome organization that faithfully reproduces a diverse range of experimental data. Furthermore, beyond static structures, OpenNucleome facilitates dynamic simulations with explicit representations of various nuclear condensates, akin to the model developed by [citet].”

Comment 3: Lastly, I would appreciate it if the authors could expand their definition of ’standardized practices’.

We apologize for any confusion caused. By ”standardized practices,” we refer to the fact that different groups often employ unique procedures for structural modeling. These procedures differ in the representation of chromosomes, the nucleus environment, and the algorithms for parameter optimization. This absence of a consensus on the optimal practices for genome modeling can be daunting for newcomers to the field.

We have revised the text to the following to avoid confusion:

“Many research groups develop their own independent software, which complicates crossvalidation and hinders the establishment of best practices for genome modeling [3–5].”

Comment 4: On page 7, the authors refer to the SI Section: Components of the whole nucleus model for further details. Could the authors provide more information on the simulated density of nuclear bodies? Is there experimental data available that details the ratio of chromatin to other nuclear components, which was used as a reference in the simulation?

We thank the reviewer for the comment. Imaging studies have provided quantitative measures about the size and number of various nuclear bodies. For example, there are 2 ∼ 5 nucleoli per nucleus, with the typical size RNo ≈ 0.5μm [6–10]. In the review by Spector and Lamond [11], the authors showed that there are 20 ∼ 50 speckles, with the typical size RSp ≈ 0.3μm. We used these numbers to guide our simulation of nuclear bodies. These information was mentioned in the Section: Chromosomes as beads on the string polymers of the supporting information.

The chromatin density is fixed by the average size of chromatin bead and the nucleus size. We chose the size of chromatin based on imaging studies as detailed in the Subsection: Mapping chromatin bead size to real unit of the supporting information. Upon fixing the bead size, the chromatin volume is determined.

Comment 5: In the statement, ’the ideal potential is only applied for beads from the same chromosome to approximate the effect of loop extrusion by Cohesin molecules for chromosome compaction and territory formation,’ it would be helpful if the authors could clarify the scope of this potential. Specifically, the code indicates that the variable ’dend ideal’ is set at 1000, suggesting an interaction along a 100Mb polymer chain at a resolution of 100Kb per bead. Could the authors elaborate on their motivation for the Cohesin complex’s activity having a significant effect over such long distances within the polymer chain?

We thank the reviewer for the insight comment. They are correct that the ideal potential was introduced to capture chromosome folding beyond the interactions between compartments, including loop extrusion. Practically, we parameterized the ideal potential such that the simulated average contact probabilities as a function of sequence separation match the experimental values. The reviewer is correct that beyond a specific value of sequence separation, one would expect the impact of loop extrusion on chromosome folding should be negligible, due to Cohesin dissociation. Correspondingly, the interaction potential should be zero at large sequence separations.

However, it is important to note that the precise separation scale cannot be known a priori. We chose 100Mb as a conservative estimation. However, as we can see from Fig. S7, our parameterization scheme indeed produced interaction parameters are mainly zero at large sequence separations. Interesting, the scale at which the potential approaches 0 (∼ 500KB), indeed agree with the estimated length traveled by Cohesin molecules before dissociation [12].

Comment 6: On pages 8 and 9, the authors discuss the optimization process. However, in reviewing the code and documentation available on the GitHub page, I could not find specific sections related to the optimization procedure described in the paper. In this context, I have a few questions: Could the authors provide more details or direct me to the parts of the documentation and the text/SI that address the optimization procedure used in their study? Additional clarification on the cost/objective function employed during the optimization process would be highly beneficial, as this was not readily apparent in the text.

We thank the reviewer for the comment. We revised the SI to include the definition of the cost function for the Adam optimizer.

“During the optimization process, our aim was to minimize the disparity between experimental findings and simulated data. To achieve this, we defined the cost function as follows:

where the index i iterates over all the constraints defined in Eq. S28.”

The detailed optimization procedure was included in the SI as quoted below

“The details of the algorithm for parameter optimization are as follows

(1) Starting with a set of values for and we performed 50 independent 3-million-step long MD simulations to obtain an ensemble of nuclear configurations. The 500K steps of each trajectory are discarded

as equilibration. We collected the configurations at every 2000 simulation steps from the rest of the simulation trajectories to compute the ensemble averages defined on the left-hand side of Eq. S13.

(2) Check the convergence of the optimization by calculating the percentage of error

defined as . The summation over i includes all the average contact probabilities defined in Eq. S28.

(3) If the error is less than a tolerance value etol, the optimization has converged, and we stop the simulations. Otherwise, we update the parameters, α, using the Adam optimizer [13]. With the new parameter values, we return to step one and restart the iteration.”

Previously, the optimization code was included as part of the analysis folder. To avoid confusion and improve readability, a separate folder named optimization has been created. This folder provides the Adam optimization of chromosome-chromosome interactions (chr-chr optimization) and chromosome-nuclear landmarks interactions (chr-NL optimization).

Comment 7: What was the motivation for choosing the Adam algorithm for optimization? Adam is designed for training on stochastic objective functions. Could the authors elucidate on the ’stochastic’ aspect of their function to be optimized? Why the Adam algorithm was considered the most appropriate choice for this application?

We thank the reviewer for the comment. As defined in Eq. R1, the cost function measures the difference between the simulated constraints with corresponding experimental values. The estimation of simulation values, by averaging over an ensemble of chromosome configurations, is inherently noisy and stochastic. Exact ensemble averages can only be achieved with unlimited samples obtained from infinite long simulations.

In the past, we have used the Newton’s method for parameterization, and the detailed algorithm can be found in the SI of Ref. 14. However, we found that Adam is more efficient as it is a first-order approximation method. The Newton’s method, on the other hand, is second-order approximation method and requires estimation of the Hessian matrix. When the number of constraints is large, as is in our case, the computational cost for estimating the Hessian matrix can be significant. Another advantage of the Adam algorithm lies in its adjustment of the learning rate along the optimization to further speedup convergence.

Comment 8: The authors mention that examples of setting up simulations, parameter optimization, and introducing new features are provided in the GitHub repository. However, I was unable to locate these examples. Could the authors guide me to these specific resources or consider adding them if they are not currently available?

We thank the reviewer for the comment. We have improved the GitHub repository and all the tutorials can be found using the links provided in Response to Comment 2.

Comment 9: Furthermore, the paper states that ’a configuration file that provides the position of individual particles in the PDB file format is needed to initialize the simulations.’ It would be beneficial for new users if the authors could elaborate on how this file is generated. And all other input files in general. Detailing the procedures for a new user to run their system using OpenNucleome would be helpful.

We thank the reviewer for the comment. The procedure for generating initial configurations was explained in the SI Section: Initial configurations for simulations and quoted below.

“We first created a total of 1000 configurations for the genome by sequentially generating the conformation of each one of the 46 chromosomes as follows. For a given chromosome, we start by placing the first bead at the center (origin) of the nucleus. The positions of the following beads, i, were determined from the (i − 1)-th bead as . v is a normalized random vector, and 0.5 was selected as the bond length between neighboring beads. To produce globular chromosome conformations, we rejected vectors, v, that led to bead positions with distance from the center larger than 4σ. Upon creating the conformation of a chromosome i, we shift its center of mass to a value ri com determined as follows. We first compute a mean radial distance, with the following equation

where Di is the average value of Lamin B DamID profile for chromosome i. Dhi and Dlo represent the highest and lowest average DamID values of all chromosomes, and 6σ and 2σ represent the upper and lower bound in radial positions for chromosomes. As shown in Fig. S6, the average Lamin B DamID profiles are highly correlated with normalized chromosome radial positions as reported by DNA MERFISH [cite], supporting their use as a proxy for estimating normalized chromosome radial positions. We then select as a uniformly distributed random variable within the range . Without loss of generality, we randomly chose the directions for shifting all 46 chromosomes.

We further relaxed the 1000 configurations to build more realistic genome structures. Following an energy minimization process, one-million-step molecular dynamics (MD) simulations were performed starting from each configuration. Simulations were performed with the following energy function

where UGenome is defined as in Eq. S7. UG-La is the excluded volume potential between chromosomes and lamina, i.e, only the second term in Eq. S24. Parameters in UGenome were from a preliminary optimization. The end configurations of the MD simulations were collected to build the final configuration ensemble (FCE).”

The tutorial for preparing initial configurations can be found at this link.

Comment 10: In the section discussing the correlation between simulated and experimental contact maps, as referenced in Figure 4A and Figure S2, the authors mention a high degree of correlation. Could the authors specify the exact value of this correlation and explain the method used for its computation? Considering that comparing two Hi-C matrices involves a large number of data points, it would be helpful to know if all data points were included in this analysis.

We have updated Fig 4A and S2 to include Pearson correlation coefficients next to the contact maps. The reviewer is correct in that all the non-redundant data points of the contact maps are included in computing the correlation coefficients.

For improved clarity, we added a new section in the supporting information to detail the calculations. The section is titled Computing Pearson correlation coefficients between experimental and simulated contact maps, and the relevant text is quoted below.

“We computed the Pearson correlation coefficients (PCC) between experimental and simulated contact maps in Fig. 4A and Fig. S2 as

xi and yi represent the experimental and simulated contact probabilities, and n is the total number of data points. Only non-redundant data points, i.e., half of the pairwise contacts, are used in the PCC calculation.”

Comment 11: In addition, the author said: ”Moreover, the simulated and experimental average contact probabilities between pairs of chromosomes agree well, and the Pearson correlation coefficient between the two datasets reaches 0.89.” How does this correlation behave when not accounting for polymer compaction or scaling? An analysis presenting the correlation as a function of genomic distance would be interesting.

Author response image 2.

Pearson correlation coefficient between experimental and simulated contact probabilities as a function of the sequence separation within specific chromosomes. For each chromosome, we first gathered a set of experimental contacts alongside a matching set of simulated ones for genomic pairs within a particular separation range. The Pearson correlation coefficient at the corresponding sequence separation was then determined using Equation R4. We limited the calculations to half of the chromosome length to ensure the availability of sufficient data.

We thank the reviewer for the comment. The analysis presenting the correlation as a function of genomic distance (sequence separation) for each chromosome is shown in Figure S12 and also included in the SI. While the correlation coefficients decreases at larger separation, the values around 0.5 is quite reasonable and comparable to results obtained using Open-Michrom.

We also computed the correlation of whole genome contact maps after excluding intra-chromosomal contacts. The PCC decreased from 0.89 to 0.4. Again, the correlation coefficient is quite reasonable considering that these contacts are purely predicted by the compartmental interactions and were not directly optimized.

Comment 12: I recommend using the web-server that is familiar to the authors to benchmark the OpenNucleome tool/model: ”3DGenBench: A Web-Server to Benchmark Computational Models for 3D Genomics.” Nucleic Acids Research, vol. 50, no. W1, July 2022, pp. W4-12.

We appreciate the reviewer’s suggestion. Unfortunately, the website is no longer active during the time of the revision. However, as detailed in Response to comment 11, we used the one of the popular metrics to exclude polymer compact effect and evaluate the agreement between simulation and experiments.

Comment 13: Regarding the comparison of simulation results with microscopy data from reference 34. Given their different resolutions and data point/space groupings, how do the authors align these datasets? Could the authors describe how they performed this comparison? How were the radial positions calculated in both the simulations and experiments? Since the data from reference 34 indicates a non-globular shape of the nucleus; how did this factor into the calculation of radial distributions?

We thank the reviewer for the comment and apologize for the confusion. First, the average properties we examined, including radial positions and interchromosomal contacts, were averaged over all genomic loci. Therefore, they are independent of data resolution.

Secondly, instead of calculating the absolute radial positions, which are subject to variations in nucleus shape and size, we defined the normalized radial positions. They measure the ratio between the distance from the nucleus center to the chromosome center and the distance from the nucleus center to the lamina. This definition was frequently used in prior imaging studies to measure chromosome radial positions.

The calculation of the simulated normalized radial positions and the experimental normalized radial positions are discussed in the Section: Computing simulated normalized chromosome radial positions

“For a given chromosome i, we first determined its center of mass position denoted as Ci. Starting from the center of the nucleus, O, we extend the the vector vOC to identify the intersection point with the nuclear lamina as Pi. The normalized chromosome radial position i is then defined as , where ||·|| represents the L2 norm.

and Section: Computing experimental normalized chromosome radial positions.

“We followed the same procedure outlined in Section: Computing simulated normalized chromosome radial positions to compute the experimental values. To determine the center of the nucleus using DNA MERFISH data, we used the algorithm, minimum volume enclosing ellipsoid (MVEE)[15], to fit an ellipsoid for each genome structure. The optimal ellipsoid defined as is obtained by optimizing subjecting to the constraint that . xi correspond to the list of chromatin positions determined experimentally.”

Comment 14: In the sentence: ”It is evident that telomeres exhibit anomalous subdiffusive motion.” I recommend mentioning the work ”Di Pierro, Michele, et al., ”Anomalous Diffusion, Spatial Coherence, and Viscoelasticity from the Energy Landscape of Human Chromosomes.” Proceedings of the National Academy of Sciences, vol. 115, no. 30, July 2018, pp. 7753-58.”.

We have revised the sentence to include the citation as follows.

“In line with previous research [cite], telomeres display anomalous subdiffusive motion. When fitted with the equation , these trajectories yield a spectrum of α values, with a peak around 0.59.”

Comment 15: Regarding the observation that ’chromosomes appear arrested and no significant changes in their radial positions are observed over timescales comparable to the cell cycle,’ could the authors provide more details on the calculations or analyses that led to this conclusion? Specifically, information on the equilibration/relaxation time of chromosome territories relative to rearrangements within a cell cycle would be interesting.

Our conclusion here was mostly based on the time trace of normalized radial positions shown in Figure 6A of the main text. Over the timescale of an entire cell cycle (24 hours), the relatively little to no changes in the radial positions supports glassy dynamics of chromosomes. We further determined the mean squared displacement (MSD) for chromosome center of masses. As shown in the left panel of Fig. S12, the MSDs are much smaller than the average size of chromosomes (see Rg values in Fig. 5A), supporting arrested dynamics.

We further computed the auto-correlation function of the normalized chromosome radial position as

where t indexes over the trajectory frames and ¯r is the mean position. As shown in Fig. S12, the positions are not completely decorrelated over 10 hours, again supporting slow dynamics. It would be interesting to examine the relaxation timescale more closely in future studies.

Comment 16: The authors also comment on the SI ”Section: Initial configurations for simulations provides more details on preparing the 1000 initial configurations.” and related to reference 34 mentioning that ”the average Lamin B DamID profiles are highly correlated with chromosome radial positions as reported by DNA MERFISH”. How do the authors account for situations where homologous chromosomes are neighbors or have an interacting interface? Ref. 34 indicates that distinguishing between these scenarios can be challenging, potentially leading to ’invalid distributions’ that are filtered out. Clarification on how such cases were handled in the simulations would be helpful.

We would like to first clarify that when comparing with experimental data, we averaged over the homologous chromosomes to obtain haploid data. We added the following text in the manuscript to emphasize this point

“Given that the majority of experimental data were analyzed for the haploid genome, we adopted a similar approach by averaging over paternal and maternal chromosomes to facilitate direct comparison. More details on data analysis can be found in the Supporting Information Section: Details of simulation data analysis.”

Furthermore, we used the processed DNA MERFISH data from the Zhuang lab, which unambiguously assigns a chromosome ID to each data point. Therefore, the issue mentioned by the reviewer is not present in the procssed data. In our simulations, since we keep track of the explicit connection between genomic segments, the trace of individual chromosomes can be determined for any configuration. Therefore, there is no ambiguity in terms of simulation data.

Comment 17: When discussing the interaction with nuclear lamina and nuclear envelop deformation, I suggest mentioning the following studies: The already cited ref 52 and ”Contessoto, Vin´ıcius G., et al. ”Interphase Chromosomes of the Aedes Aegypti Mosquito Are Liquid Crystalline and Can Sense Mechanical Cues.” Nature Communications, vol. 14, no. 1, Jan. 2023, p. 326.”

We updated the text to include the suggested reference.

“Numerous studies have highlighted the remarkable influence of nuclear shape on the positioning of chromosomes and the regulation of gene expression [16, 17].”

Comment 18: The authors state that ’Tutorials in the format of Python Scripts with extensive documentation are provided to facilitate the adoption of the model by the community.’ However, as I mentioned, the documentation appears to be limited, and the available tutorials could benefit from further expansion. I suggest that the authors consider enhancing these resources to better assist users in adopting and understanding the model.

As detailed in the Response to Comment 2, we have updated the GitHub repository to better document the included Jupyter notebooks and tutorials.

Comment 19: In the Methods section, the authors discuss using Langevin dynamics for certain simulations and Brownian dynamics for others. Could the authors provide more detailed reasoning behind the choice of these different dynamics for different aspects of the simulation? Furthermore, it would be insightful to know how the results might vary if only one of these dynamics was utilized throughout the study. Such clarification would help in understanding the implications of these methodological choices on the outcomes of the simulations.

We thank the reviewer for the comment. As detailed in the supporting information Section: Mapping the Reduced Time Unit to Real Time, the Brownian dynamics simulations provide a rigorous mapping to the biological timescale. By choosing a specific value for the nucleoplasmic viscosity, we determined the time unit in simulations as τ = 0.65s. With this time conversion, the simulated diffusion coefficients of telomeres match well with experimental values. Therefore, Brownian dynamics simulations are recommended for computing time dependent quantities and the large damping coefficients mimics the complex nuclear environment well.

On the other hand, the large damping coefficient slows down the configuration relaxation of the system significantly. For computing equilibrium statistical properties, it is useful to use a small coefficient and the Langevin integrator with large time steps to facilitate conformational relaxation.

References

[1] Rao, S. S.; Huntley, M. H.; Durand, N. C.; Stamenova, E. K.; Bochkov, I. D.; Robinson, J. T.; Sanborn, A. L.; Machol, I.; Omer, A. D.; Lander, E. S.; others A 3D map of the human genome at kilobase resolution reveals principles of chromatin looping. Cell 2014, 159, 1665–1680.

[2] Qi, Y.; Zhang, B. Predicting three-dimensional genome organization with chromatin states. PLoS computational biology 2019, 15, e1007024.

[3] Yildirim, A.; Hua, N.; Boninsegna, L.; Zhan, Y.; Polles, G.; Gong, K.; Hao, S.; Li, W.; Zhou, X. J.; Alber, F. Evaluating the role of the nuclear microenvironment in gene function by population-based modeling. Nature Structural & Molecular Biology 2023, 1–14.

[4] Junior, A. B. O.; Contessoto, V. G.; Mello, M. F.; Onuchic, J. N. A scalable computational approach for simulating complexes of multiple chromosomes. Journal of molecular biology 2021, 433, 166700.

[5] Fujishiro, S.; Sasai, M. Generation of dynamic three-dimensional genome structure through phase separation of chromatin. Proceedings of the National Academy of Sciences 2022, 119, e2109838119.

[6] Caragine, C. M.; Haley, S. C.; Zidovska, A. Nucleolar dynamics and interactions with nucleoplasm in living cells. Elife 2019, 8, e47533.

[7] Brangwynne, C. P.; Mitchison, T. J.; Hyman, A. A. Active liquid-like behavior of nucleoli determines their size and shape in Xenopus laevis oocytes. Proceedings of the National Academy of Sciences 2011, 108, 4334–4339.

[8] Farley, K. I.; Surovtseva, Y.; Merkel, J.; Baserga, S. J. Determinants of mammalian nucleolar architecture. Chromosoma 2015, 124, 323–331.

[9] Qi, Y.; Zhang, B. Chromatin network retards nucleoli coalescence. Nature Communications 2021, 12, 6824.

[10] Caragine, C. M.; Haley, S. C.; Zidovska, A. Surface fluctuations and coalescence of nucleolar droplets in the human cell nucleus. Physical review letters 2018, 121, 148101.

[11] Spector, D. L.; Lamond, A. I. Nuclear speckles. Cold Spring Harbor perspectives in biology 2011, 3, a000646.

[12] Banigan, E. J.; Mirny, L. A. Loop extrusion: theory meets single-molecule experiments. Current opinion in cell biology 2020, 64, 124–138.

[13] Kingma, D. P.; Ba, J. Adam: A method for stochastic optimization. arXiv preprint arXiv:1412.6980 2014,

[14] Zhang, B.; Wolynes, P. G. Topology, structures, and energy landscapes of human chromosomes. Proceedings of the National Academy of Sciences 2015, 112, 6062–6067.

[15] Moshtagh, N.; others Minimum volume enclosing ellipsoid. Convex optimization 2005, 111, 1–9.

[16] Brahmachari, S.; Contessoto, V. G.; Di Pierro, M.; Onuchic, J. N. Shaping the genome via lengthwise compaction, phase separation, and lamina adhesion. Nucleic Acids Res. 2022, 50, 1–14.

[17] Contessoto, V. G.; Dudchenko, O.; Aiden, E. L.; Wolynes, P. G.; Onuchic, J. N.; Di Pierro, M. Interphase chromosomes of the Aedes aegypti mosquito are liquid crystalline and can sense mechanical cues. Nature Communications 2023, 14, 326.

https://doi.org/10.7554/eLife.93223.2.sa4

OpenNucleome for high resolution nuclear structural and dynamical modeling

Significance of findings

Strength of evidence

Abstract

Introduction

Results

Non-equilibrium nucleus model at 100 KB resolution

Computer model of the human nucleus for structural and dynamical characterizations.

Optimization of model parameters with experimental data

Overview of the iterative algorithm for parameterizing the nucleus model with experimental data.

Molecular dynamics simulations with GPU acceleration

OpenNucleome facilitates GPU-accelerated simulations of the human nucleus.

Simulations reproduce and predict diverse experimental data

Simulated structures reproduce contact frequencies between chromosomes and between chromosomes and nuclear landmarks.

Structural and dynamical predictions of the nucleus model match results from microscopy imaging.

Heterogeneity and conserved features of nuclear organizations.

Heterogeneity and robustness of the simulated conformational ensemble

Nuclear deformation preserves chromosome-nuclear body contacts

Nuclear deformations influence genome organization while preserving chromatin-speckle contacts.

Conclusions and Discussion

Methods

Molecular dynamics simulation details

Experimental data processing and analysis

Acknowledgements

Competing interests

Data and materials availability

References

Article and author information

Author information

Zhuohan Lao

Kartik Kamat

Zhongling Jiang

Bin Zhang

Version history

Copyright

Peer review process

Editors

Be the first to read new articles from eLife