Multiple Chains: Chromatin Dynamics Simulations on Chromosome 10 and Chromosome 11 of GM12878 Cell Line

This tutorial should take between 5 to 15 minutes of reading and performing simulations.

The first step is to import the OpenMiChroM module

from OpenMiChroM.ChromDynamics import MiChroM
from OpenMiChroM.CndbTools import cndbTools

OpenMiChroM allows the simulation of multiples chromosomes. In this tutorial, the multiple chain simulation will be performed using the chromosomes 10 and 11 of the human GM12878 cell line. The system is generated based on the collapsed structure of each chromosome. This collapse step is similar to the one presented in the single chromosome simulation. For each chromosome, the collapse simulation should be performed.

Chromosome 10 Collapse Simulation

sim_chr10 = MiChroM(name="chr10", temperature=1.0, time_step=0.01)
sim_chr10.setup(platform="opencl")
sim_chr10.saveFolder('output_chr10_chr11')
Chrom10 = sim_chr10.create_springSpiral(ChromSeq='inputs/chr10_beads.txt')
sim_chr10.loadStructure(Chrom10, center=True)

sim_chr10.addFENEBonds(kfb=30.0)
sim_chr10.addAngles(ka=2.0)
sim_chr10.addRepulsiveSoftCore(Ecut=4.0)
sim_chr10.addFlatBottomHarmonic()

sim_chr10.addTypetoType(mu=3.22, rc = 1.78)
sim_chr10.addIdealChromosome(mu=3.22, rc = 1.78, dinit=3, dend=500)

print("Performing chr10 simulation...")
for _ in range(1000):
    sim_chr10.runSimBlock(500)

sim_chr10.saveStructure(filename="chr10" ,mode="ndb")
del sim_chr10

Chromosome 11 Collapse Simulation

sim_chr11 = MiChroM(name="chr11", temperature=1.0, time_step=0.01)
sim_chr11.setup(platform="opencl")
sim_chr11.saveFolder('output_chr10_chr11')
Chrom11 = sim_chr11.create_springSpiral(ChromSeq='inputs/chr11_beads.txt')
sim_chr11.loadStructure(Chrom11, center=True)

sim_chr11.addFENEBonds(kfb=30.0)
sim_chr11.addAngles(ka=2.0)
sim_chr11.addRepulsiveSoftCore(Ecut=4.0)
sim_chr11.addFlatBottomHarmonic()

sim_chr11.addTypetoType(mu=3.22, rc = 1.78)
sim_chr11.addIdealChromosome(mu=3.22, rc = 1.78, dinit=3, dend=500)

print("Perform chr11 simulation...")
for _ in range(1000):
    sim_chr11.runSimBlock(500)

sim_chr11.saveStructure(filename="chr11", mode="ndb")
del sim_chr11

Once the collapse simulations are done, the collapsed structures of each chromosome should be included in the same simulation system.

sim_chr10_chr11 = MiChroM(name="chr10_chr11", temperature=1.0, time_step=0.01)

sim_chr10_chr11.setup(platform="opencl")

sim_chr10_chr11.saveFolder('output_chr10_chr11')

The function loadNDB receives a list of files and saves the chromosome collapse structures in the variable Struc_chr10_chr11.

Struc_chr10_chr11 = sim_chr10_chr11.loadNDB(NDBfiles=['output_chr10_chr11/chr10_0_block1000.ndb',
                                   'output_chr10_chr11/chr11_0_block1000.ndb'])

Struc_chr10_chr11 contains the position of all beads for each chromosome. The variable chains shows information of each chromosome. (start,end,is_ring), is_ring=0 represents a open chromosome chain and is_ring=1 represents a circular polymer (used for simulating bacteria genome).

Before loading the structures in the simulation context, it is necessary to distribute the chromosome chains.

• This step is essential for randomizing the initial condition when simulating different replicas.

Struc_chr10_chr11 = sim_chr10_chr11.setFibPosition(Struc_chr10_chr11, dist=(1.5,3.0))

Load the chromosomes in the simulation context.

sim_chr10_chr11.loadStructure(Struc_chr10_chr11, center=True)

The initial distribution of the chromosome structures can be saved in .ndb file format. The file is stored in the path given in saveFolder.

• In the case of having multiple chains in the simulation context, the saveStructure function will save each structure separately, and each saved structure will have an index starting from 0 . In this tutorial, 0 is for chr10, and 1 is for chr11

sim_chr10_chr11.saveStructure(mode='ndb')

The next step is to add the force field in the simulation object sim_chr10_chr11.

In this tutorial, the forces are added in two sets:

MiChroM Homopolymer (Bonded) Potentials

sim_chr10_chr11.addFENEBonds(kfb=30.0)
sim_chr10_chr11.addAngles(ka=2.0)
sim_chr10_chr11.addRepulsiveSoftCore(Ecut=4.0)
sim_chr10_chr11.addFlatBottomHarmonic(n_rad=20)

MiChroM non-Bonded Potentials Differently from the single chromosome tutorial, here, each force is added separately

The addTypetoType interaction is independent of the chromosome chain and only depends on the loci types. This potential is added for both chains simultaneously.

sim_chr10_chr11.addTypetoType(mu=3.22, rc = 1.78)

On the other hand, in the Ideal Chromossome potential, the forces need to be added for each chromosome separately.

The function addMultiChainIC receives the chromosome chain and adds the IC potential.

• The chromosome chain information can be obtained from the function chains. In this case, the data is retrieved from the object sim_chr10_chr11.chains.

• The Ideal Chromosome potential was applied from the genomic distance d = 3 to d=500. The range values can be adjusted based on the chromosome length.

sim_chr10_chr11.addMultiChainIC(chains=sim_chr10_chr11.chains[0], mu=3.22, rc = 1.78, dinit=3, dend=500)
sim_chr10_chr11.addMultiChainIC(chains=sim_chr10_chr11.chains[1], mu=3.22, rc = 1.78, dinit=3, dend=500)

The simulation setup is complete!

To run the production simulation, it is necessary to initialize the .cndb file to save the chromatin dynamics trajectory.

The files will be saved in the output folder defined using the saveFolder function. Each chromosome will be stored separately. The same procedure is adopted when using the saveStructure function.

sim_chr10_chr11.initStorage('traj_chr10_chr11', mode='w')

Sets the parameters of the production simulations:

\(block = 5\times10^2\) \(n\_blocks = 2\times10^3\)

block = 5*10**2
n_blocks = 5*10**2

for _ in range(n_blocks):
    sim_chr10_chr11.runSimBlock(block)
    sim_chr10_chr11.saveStructure()

Once the simulation is performed, it is necessary to close the .cndb file to avoid losing the trajectory data.

sim_chr10_chr11.storage[0].close() #close the cndb file for chr10
sim_chr10_chr11.storage[1].close() #close the cndb file for chr11

To visualize the chromosome’s 3D structures in the standard visualization software for macromolecules, there are available scripts for converting the ndb/cndb file format to .pdb and .gro. For details, please check the Nucleome Data Bank.

The ndb plugin for visualizing the chromatin dynamics trajectories in VMD/Chimera/Pymol is under development.

sim_chr10_chr11.saveStructure(mode="ndb")