[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]

[IMP-users] BuildModel 2 BuildSystem



Dear Imp-Users,

I am trying to model a complex using IMP however, I can't use the older IMP on a server. My problem is that when I try to change to using BuildSystem, my simulation fails.

I attach, the original modeling.py and the new modeling.py with the two inputs.

Can you help me how to fix this?

Best regards,
Marcell Miski
Ph.D. Student
PPCU FIT BP HU
|directories|
|pdb_dir|./|
|fasta_dir|./|
|gmm_dir|./|

|topology_dictionary|
|molecule_name|color|fasta_fn|fasta_id|pdb_fn|chain|residue_range|pdb_offset|bead_size|em_residues_per_gaussian|rigid_body | super_rigid_body | chain_of_super_rigid_bodies |
|PSD-95|blue|  2xkx.fasta.txt|2XKX:A|  2xkx.pdb|A|   62,148|0|300|0|1|1||
|PSD-95|blue|  2xkx.fasta.txt|2XKX:A|  BEADS   |A|  149,156|0| 10|0|0|1||
|PSD-95|blue| 2xkx.fasta.txt|2XKX:A|  2xkx.pdb|A|  157,243|0|300|0|2|1||
|PSD-95|blue|  2xkx.fasta.txt|2XKX:A|  BEADS   |A|  244,309|0| 10|0|0|1||
|PSD-95|blue|  2xkx.fasta.txt|2XKX:A|  2xkx.pdb|A|  310,390|0|300|0|3|1||
|PSD-95|blue| 2xkx.fasta.txt|2XKX:A|  BEADS   |A|  391,424|0| 10|0|0|1||
|PSD-95|blue| 2xkx.fasta.txt|2XKX:A|  2xkx.pdb|A|  425,495|0|300|0|4|1||
|PSD-95|blue|  2xkx.fasta.txt|2XKX:A|  BEADS   |A|  496,530|0| 10|0|0|1||
|PSD-95|blue|  2xkx.fasta.txt|2XKX:A|  2xkx.pdb|A|  531,590|0|300|0|5|1||
|Kir2.1-I|red|  2xky.fasta.txt|2XKY:I|  2xky.pdb|I|   69,309|0|300|0|6|2||
|Kir2.1-J|red|  2xky.fasta.txt|2XKY:J|  2xky.pdb|J|   69,309|0|300|0|7|2||
|Kir2.1-K|red|  2xky.fasta.txt|2XKY:K|  2xky.pdb|K|   69,309|0|300|0|8|2||
|Kir2.1-L|red|  2xky.fasta.txt|2XKY:L|  2xky.pdb|L|   69,309|0|300|0|9|2||
|Kir2.1-2-I|green|  2xky.fasta.txt|2XKY:I|  2xky_2_6.pdb|I|   69,309|0|300|0|10|3||
|Kir2.1-2-J|green|  2xky.fasta.txt|2XKY:J|  2xky_2_6.pdb|J|   69,309|0|300|0|11|3||
|Kir2.1-2-K|green|  2xky.fasta.txt|2XKY:K|  2xky_2_6.pdb|K|   69,309|0|300|0|12|3||
|Kir2.1-2-L|green|  2xky.fasta.txt|2XKY:L|  2xky_2_6.pdb|L|   69,309|0|300|0|13|3||
import datetime
import IMP
import IMP.core
import IMP.algebra
import IMP.atom
import IMP.container

import IMP.pmi.restraints.crosslinking
import IMP.pmi.restraints.stereochemistry
import IMP.pmi.restraints.em
import IMP.pmi.restraints.basic
import IMP.pmi.representation
import IMP.pmi.tools
import IMP.pmi.samplers
import IMP.pmi.output
import IMP.pmi.macros
import IMP.pmi.topology

import os
import sys

datadirectory = "../data/"
topology_file = datadirectory+"topology.txt"


num_frames = 500
num_mc_steps = 50
bead_max_trans = 5.00
rb_max_trans = 5.00	
rb_max_rot = 0

#van szerkezet es mozog
#homer coiled coil egy rigid body h egyutt mozogjon (3)
rigid_bodies = [["PDZ_1"],["PDZ_2"],["PDZ_3"],["SH3"],["GK"],
				]
#egy molekulaban levo domainek
#homer tobb molekulat akarunk egyutt mogatni, ezeket lanconkent deklaralhattuk, szal itt a lancokat egynek kell venni			
super_rigid_bodies = [ ["PDZ_1","PDZ_2","PDZ_3","SH3","GK"],
					  ]

starttime = datetime.datetime.now()	
				  
m = IMP.Model()
topology = IMP.pmi.topology.TopologyReader(topology_file)
domains = topology.get_components()

bm = IMP.pmi.macros.BuildModel(m,
                    component_topologies=domains,
                    list_of_rigid_bodies=rigid_bodies,
					list_of_super_rigid_bodies=super_rigid_bodies)

					
representation = bm.get_representation()

for nc,component in enumerate(domains):
    name = component.name
    sel = IMP.atom.Selection(representation.prot,molecule=name)
    ps = sel.get_selected_particles()
    clr = IMP.display.get_rgb_color(float(nc)/len(domains))
    for p in ps:
        if not IMP.display.Colored.get_is_setup(p):
            IMP.display.Colored.setup_particle(p,clr)
        else:
            IMP.display.Colored(p).set_color(clr)



#representation.shuffle_configuration(50)

representation.set_rigid_bodies_max_rot(rb_max_rot)
representation.set_floppy_bodies_max_trans(bead_max_trans)
representation.set_rigid_bodies_max_trans(rb_max_trans)


outputobjects = [] # reporter objects (for stat files)
sampleobjects = [] # sampling objects





# add colors to the components
for nc,component in enumerate(domains):
    name = component.name
    sel = IMP.atom.Selection(representation.prot,molecule=name)
    ps = sel.get_selected_particles()
    clr = IMP.display.get_rgb_color(float(nc)/len(domains))
    for p in ps:
        if not IMP.display.Colored.get_is_setup(p):
            IMP.display.Colored.setup_particle(p,clr)
        else:
            IMP.display.Colored(p).set_color(clr)




outputobjects = [] # reporter objects (for stat files)
sampleobjects = [] # sampling objects

outputobjects.append(representation)
sampleobjects.append(representation)


ev = IMP.pmi.restraints.stereochemistry.ExcludedVolumeSphere(representation, 
															resolution=8,
															kappa=1.0)
ev.add_to_model()
outputobjects.append(ev)


# Crosslinks - dataset 1 fejlecei az interactions.csvnek
columnmap={}
columnmap["Protein1"]="prot1"
columnmap["Protein2"]="prot2"
columnmap["Residue1"]="res1"
columnmap["Residue2"]="res2"
columnmap["IDScore"]=None


xl1 = IMP.pmi.restraints.crosslinking.ISDCrossLinkMS(representation,
                                   datadirectory+'interactions.csv', #most interactions.csv
                                   length=15.0, #interakcio kornyezete
                                   slope=0.05, #kozeledes gyorsasaga
                                   columnmapping=columnmap,
                                   resolution=10.0, #beam
                                   label="Chen", #interaction type?
                                   csvfile=True)
xl1.add_to_model()

sampleobjects.append(xl1)
outputobjects.append(xl1)


#floppy_bodies flexibilis szerkezetnelkuli csak szekvencia van
representation.optimize_floppy_bodies(10)


#--------------------------
# Monte-Carlo Sampling
#--------------------------
# This object defines all components to be sampled as well as the sampling protocol
mc1=IMP.pmi.macros.ReplicaExchange0(m,
                                    representation,
                                    monte_carlo_sample_objects=sampleobjects,
                                    output_objects=outputobjects,
									crosslink_restraints=[xl1], #chimera kimenethez csv file nev
                                    monte_carlo_temperature=1.0,
                                    simulated_annealing=True,
                                    simulated_annealing_minimum_temperature=1.0,
                                    simulated_annealing_maximum_temperature=2.5,
                                    simulated_annealing_minimum_temperature_nframes=200,
                                    simulated_annealing_maximum_temperature_nframes=20,
                                    replica_exchange_minimum_temperature=1.0,
                                    replica_exchange_maximum_temperature=2.5,
                                    number_of_best_scoring_models=20,
                                    monte_carlo_steps=num_mc_steps,
                                    number_of_frames=num_frames,
                                    global_output_directory="output")

# Start Sampling
mc1.execute_macro()

endtime = datetime.datetime.now()

file = open("simulation_data.txt","w") 


file.write("\n")
file.write("Start: ")
file.write(str(starttime))
file.write("End: ")
file.write(str(endtime))
file.write("Elapsed: ")
file.write(str(endtime-starttime))

file.close()

print("\n")
print("Start: ")
print(starttime)
print("End: ")
print(endtime)
print("Elapsed: ")
print(endtime-starttime)
|directories|
|pdb_dir|./|
|fasta_dir|./|
|gmm_dir|./|

|topology_dictionary|
|component_name|domain_name|fasta_fn|fasta_id|pdb_fn|chain|residue_range|pdb_offset|bead_size|em_residues_per_gaussian|
|PSD-95|  PDZ_1|  2xkx.fasta.txt|2XKX:A|  2xkx.pdb|A|   62,148|0|300|0|
|PSD-95| LINK_1|  2xkx.fasta.txt|2XKX:A|  2xkx.pdb|A|  149,156|0| 10|0|
|PSD-95|  PDZ_2|  2xkx.fasta.txt|2XKX:A|  2xkx.pdb|A|  157,243|0|300|0|
|PSD-95| LINK_2|  2xkx.fasta.txt|2XKX:A|  2xkx.pdb|A|  244,309|0| 10|0|
|PSD-95|  PDZ_3|  2xkx.fasta.txt|2XKX:A|  2xkx.pdb|A|  310,390|0|300|0|
|PSD-95| LINK_3|  2xkx.fasta.txt|2XKX:A|  2xkx.pdb|A|  391,424|0| 10|0|
|PSD-95|    SH3|  2xkx.fasta.txt|2XKX:A|  2xkx.pdb|A|  425,495|0|300|0|
|PSD-95| LINK_4|  2xkx.fasta.txt|2XKX:A|  2xkx.pdb|A|  496,530|0| 10|0|
|PSD-95|     GK|  2xkx.fasta.txt|2XKX:A|  2xkx.pdb|A|  531,590|0|300|0|
| KIR_I| ANCH_I|  2xky.fasta.txt|2XKY:I|  2xky.pdb|I|   69,309|0|300|0|
| KIR_J| ANCH_J|  2xky.fasta.txt|2XKY:J|  2xky.pdb|J|   69,309|0|300|0|
| KIR_K| ANCH_K|  2xky.fasta.txt|2XKY:K|  2xky.pdb|K|   69,309|0|300|0|
| KIR_L| ANCH_L|  2xky.fasta.txt|2XKY:L|  2xky.pdb|L|   69,309|0|300|0|
| KIR_I_2| ANCH_I_2|  2xky.fasta.txt|2XKY:I|  2xky_2_6.pdb|I|   69,309|0|300|0|
| KIR_J_2| ANCH_J_2|  2xky.fasta.txt|2XKY:J|  2xky_2_6.pdb|J|   69,309|0|300|0|
| KIR_K_2| ANCH_K_2|  2xky.fasta.txt|2XKY:K|  2xky_2_6.pdb|K|   69,309|0|300|0|
| KIR_L_2| ANCH_L_2|  2xky.fasta.txt|2XKY:L|  2xky_2_6.pdb|L|   69,309|0|300|0|
prot1,res1,prot2,res2
PSD-95,100,KIR_I,308
PSD-95,200,KIR_J_2,308
import datetime
import IMP
import IMP.core
import IMP.algebra
import IMP.atom
import IMP.container

import IMP.pmi.restraints.crosslinking
import IMP.pmi.restraints.stereochemistry
import IMP.pmi.restraints.em
import IMP.pmi.restraints.basic
import IMP.pmi.representation
import IMP.pmi.tools
import IMP.pmi.samplers
import IMP.pmi.output
import IMP.pmi.macros
import IMP.pmi.topology

import os
import sys

#---------------------------
# Define Input Files
#---------------------------
datadirectory = "../data/"
topology_file = datadirectory+"topology2.txt"


#--------------------------
# Monte-Carlo Sampling
#--------------------------

#--------------------------
# Set MC Sampling Parameters
#--------------------------
num_frames = 20000
if '--test' in sys.argv: num_frames=100
num_mc_steps = 10

#---------------------------
# Define Input Degrees of Freedom
#---------------------------
bead_max_trans = 5.00
rb_max_trans = 5.00
rb_max_rot = 0


# Initialize model
m = IMP.Model()
# Read in the topology file.
# Specify the directory wheere the PDB files, fasta files and GMM files are
topology = IMP.pmi.topology.TopologyReader(topology_file,
                                  pdb_dir=datadirectory,
                                  fasta_dir=datadirectory,
                                  gmm_dir=datadirectory)
# Use the BuildSystem macro to build states from the topology file
bs = IMP.pmi.macros.BuildSystem(m)
# Each state can be specified by a topology file.
bs.add_state(topology)


starttime = datetime.datetime.now()

#Root-hierarchy and degrees of freedom
root_hier, dof = bs.execute_macro(max_rb_trans=rb_max_trans,
                                  max_rb_rot=rb_max_rot,
                                  max_bead_trans=bead_max_trans,
                                  max_srb_trans=4.0,
                                  max_srb_rot=0.3)

# Fix all rigid bodies but not Rpb4 and Rpb7 (the stalk)
# First select and gather all particles to fix.
fixed_particles=[]
for prot in ["Kir-I","Kir-J","Kir-K","Kir-L","Kir2-I","Kir2-J","Kir2-K","Kir2-L"]:
    fixed_particles+=IMP.atom.Selection(root_hier,molecule=prot).get_selected_particles()
# Fix the Corresponding Rigid movers and Super Rigid Body movers using dof
# The flexible beads will still be flexible (fixed_beads is an empty list)!
fixed_beads,fixed_rbs=dof.disable_movers(fixed_particles,
                                         [IMP.core.RigidBodyMover,
                                          IMP.pmi.TransformMover])
"""
# Randomize the initial configuration before sampling, of only the molecules
# we are interested in (Rpb4 and Rpb7)
IMP.pmi.tools.shuffle_configuration(root_hier,
                                    excluded_rigid_bodies=fixed_rbs,
                                    max_translation=50,
                                    verbose=False,
                                    cutoff=5.0,
                                    niterations=100)
"""
# Connectivity keeps things connected along the backbone (ignores if inside
# same rigid body)
outputobjects=[] # reporter objects (for stat files)
mols = IMP.pmi.tools.get_molecules(root_hier)
for mol in mols:
    molname=mol.get_name()
    IMP.pmi.tools.display_bonds(mol)
    cr = IMP.pmi.restraints.stereochemistry.ConnectivityRestraint(mol,scale=2.0)
    cr.add_to_model()
    cr.set_label(molname)
    outputobjects.append(cr)

ev = IMP.pmi.restraints.stereochemistry.ExcludedVolumeSphere(
                                         included_objects=root_hier,
                                         resolution=10)
ev.add_to_model()
outputobjects.append(ev)

#---------------------------
# Give Cross-Links
#---------------------------

# Crosslinks - dataset 2
#  We can easily add a second set of crosslinks.
#  These have a different format and label, but other settings are the same
sampleobjects=[]
xldbkwc = IMP.pmi.io.crosslink.CrossLinkDataBaseKeywordsConverter()
xldbkwc.set_protein1_key("prot1")
xldbkwc.set_protein2_key("prot2")
xldbkwc.set_residue1_key("res1")
xldbkwc.set_residue2_key("res2")

xl2db = IMP.pmi.io.crosslink.CrossLinkDataBase(xldbkwc)
xl2db.create_set_from_file(datadirectory+'interactions2.csv')

xl2 = IMP.pmi.restraints.crosslinking.CrossLinkingMassSpectrometryRestraint(
                                   root_hier=root_hier,
                                   CrossLinkDataBase=xl2db,
                                   length=21.0,
                                   slope=0.01,
                                   resolution=1.0,
                                   label="Chen",
                                   weight=1.)
xl2.add_to_model()
outputobjects.append(xl2)
sampleobjects.append(xl2)
#---------------------------
# Give EM restraints - NO
#---------------------------

#---------------------------
# SAMPLING
#---------------------------

mc1=IMP.pmi.macros.ReplicaExchange0(m,
                                    root_hier=root_hier,
                                    monte_carlo_sample_objects=sampleobjects,#dof.get_movers(),
                                    output_objects=outputobjects,
                                    crosslink_restraints=xl2,
                                    monte_carlo_temperature=1.0,
                                    simulated_annealing=True,
                                    simulated_annealing_minimum_temperature=1.0,
                                    simulated_annealing_maximum_temperature=2.5,
                                    simulated_annealing_minimum_temperature_nframes=200,
                                    simulated_annealing_maximum_temperature_nframes=20,
                                    replica_exchange_minimum_temperature=1.0,
                                    replica_exchange_maximum_temperature=2.5,
                                    number_of_best_scoring_models=100,
                                    monte_carlo_steps=num_mc_steps,
                                    number_of_frames=num_frames,
                                    global_output_directory="output")

mc1.execute_macro()

endtime = datetime.datetime.now()

file = open("simulation_data.txt","w")


file.write("\n")
file.write("Start: ")
file.write(str(starttime))
file.write("End: ")
file.write(str(endtime))
file.write("Elapsed: ")
file.write(str(endtime-starttime))

file.close()

print("\n")
print("Start: ")
print(starttime)
print("End: ")
print(endtime)
print("Elapsed: ")
print(endtime-starttime)
prot1,res1,prot2,res2
PSD-95,100,Kir-I,308
PSD-95,200,Kir2-J,308