3 """@namespace IMP.pmi.representation
4 Representation of the system.
7 from __future__
import print_function
20 from math
import pi, sqrt
21 from operator
import itemgetter
24 class Representation(object):
28 Set up the representation of all proteins and nucleic acid macromolecules.
30 Create the molecular hierarchies, representation,
31 sequence connectivity for the various involved proteins and
32 nucleic acid macromolecules:
34 Create a protein, DNA or RNA, represent it as a set of connected balls of appropriate
35 radii and number of residues, PDB at given resolution(s), or ideal helices.
37 How to use the SimplifiedModel class (typical use):
39 see test/test_hierarchy_contruction.py
43 1) Create a chain of helices and flexible parts
45 c_1_119 =self.add_component_necklace("prot1",1,119,20)
46 c_120_131 =self.add_component_ideal_helix("prot1",resolutions=[1,10],resrange=(120,131))
47 c_132_138 =self.add_component_beads("prot1",[(132,138)])
48 c_139_156 =self.add_component_ideal_helix("prot1",resolutions=[1,10],resrange=(139,156))
49 c_157_174 =self.add_component_beads("prot1",[(157,174)])
50 c_175_182 =self.add_component_ideal_helix("prot1",resolutions=[1,10],resrange=(175,182))
51 c_183_194 =self.add_component_beads("prot1",[(183,194)])
52 c_195_216 =self.add_component_ideal_helix("prot1",resolutions=[1,10],resrange=(195,216))
53 c_217_250 =self.add_component_beads("prot1",[(217,250)])
56 self.set_rigid_body_from_hierarchies(c_120_131)
57 self.set_rigid_body_from_hierarchies(c_139_156)
58 self.set_rigid_body_from_hierarchies(c_175_182)
59 self.set_rigid_body_from_hierarchies(c_195_216)
61 clist=[c_1_119,c_120_131,c_132_138,c_139_156,c_157_174,c_175_182,c_183_194,c_195_216,
64 self.set_chain_of_super_rigid_bodies(clist,2,3)
66 self.set_super_rigid_bodies(["prot1"])
70 def __init__(self, m, upperharmonic=True, disorderedlength=True):
73 @param upperharmonic This flag uses either harmonic (False)
74 or upperharmonic (True) in the intra-pair
75 connectivity restraint.
76 @param disorderedlength This flag uses either disordered length
77 calculated for random coil peptides (True) or zero
78 surface-to-surface distance between beads (False)
79 as optimal distance for the sequence connectivity
90 self.upperharmonic = upperharmonic
91 self.disorderedlength = disorderedlength
92 self.rigid_bodies = []
93 self.fixed_rigid_bodies = []
94 self.fixed_floppy_bodies = []
95 self.floppy_bodies = []
101 self.super_rigid_bodies = []
102 self.rigid_body_symmetries = []
103 self.output_level =
"low"
106 self.maxtrans_rb = 2.0
107 self.maxrot_rb = 0.04
108 self.maxtrans_srb = 2.0
109 self.maxrot_srb = 0.2
110 self.rigidbodiesarefixed =
False
111 self.floppybodiesarefixed =
False
112 self.maxtrans_fb = 3.0
113 self.resolution = 10.0
114 self.bblenght = 100.0
118 self.representation_is_modified =
False
119 self.unmodeledregions_cr_dict = {}
120 self.sortedsegments_cr_dict = {}
122 self.connected_intra_pairs = []
125 self.sequence_dict = {}
126 self.hier_geometry_pairs = {}
132 self.hier_representation = {}
133 self.hier_resolution = {}
136 self.reference_structures = {}
139 self.linker_restraints_dict = {}
140 self.threetoone = {
'ALA':
'A',
'ARG':
'R', 'ASN': 'N', 'ASP': 'D',
141 'CYS':
'C',
'GLU':
'E',
'GLN':
'Q',
'GLY':
'G',
142 'HIS':
'H',
'ILE':
'I',
'LEU':
'L',
'LYS':
'K',
143 'MET':
'M',
'PHE':
'F',
'PRO':
'P',
'SER':
'S',
144 'THR':
'T',
'TRP':
'W',
'TYR':
'Y',
'VAL':
'V',
'UNK':
'X'}
146 self.onetothree = dict((v, k)
for k, v
in self.threetoone.items())
150 def set_label(self, label):
153 def create_component(self, name, color=0.0):
155 protein_h.set_name(name)
156 self.hier_dict[name] = protein_h
157 self.hier_representation[name] = {}
158 self.hier_db.add_name(name)
159 self.prot.add_child(protein_h)
160 self.color_dict[name] = color
161 self.elements[name] = []
166 def add_component_name(self, *args, **kwargs):
167 self.create_component(*args, **kwargs)
169 def get_component_names(self):
170 return list(self.hier_dict.keys())
175 Add the primary sequence for a single component.
177 @param name Human-readable name of the component
178 @param filename Name of the FASTA file
179 @param id Identifier of the sequence in the FASTA file header
180 (if not provided, use `name` instead)
185 if id
not in record_dict:
186 raise KeyError(
"id %s not found in fasta file" % id)
187 length = len(record_dict[id])
188 self.sequence_dict[name] = str(record_dict[id])
191 self.sequence_dict[name]=offs_str+self.sequence_dict[name]
193 self.elements[name].append((length, length,
" ",
"end"))
195 def autobuild_model(self, name, pdbname, chain,
196 resolutions=
None, resrange=
None,
198 color=
None, pdbresrange=
None, offset=0,
199 show=
False, isnucleicacid=
False,
202 self.representation_is_modified =
True
206 color = self.color_dict[name]
208 self.color_dict[name] = color
210 if resolutions
is None:
212 print(
"autobuild_model: constructing %s from pdb %s and chain %s" % (name, pdbname, str(chain)))
221 t.get_children()[0].get_children()[0]).
get_index()
223 t.get_children()[0].get_children()[-1]).
get_index()
229 if name
in self.sequence_dict:
230 resrange = (1, len(self.sequence_dict[name]))
232 resrange = (start + offset, end + offset)
235 resrange = (resrange[0],end)
236 start = resrange[0] - offset
237 end = resrange[1] - offset
254 for n, g
in enumerate(gaps):
258 print(
"autobuild_model: constructing fragment %s from pdb" % (str((first, last))))
260 chain, resolutions=resolutions,
261 color=color, cacenters=
True,
262 resrange=(first, last),
263 offset=offset, isnucleicacid=isnucleicacid)
264 elif g[2] ==
"gap" and n > 0:
265 print(
"autobuild_model: constructing fragment %s as a bead" % (str((first, last))))
266 parts = self.hier_db.get_particles_at_closest_resolution(name,
271 first+offset, last+offset, missingbeadsize, incoord=xyz)
273 elif g[2] ==
"gap" and n == 0:
275 print(
"autobuild_model: constructing fragment %s as a bead" % (str((first, last))))
277 first+offset, last+offset, missingbeadsize, incoord=xyznter)
284 def autobuild_pdb_and_intervening_beads(self, *args, **kwargs):
285 r = self.autobuild_model(*args, **kwargs)
289 resrange=
None, offset=0, cacenters=
True, show=
False,
290 isnucleicacid=
False, readnonwateratoms=
False,
291 read_ca_cb_only=
False):
293 Add a component that has an associated 3D structure in a PDB file.
295 Reads the PDB, and constructs the fragments corresponding to contiguous
298 @return a list of hierarchies.
300 @param name (string) the name of the component
301 @param pdbname (string) the name of the PDB file
302 @param chain (string or integer) can be either a string (eg, "A")
303 or an integer (eg, 0 or 1) in case you want
304 to get the corresponding chain number in the PDB.
305 @param resolutions (integers) a list of integers that corresponds
306 to the resolutions that have to be generated
307 @param color (float from 0 to 1) the color applied to the
308 hierarchies generated
309 @param resrange (tuple of integers): the residue range to extract
310 from the PDB. It is a tuple (beg,end). If not specified,
311 all residues belonging to the specified chain are read.
312 @param offset (integer) specifies the residue index offset to be
313 applied when reading the PDB (the FASTA sequence is
314 assumed to start from residue 1, so use this parameter
315 if the PDB file does not start at residue 1)
316 @param cacenters (boolean) if True generates resolution=1 beads
317 centered on C-alpha atoms.
318 @param show (boolean) print out the molecular hierarchy at the end.
319 @param isnucleicacid (boolean) use True if you're reading a PDB
321 @param readnonwateratoms (boolean) if True fixes some pathological PDB.
322 @param read_ca_cb_only (boolean) if True, only reads CA/CB
325 self.representation_is_modified =
True
328 color = self.color_dict[name]
329 protein_h = self.hier_dict[name]
339 if type(chain) == str:
347 t.get_children()[0].get_children()[0]).
get_index()
349 t.get_children()[0].get_children()[-1]).
get_index()
350 c =
IMP.atom.Chain(IMP.atom.get_by_type(t, IMP.atom.CHAIN_TYPE)[0])
354 IMP.atom.get_by_type(t, IMP.atom.CHAIN_TYPE)[chain])
361 if not resrange
is None:
362 if resrange[0] > start:
364 if resrange[1] < end:
367 if not isnucleicacid:
371 residue_indexes=list(range(
374 atom_type=IMP.atom.AT_CA)
381 residue_indexes=list(range(
384 atom_type=IMP.atom.AT_P)
386 ps = sel.get_selected_particles()
388 raise ValueError(
"%s no residue found in pdb %s chain %s that overlaps with the queried stretch %s-%s" \
389 % (name, pdbname, str(chain), str(resrange[0]), str(resrange[1])))
397 start = start + offset
400 self.elements[name].append(
401 (start, end, pdbname.split(
"/")[-1] +
":" + chain,
"pdb"))
404 resolutions, isnucleicacid, c0, protein_h,
"pdb", color)
415 def add_component_ideal_helix(
423 self.representation_is_modified =
True
424 from math
import pi, cos, sin
426 protein_h = self.hier_dict[name]
429 color = self.color_dict[name]
433 self.elements[name].append((start, end,
" ",
"helix"))
435 for n, res
in enumerate(range(start, end + 1)):
436 if name
in self.sequence_dict:
438 rtstr = self.onetothree[
439 self.sequence_dict[name][res-1]]
459 x = 2.3 * cos(n * 2 * pi / 3.6)
460 y = 2.3 * sin(n * 2 * pi / 3.6)
461 z = 6.2 / 3.6 / 2 * n * 2 * pi / 3.6
470 resolutions,
False, c0, protein_h,
"helix", color)
479 """ add beads to the representation
480 @param name the component name
481 @param ds a list of tuples corresponding to the residue ranges
483 @param colors a list of colors associated to the beads
484 @param incoord the coordinate tuple corresponding to the position
489 self.representation_is_modified =
True
491 protein_h = self.hier_dict[name]
494 colors = [self.color_dict[name]]
497 for n, dss
in enumerate(ds):
498 ds_frag = (dss[0], dss[1])
499 self.elements[name].append((dss[0], dss[1],
" ",
"bead"))
501 if ds_frag[0] == ds_frag[1]:
503 if name
in self.sequence_dict:
505 rtstr = self.onetothree[
506 self.sequence_dict[name][ds_frag[0]-1]]
513 h.set_name(name +
'_%i_bead' % (ds_frag[0]))
514 prt.set_name(name +
'_%i_bead' % (ds_frag[0]))
518 h.set_name(name +
'_%i-%i_bead' % (ds_frag[0], ds_frag[1]))
519 prt.set_name(name +
'_%i-%i_bead' % (ds_frag[0], ds_frag[1]))
520 h.set_residue_indexes(list(range(ds_frag[0], ds_frag[1] + 1)))
521 resolution = len(h.get_residue_indexes())
522 if "Beads" not in self.hier_representation[name]:
524 root.set_name(
"Beads")
525 self.hier_representation[name][
"Beads"] = root
526 protein_h.add_child(root)
527 self.hier_representation[name][
"Beads"].add_child(h)
529 for kk
in range(ds_frag[0], ds_frag[1] + 1):
530 self.hier_db.add_particles(name, kk, resolution, [h])
553 ptem.set_radius(radius)
556 radius = 0.8 * (3.0 / 4.0 / pi * volume) ** (1.0 / 3.0)
558 ptem.set_radius(radius)
560 if not tuple(incoord)
is None:
561 ptem.set_coordinates(incoord)
566 self.floppy_bodies.append(prt)
575 Generates a string of beads with given length.
577 self.representation_is_modified =
True
585 [(chunk[0], chunk[-1])], colors=colors,incoord=incoord)
589 self, name, hierarchies=
None, selection_tuples=
None,
591 resolution=0.0, num_components=10,
592 inputfile=
None, outputfile=
None,
595 covariance_type=
'full', voxel_size=1.0,
597 sampled_points=1000000, num_iter=100,
599 multiply_by_total_mass=
True,
601 intermediate_map_fn=
None,
602 density_ps_to_copy=
None,
603 use_precomputed_gaussians=
False):
605 Sets up a Gaussian Mixture Model for this component.
606 Can specify input GMM file or it will be computed.
607 @param name component name
608 @param hierarchies set up GMM for some hierarchies
609 @param selection_tuples (list of tuples) example (first_residue,last_residue,component_name)
610 @param particles set up GMM for particles directly
611 @param resolution usual PMI resolution for selecting particles from the hierarchies
612 @param inputfile read the GMM from this file
613 @param outputfile fit and write the GMM to this file (text)
614 @param outputmap after fitting, create GMM density file (mrc)
615 @param kernel_type for creating the intermediate density (points are sampled to make GMM). Options are IMP.em.GAUSSIAN, IMP.em.SPHERE, and IMP.em.BINARIZED_SPHERE
616 @param covariance_type for fitting the GMM. options are 'full', 'diagonal' and 'spherical'
617 @param voxel_size for creating the intermediate density map and output map.
618 lower number increases accuracy but also rasterizing time grows
619 @param out_hier_name name of the output density hierarchy
620 @param sampled_points number of points to sample. more will increase accuracy and fitting time
621 @param num_iter num GMM iterations. more will increase accuracy and fitting time
622 @param multiply_by_total_mass multiply the weights of the GMM by this value (only works on creation!)
623 @param transform for input file only, apply a transformation (eg for multiple copies same GMM)
624 @param intermediate_map_fn for debugging, this will write the intermediate (simulated) map
625 @param density_ps_to_copy in case you already created the appropriate GMM (eg, for beads)
626 @param use_precomputed_gaussians Set this flag and pass fragments - will use roughly spherical Gaussian setup
634 self.representation_is_modified =
True
636 protein_h = self.hier_dict[name]
637 if "Densities" not in self.hier_representation[name]:
639 root.set_name(
"Densities")
640 self.hier_representation[name][
"Densities"] = root
641 protein_h.add_child(root)
644 fragment_particles = []
645 if not particles
is None:
646 fragment_particles += particles
647 if not hierarchies
is None:
649 self, resolution=resolution,
650 hierarchies=hierarchies)
651 if not selection_tuples
is None:
652 for st
in selection_tuples:
653 fragment_particles += IMP.pmi.tools.select_by_tuple(
654 self, tupleselection=st,
655 resolution=resolution,
656 name_is_ambiguous=
False)
659 density_particles = []
662 inputfile, density_particles,
664 elif density_ps_to_copy:
665 for ip
in density_ps_to_copy:
671 density_particles.append(p)
672 elif use_precomputed_gaussians:
673 if len(fragment_particles) == 0:
674 print(
"add_component_density: no particle was selected")
676 for p
in fragment_particles:
679 raise Exception(
"The particles you selected must be Fragments and XYZs")
681 pos=
IMP.core.XYZ(self.m,p.get_particle_index()).get_coordinates()
686 raise Exception(
"We haven't created a bead file for",nres,
"residues")
693 if len(fragment_particles) == 0:
694 print(
"add_component_density: no particle was selected")
697 density_particles = IMP.isd.gmm_tools.sample_and_fit_to_particles(
706 multiply_by_total_mass,
712 s0.set_name(out_hier_name)
713 self.hier_representation[name][
"Densities"].add_child(s0)
715 for nps, p
in enumerate(density_particles):
717 p.set_name(s0.get_name() +
'_gaussian_%i' % nps)
720 def get_component_density(self, name):
721 return self.hier_representation[name][
"Densities"]
724 selection_tuples=
None,
729 '''Decorates all specified particles as Gaussians directly.
730 @param name component name
731 @param hierarchies set up GMM for some hierarchies
732 @param selection_tuples (list of tuples) example (first_residue,last_residue,component_name)
733 @param particles set up GMM for particles directly
734 @param resolution usual PMI resolution for selecting particles from the hierarchies
740 from math
import sqrt
741 self.representation_is_modified =
True
743 if particles
is None:
744 fragment_particles = []
746 fragment_particles = particles
748 if not hierarchies
is None:
750 self, resolution=resolution,
751 hierarchies=hierarchies)
753 if not selection_tuples
is None:
754 for st
in selection_tuples:
755 fragment_particles += IMP.pmi.tools.select_by_tuple(
756 self, tupleselection=st,
757 resolution=resolution,
758 name_is_ambiguous=
False)
760 if len(fragment_particles) == 0:
761 print(
"add all atom densities: no particle was selected")
765 print(
'setting up all atom gaussians num_particles',len(fragment_particles))
766 for n,p
in enumerate(fragment_particles):
774 print(
'setting up particle',p.get_name(),
" as individual gaussian particle")
776 if not output_map
is None:
777 print(
'writing map to', output_map)
785 Make a copy of a hierarchy and append it to a component.
788 self.representation_is_modified =
True
789 protein_h = self.hier_dict[name]
790 hierclone = IMP.atom.create_clone(hierarchy)
791 hierclone.set_name(hierclone.get_name() +
"_clone")
792 protein_h.add_child(hierclone)
793 outhier.append(hierclone)
799 for n, pmain
in enumerate(psmain):
805 self.hier_db.add_particles(
821 def dump_particle_descriptors(self):
827 particles_attributes={}
828 floppy_body_attributes={}
834 hparent=h.get_parent()
835 while hparent != hroot:
836 hparent=h.get_parent()
837 name+=
"|"+hparent.get_name()
839 particles_attributes[name]={
"COORDINATES":numpy.array(
IMP.core.XYZR(leaf.get_particle()).get_coordinates()),
845 rigid_body_attributes={}
846 for rb
in self.rigid_bodies:
848 rf=rb.get_reference_frame()
849 t=rf.get_transformation_to()
850 trans=t.get_translation()
852 rigid_body_attributes[name]={
"TRANSLATION":numpy.array(trans),
853 "ROTATION":numpy.array(rot.get_quaternion()),
854 "COORDINATES_NONRIGID_MEMBER":{},
855 "COORDINATES_RIGID_MEMBER":{}}
856 for mi
in rb.get_member_indexes():
857 rm=self.m.get_particle(mi)
859 name_part=rm.get_name()
861 rigid_body_attributes[name][
"COORDINATES_NONRIGID_MEMBER"][name_part]=numpy.array(xyz)
863 name_part=rm.get_name()
865 rigid_body_attributes[name][
"COORDINATES_RIGID_MEMBER"][name_part]=numpy.array(xyz)
869 pickle.dump(particles_attributes,
870 open(
"particles_attributes.pkl",
"wb"))
871 pickle.dump(rigid_body_attributes,
872 open(
"rigid_body_attributes.pkl",
"wb"))
876 def load_particle_descriptors(self):
882 particles_attributes = pickle.load(open(
"particles_attributes.pkl",
884 rigid_body_attributes = pickle.load(open(
"rigid_body_attributes.pkl",
894 hparent=h.get_parent()
895 while hparent != hroot:
896 hparent=h.get_parent()
897 name+=
"|"+hparent.get_name()
901 xyzr.set_coordinates(particles_attributes[name][
"COORDINATES"])
907 for rb
in self.rigid_bodies:
909 trans=rigid_body_attributes[name][
"TRANSLATION"]
910 rot=rigid_body_attributes[name][
"ROTATION"]
913 rb.set_reference_frame(rf)
914 coor_nrm_ref=rigid_body_attributes[name][
"COORDINATES_NONRIGID_MEMBER"]
915 coor_rm_ref_dict=rigid_body_attributes[name][
"COORDINATES_RIGID_MEMBER"]
918 for mi
in rb.get_member_indexes():
919 rm=self.m.get_particle(mi)
921 name_part=rm.get_name()
922 xyz=coor_nrm_ref[name_part]
924 self.m.set_attribute(fk, rm,xyz[n])
926 name_part=rm.get_name()
928 coor_rm_model.append(
IMP.core.XYZ(rm).get_coordinates())
929 if len(coor_rm_model)==0:
continue
937 rmf_component_name=
None,
938 check_number_particles=
True,
939 representation_name_to_rmf_name_map=
None,
942 '''Read and replace coordinates from an RMF file.
943 Replace the coordinates of particles with the same name.
944 It assumes that the RMF and the representation have the particles
946 @param component_name Component name
947 @param rmf_component_name Name of the component in the RMF file
948 (if not specified, use `component_name`)
949 @param representation_name_to_rmf_name_map a dictionary that map
950 the original rmf particle name to the recipient particle component name
951 @param save_file: save a file with the names of particles of the component
956 prots = IMP.pmi.analysis.get_hiers_from_rmf(
962 raise ValueError(
"cannot read hierarchy from rmf")
975 p, self.hier_dict.keys())
976 if not rmf_component_name
is None and protname == rmf_component_name:
978 elif rmf_component_name
is None and protname == component_name:
984 reprnames=[p.get_name()
for p
in psrepr]
985 rmfnames=[p.get_name()
for p
in psrmf]
988 fl=open(component_name+
".txt",
"w")
989 for i
in itertools.izip_longest(reprnames,rmfnames): fl.write(str(i[0])+
","+str(i[1])+
"\n")
992 if check_number_particles
and not representation_name_to_rmf_name_map:
993 if len(psrmf) != len(psrepr):
994 fl=open(component_name+
".txt",
"w")
995 for i
in itertools.izip_longest(reprnames,rmfnames): fl.write(str(i[0])+
","+str(i[1])+
"\n")
996 raise ValueError(
"%s cannot proceed the rmf and the representation don't have the same number of particles; "
997 "particles in rmf: %s particles in the representation: %s" % (str(component_name), str(len(psrmf)), str(len(psrepr))))
1000 if not representation_name_to_rmf_name_map:
1001 for n, prmf
in enumerate(psrmf):
1003 prmfname = prmf.get_name()
1004 preprname = psrepr[n].get_name()
1006 raise ValueError(
"component %s cannot proceed if rigid bodies were initialized. Use the function before defining the rigid bodies" % component_name)
1008 raise ValueError(
"component %s cannot proceed if rigid bodies were initialized. Use the function before defining the rigid bodies" % component_name)
1010 if prmfname != preprname:
1011 print(
"set_coordinates_from_rmf: WARNING rmf particle and representation particles have not the same name %s %s " % (prmfname, preprname))
1021 g=gprmf.get_gaussian()
1022 grepr.set_gaussian(g)
1026 repr_name_particle_map={}
1027 rmf_name_particle_map={}
1029 rmf_name_particle_map[p.get_name()]=p
1033 for prepr
in psrepr:
1035 prmf=rmf_name_particle_map[representation_name_to_rmf_name_map[prepr.get_name()]]
1037 print(
"set_coordinates_from_rmf: WARNING missing particle name in representation_name_to_rmf_name_map, skipping...")
1045 If the root hierarchy does not exist, construct it.
1047 if "Res:" + str(int(resolution))
not in self.hier_representation[name]:
1049 root.set_name(name +
"_Res:" + str(int(resolution)))
1050 self.hier_representation[name][
1051 "Res:" + str(int(resolution))] = root
1052 protein_h.add_child(root)
1055 input_hierarchy, protein_h, type, color):
1057 Generate all needed coarse grained layers.
1059 @param name name of the protein
1060 @param resolutions list of resolutions
1061 @param protein_h root hierarchy
1062 @param input_hierarchy hierarchy to coarse grain
1063 @param type a string, typically "pdb" or "helix"
1067 if (1
in resolutions)
or (0
in resolutions):
1070 if 1
in resolutions:
1073 s1.set_name(
'%s_%i-%i_%s' % (name, start, end, type))
1075 self.hier_representation[name][
"Res:1"].add_child(s1)
1077 if 0
in resolutions:
1080 s0.set_name(
'%s_%i-%i_%s' % (name, start, end, type))
1082 self.hier_representation[name][
"Res:0"].add_child(s0)
1085 if not isnucleicacid:
1088 atom_type=IMP.atom.AT_CA)
1092 atom_type=IMP.atom.AT_P)
1094 for p
in sel.get_selected_particles():
1096 if 0
in resolutions:
1098 resclone0 = IMP.atom.create_clone(resobject)
1100 s0.add_child(resclone0)
1101 self.hier_db.add_particles(
1105 resclone0.get_children())
1107 chil = resclone0.get_children()
1116 if 1
in resolutions:
1118 resclone1 = IMP.atom.create_clone_one(resobject)
1120 s1.add_child(resclone1)
1121 self.hier_db.add_particles(name, resindex, 1, [resclone1])
1125 prt = resclone1.get_particle()
1126 prt.set_name(
'%s_%i_%s' % (name, resindex, type))
1148 for r
in resolutions:
1149 if r != 0
and r != 1:
1155 chil = s.get_children()
1157 s0.set_name(
'%s_%i-%i_%s' % (name, start, end, type))
1161 self.hier_representation[name][
1162 "Res:" + str(int(r))].add_child(s0)
1171 prt.set_name(
'%s_%i_%s' % (name, first, type))
1173 prt.set_name(
'%s_%i_%i_%s' % (name, first, last, type))
1175 self.hier_db.add_particles(name, kk, r, [prt])
1196 Get the hierarchies at the given resolution.
1198 The map between resolution and hierarchies is cached to accelerate
1199 the selection algorithm. The cache is invalidated when the
1200 representation was changed by any add_component_xxx.
1203 if self.representation_is_modified:
1204 rhbr = self.hier_db.get_all_root_hierarchies_by_resolution(
1206 self.hier_resolution[resolution] = rhbr
1207 self.representation_is_modified =
False
1210 if resolution
in self.hier_resolution:
1211 return self.hier_resolution[resolution]
1213 rhbr = self.hier_db.get_all_root_hierarchies_by_resolution(
1215 self.hier_resolution[resolution] = rhbr
1219 self, max_translation=300., max_rotation=2.0 * pi,
1220 avoidcollision=
True, cutoff=10.0, niterations=100,
1222 excluded_rigid_bodies=
None,
1223 hierarchies_excluded_from_collision=
None):
1225 Shuffle configuration; used to restart the optimization.
1227 The configuration of the system is initialized by placing each
1228 rigid body and each bead randomly in a box with a side of
1229 `max_translation` angstroms, and far enough from each other to
1230 prevent any steric clashes. The rigid bodies are also randomly rotated.
1232 @param avoidcollision check if the particle/rigid body was
1233 placed close to another particle; uses the optional
1234 arguments cutoff and niterations
1235 @param bounding_box defined by ((x1,y1,z1),(x2,y2,z2))
1238 if excluded_rigid_bodies
is None:
1239 excluded_rigid_bodies = []
1240 if hierarchies_excluded_from_collision
is None:
1241 hierarchies_excluded_from_collision = []
1243 if len(self.rigid_bodies) == 0:
1244 print(
"shuffle_configuration: rigid bodies were not intialized")
1247 gcpf.set_distance(cutoff)
1249 hierarchies_excluded_from_collision_indexes = []
1258 if not bounding_box
is None:
1259 ((x1, y1, z1), (x2, y2, z2)) = bounding_box
1264 for h
in hierarchies_excluded_from_collision:
1268 allparticleindexes = list(
1269 set(allparticleindexes) - set(hierarchies_excluded_from_collision_indexes))
1271 print(hierarchies_excluded_from_collision)
1272 print(len(allparticleindexes),len(hierarchies_excluded_from_collision_indexes))
1274 for rb
in self.rigid_bodies:
1275 if rb
not in excluded_rigid_bodies:
1278 rbindexes = rb.get_member_particle_indexes()
1281 set(rbindexes) - set(hierarchies_excluded_from_collision_indexes))
1282 otherparticleindexes = list(
1283 set(allparticleindexes) - set(rbindexes))
1285 if len(otherparticleindexes)
is None:
1289 while niter < niterations:
1290 rbxyz = (rb.get_x(), rb.get_y(), rb.get_z())
1292 if not bounding_box
is None:
1299 translation-old_coord)
1311 gcpf.get_close_pairs(
1313 otherparticleindexes,
1319 print(
"shuffle_configuration: rigid body placed close to other %d particles, trying again..." % npairs)
1320 print(
"shuffle_configuration: rigid body name: " + rb.get_name())
1321 if niter == niterations:
1322 raise ValueError(
"tried the maximum number of iterations to avoid collisions, increase the distance cutoff")
1326 print(
'shuffling', len(self.floppy_bodies),
'floppy bodies')
1327 for fb
in self.floppy_bodies:
1333 otherparticleindexes = list(
1334 set(allparticleindexes) - set(fbindexes))
1335 if len(otherparticleindexes)
is None:
1341 while niter < niterations:
1343 if not bounding_box
is None:
1365 gcpf.get_close_pairs(
1367 otherparticleindexes,
1373 print(
"shuffle_configuration: floppy body placed close to other %d particles, trying again..." % npairs)
1374 if niter == niterations:
1375 raise ValueError(
"tried the maximum number of iterations to avoid collisions, increase the distance cutoff")
1379 def set_current_coordinates_as_reference_for_rmsd(self, label="None"):
1383 self.reference_structures[label] = (
1387 def get_all_rmsds(self):
1390 for label
in self.reference_structures:
1392 current_coordinates = [
IMP.core.XYZ(p).get_coordinates()
1393 for p
in self.reference_structures[label][1]]
1394 reference_coordinates = self.reference_structures[label][0]
1395 if len(reference_coordinates) != len(current_coordinates):
1396 print(
"calculate_all_rmsds: reference and actual coordinates are not the same")
1399 current_coordinates,
1400 reference_coordinates)
1402 reference_coordinates,
1403 current_coordinates)
1407 reference_coordinates,
1408 current_coordinates)
1409 rmsds[label +
"_GlobalRMSD"] = rmsd_global
1410 rmsds[label +
"_RelativeDRMS"] = rmsd_relative
1413 def setup_component_geometry(self, name, color=None, resolution=1.0):
1415 color = self.color_dict[name]
1419 self.hier_geometry_pairs[name] = []
1420 protein_h = self.hier_dict[name]
1422 pbr = IMP.pmi.tools.sort_by_residues(pbr)
1424 for n
in range(len(pbr) - 1):
1425 self.hier_geometry_pairs[name].append((pbr[n], pbr[n + 1], color))
1428 self, name, resolution=10, scale=1.0):
1430 Generate restraints between contiguous fragments in the hierarchy.
1431 The linkers are generated at resolution 10 by default.
1438 protein_h = self.hier_dict[name]
1440 frs = self.hier_db.get_preroot_fragments_by_resolution(
1446 start = fr.get_children()[0]
1451 end = fr.get_children()[-1]
1457 SortedSegments.append((start, end, startres))
1458 SortedSegments = sorted(SortedSegments, key=itemgetter(2))
1461 for x
in range(len(SortedSegments) - 1):
1462 last = SortedSegments[x][1]
1463 first = SortedSegments[x + 1][0]
1470 residuegap = firstresn - lastresn - 1
1471 if self.disorderedlength
and (nreslast / 2 + nresfirst / 2 + residuegap) > 20.0:
1474 optdist = sqrt(5 / 3) * 1.93 * \
1475 (nreslast / 2 + nresfirst / 2 + residuegap) ** 0.6
1477 if self.upperharmonic:
1483 optdist = (0.0 + (float(residuegap) + 1.0) * 3.6) * scale
1484 if self.upperharmonic:
1490 pt0 = last.get_particle()
1491 pt1 = first.get_particle()
1494 print(
"Adding sequence connectivity restraint between", pt0.get_name(),
" and ", pt1.get_name(),
'of distance', optdist)
1495 sortedsegments_cr.add_restraint(r)
1496 self.linker_restraints_dict[
1497 "LinkerRestraint-" + pt0.get_name() +
"-" + pt1.get_name()] = r
1498 self.connected_intra_pairs.append((pt0, pt1))
1499 self.connected_intra_pairs.append((pt1, pt0))
1501 self.linker_restraints.add_restraint(sortedsegments_cr)
1502 self.linker_restraints.add_restraint(unmodeledregions_cr)
1504 self.sortedsegments_cr_dict[name] = sortedsegments_cr
1505 self.unmodeledregions_cr_dict[name] = unmodeledregions_cr
1507 def optimize_floppy_bodies(self, nsteps, temperature=1.0):
1509 pts = IMP.pmi.tools.ParticleToSampleList()
1510 for n, fb
in enumerate(self.floppy_bodies):
1511 pts.add_particle(fb,
"Floppy_Bodies", 1.0,
"Floppy_Body_" + str(n))
1512 if len(pts.get_particles_to_sample()) > 0:
1514 print(
"optimize_floppy_bodies: optimizing %i floppy bodies" % len(self.floppy_bodies))
1517 print(
"optimize_floppy_bodies: no particle to optimize")
1519 def create_rotational_symmetry(self, maincopy, copies):
1521 self.representation_is_modified =
True
1522 ncopies = len(copies) + 1
1525 mainparticles = sel.get_selected_particles()
1527 for k
in range(len(copies)):
1533 copyparticles = sel.get_selected_particles()
1537 for n, p
in enumerate(mainparticles):
1539 pc = copyparticles[n]
1541 mainpurged.append(p)
1544 copypurged.append(pc)
1548 for n, p
in enumerate(mainpurged):
1551 print(
"setting " + p.get_name() +
" as reference for " + pc.get_name())
1554 lc.add(pc.get_index())
1557 self.m.add_score_state(c)
1562 def create_rigid_body_symmetry(self, particles_reference, particles_copy,label="None",
1565 self.representation_is_modified =
True
1567 mainparticles = particles_reference
1569 t=initial_transformation
1571 p.set_name(
"RigidBody_Symmetry")
1577 copyparticles = particles_copy
1581 for n, p
in enumerate(mainparticles):
1583 pc = copyparticles[n]
1585 mainpurged.append(p)
1591 copypurged.append(pc)
1598 for n, p
in enumerate(mainpurged):
1601 print(
"setting " + p.get_name() +
" as reference for " + pc.get_name())
1604 lc.add(pc.get_index())
1607 self.m.add_score_state(c)
1610 self.rigid_bodies.append(rb)
1611 self.rigid_body_symmetries.append(rb)
1612 rb.set_name(label+
".rigid_body_symmetry."+str(len(self.rigid_body_symmetries)))
1615 def create_amyloid_fibril_symmetry(self, maincopy, axial_copies,
1616 longitudinal_copies, axis=(0, 0, 1), translation_value=4.8):
1619 self.representation_is_modified =
True
1622 protein_h = self.hier_dict[maincopy]
1625 for ilc
in range(-longitudinal_copies, longitudinal_copies + 1):
1626 for iac
in range(axial_copies):
1627 copyname = maincopy +
"_a" + str(ilc) +
"_l" + str(iac)
1628 self.create_component(copyname, 0.0)
1629 for hier
in protein_h.get_children():
1631 copyhier = self.hier_dict[copyname]
1632 outhiers.append(copyhier)
1636 2 * pi / axial_copies * (float(iac)))
1637 translation_vector = tuple(
1638 [translation_value * float(ilc) * x
for x
in axis])
1639 print(translation_vector)
1644 for n, p
in enumerate(mainparts):
1651 lc.add(pc.get_index())
1653 self.m.add_score_state(c)
1659 Load coordinates in the current representation.
1660 This should be done only if the hierarchy self.prot is identical
1661 to the one as stored in the rmf i.e. all components were added.
1666 rh = RMF.open_rmf_file_read_only(rmfname)
1674 create the representation (i.e. hierarchies) from the rmf file.
1675 it will be stored in self.prot, which will be overwritten.
1676 load the coordinates from the rmf file at frameindex.
1678 rh = RMF.open_rmf_file(rmfname)
1684 for p
in self.prot.get_children():
1685 self.create_component(p.get_name())
1686 self.hier_dict[p.get_name()] = p
1688 still missing: save rigid bodies contained in the rmf in self.rigid_bodies
1689 save floppy bodies in self.floppy_bodies
1690 get the connectivity restraints
1695 Construct a rigid body from hierarchies (and optionally particles).
1697 @param hiers list of hierarchies to make rigid
1698 @param particles list of particles to add to the rigid body
1701 if particles
is None:
1704 rigid_parts = set(particles)
1707 print(
"set_rigid_body_from_hierarchies> setting up a new rigid body")
1714 print(
"set_rigid_body_from_hierarchies> WARNING particle %s already belongs to rigid body %s" % (p.get_name(), rb.get_name()))
1717 name += hier.get_name() +
"-"
1718 print(
"set_rigid_body_from_hierarchies> adding %s to the rigid body" % hier.get_name())
1720 if len(list(rigid_parts)) != 0:
1722 rb.set_coordinates_are_optimized(
True)
1723 rb.set_name(name +
"rigid_body")
1724 self.rigid_bodies.append(rb)
1728 print(
"set_rigid_body_from_hierarchies> rigid body could not be setup")
1732 Construct a rigid body from a list of subunits.
1734 Each subunit is a tuple that identifies the residue ranges and the
1735 component name (as used in create_component()).
1737 subunits: [(name_1,(first_residue_1,last_residue_1)),(name_2,(first_residue_2,last_residue_2)),.....]
1739 [name_1,name_2,(name_3,(first_residue_3,last_residue_3)),.....]
1741 example: ["prot1","prot2",("prot3",(1,10))]
1743 sometimes, we know about structure of an interaction
1744 and here we make such PPIs rigid
1749 if type(s) == type(tuple())
and len(s) == 2:
1753 residue_indexes=list(range(s[1][0],
1755 if len(sel.get_selected_particles()) == 0:
1756 print(
"set_rigid_bodies: selected particle does not exist")
1757 for p
in sel.get_selected_particles():
1761 print(
"set_rigid_body_from_hierarchies> WARNING particle %s already belongs to rigid body %s" % (p.get_name(), rb.get_name()))
1765 elif type(s) == type(str()):
1767 if len(sel.get_selected_particles()) == 0:
1768 print(
"set_rigid_bodies: selected particle does not exist")
1769 for p
in sel.get_selected_particles():
1773 print(
"set_rigid_body_from_hierarchies> WARNING particle %s already belongs to rigid body %s" % (p.get_name(), rb.get_name()))
1778 rb.set_coordinates_are_optimized(
True)
1779 rb.set_name(
''.join(str(subunits)) +
"_rigid_body")
1780 self.rigid_bodies.append(rb)
1783 def set_super_rigid_body_from_hierarchies(
1789 super_rigid_xyzs = set()
1790 super_rigid_rbs = set()
1792 print(
"set_super_rigid_body_from_hierarchies> setting up a new SUPER rigid body")
1799 super_rigid_rbs.add(rb)
1801 super_rigid_xyzs.add(p)
1802 print(
"set_rigid_body_from_hierarchies> adding %s to the rigid body" % hier.get_name())
1803 if len(super_rigid_rbs) < min_size:
1806 self.super_rigid_bodies.append((super_rigid_xyzs, super_rigid_rbs))
1809 self.super_rigid_bodies.append(
1810 (super_rigid_xyzs, super_rigid_rbs, axis))
1812 def fix_rigid_bodies(self, rigid_bodies):
1813 self.fixed_rigid_bodies += rigid_bodies
1817 self, hiers, min_length=
None, max_length=
None, axis=
None):
1819 Make a chain of super rigid bodies from a list of hierarchies.
1821 Takes a linear list of hierarchies (they are supposed to be
1822 sequence-contiguous) and produces a chain of super rigid bodies
1823 with given length range, specified by min_length and max_length.
1826 hiers = IMP.pmi.tools.flatten_list(hiers)
1830 self.set_super_rigid_body_from_hierarchies(hs, axis, min_length)
1832 def set_super_rigid_bodies(self, subunits, coords=None):
1833 super_rigid_xyzs = set()
1834 super_rigid_rbs = set()
1837 if type(s) == type(tuple())
and len(s) == 3:
1841 residue_indexes=list(range(s[0],
1843 if len(sel.get_selected_particles()) == 0:
1844 print(
"set_rigid_bodies: selected particle does not exist")
1845 for p
in sel.get_selected_particles():
1848 super_rigid_rbs.add(rb)
1850 super_rigid_xyzs.add(p)
1851 elif type(s) == type(str()):
1853 if len(sel.get_selected_particles()) == 0:
1854 print(
"set_rigid_bodies: selected particle does not exist")
1855 for p
in sel.get_selected_particles():
1859 super_rigid_rbs.add(rb)
1861 super_rigid_xyzs.add(p)
1862 self.super_rigid_bodies.append((super_rigid_xyzs, super_rigid_rbs))
1866 Remove leaves of hierarchies from the floppy bodies list based
1867 on the component name
1870 ps=[h.get_particle()
for h
in hs]
1871 for p
in self.floppy_bodies:
1873 if p
in ps: self.floppy_bodies.remove(p)
1874 if p
in hs: self.floppy_bodies.remove(p)
1881 Remove leaves of hierarchies from the floppy bodies list.
1883 Given a list of hierarchies, get the leaves and remove the
1884 corresponding particles from the floppy bodies list. We need this
1885 function because sometimes
1886 we want to constrain the floppy bodies in a rigid body - for instance
1887 when you want to associate a bead with a density particle.
1889 for h
in hierarchies:
1892 if p
in self.floppy_bodies:
1893 print(
"remove_floppy_bodies: removing %s from floppy body list" % p.get_name())
1894 self.floppy_bodies.remove(p)
1897 def set_floppy_bodies(self):
1898 for p
in self.floppy_bodies:
1900 p.set_name(name +
"_floppy_body")
1902 print(
"I'm trying to make this particle flexible although it was assigned to a rigid body", p.get_name())
1905 rb.set_is_rigid_member(p.get_index(),
False)
1908 rb.set_is_rigid_member(p.get_particle_index(),
False)
1909 p.set_name(p.get_name() +
"_rigid_body_member")
1911 def set_floppy_bodies_from_hierarchies(self, hiers):
1916 self.floppy_bodies.append(p)
1923 selection tuples must be [(r1,r2,"name1"),(r1,r2,"name2"),....]
1924 @return the particles
1927 print(selection_tuples)
1928 for s
in selection_tuples:
1929 ps = IMP.pmi.tools.select_by_tuple(
1930 representation=self, tupleselection=tuple(s),
1931 resolution=
None, name_is_ambiguous=
False)
1935 def get_connected_intra_pairs(self):
1936 return self.connected_intra_pairs
1938 def set_rigid_bodies_max_trans(self, maxtrans):
1939 self.maxtrans_rb = maxtrans
1941 def set_rigid_bodies_max_rot(self, maxrot):
1942 self.maxrot_rb = maxrot
1944 def set_super_rigid_bodies_max_trans(self, maxtrans):
1945 self.maxtrans_srb = maxtrans
1947 def set_super_rigid_bodies_max_rot(self, maxrot):
1948 self.maxrot_srb = maxrot
1950 def set_floppy_bodies_max_trans(self, maxtrans):
1951 self.maxtrans_fb = maxtrans
1955 Fix rigid bodies in their actual position.
1956 The get_particles_to_sample() function will return
1957 just the floppy bodies (if they are not fixed).
1959 self.rigidbodiesarefixed = rigidbodiesarefixed
1963 Fix floppy bodies in their actual position.
1964 The get_particles_to_sample() function will return
1965 just the rigid bodies (if they are not fixed).
1967 self.floppybodiesarefixed=floppybodiesarefixed
1969 def draw_hierarchy_graph(self):
1971 print(
"Drawing hierarchy graph for " + c.get_name())
1972 IMP.pmi.output.get_graph_from_hierarchy(c)
1974 def get_geometries(self):
1977 for name
in self.hier_geometry_pairs:
1978 for pt
in self.hier_geometry_pairs[name]:
1992 def setup_bonds(self):
1995 for name
in self.hier_geometry_pairs:
1996 for pt
in self.hier_geometry_pairs[name]:
2011 def show_component_table(self, name):
2012 if name
in self.sequence_dict:
2013 lastresn = len(self.sequence_dict[name])
2016 residues = self.hier_db.get_residue_numbers(name)
2017 firstresn = min(residues)
2018 lastresn = max(residues)
2020 for nres
in range(firstresn, lastresn + 1):
2022 resolutions = self.hier_db.get_residue_resolutions(name, nres)
2024 for r
in resolutions:
2025 ps += self.hier_db.get_particles(name, nres, r)
2026 print(
"%20s %7s" % (name, nres),
" ".join([
"%20s %7s" % (str(p.get_name()),
2029 print(
"%20s %20s" % (name, nres),
"**** not represented ****")
2031 def draw_hierarchy_composition(self):
2033 ks = list(self.elements.keys())
2037 for k
in self.elements:
2038 for l
in self.elements[k]:
2043 self.draw_component_composition(k, max)
2045 def draw_component_composition(self, name, max=1000, draw_pdb_names=False):
2046 from matplotlib
import pyplot
2047 import matplotlib
as mpl
2049 tmplist = sorted(self.elements[k], key=itemgetter(0))
2051 endres = tmplist[-1][1]
2053 print(
"draw_component_composition: missing information for component %s" % name)
2055 fig = pyplot.figure(figsize=(26.0 * float(endres) / max + 2, 2))
2056 ax = fig.add_axes([0.05, 0.475, 0.9, 0.15])
2061 norm = mpl.colors.Normalize(vmin=5, vmax=10)
2065 for n, l
in enumerate(tmplist):
2069 if bounds[-1] != l[0]:
2070 colors.append(
"white")
2073 colors.append(
"#99CCFF")
2075 colors.append(
"#FFFF99")
2077 colors.append(
"#33CCCC")
2079 bounds.append(l[1] + 1)
2082 colors.append(
"#99CCFF")
2084 colors.append(
"#FFFF99")
2086 colors.append(
"#33CCCC")
2088 bounds.append(l[1] + 1)
2090 if bounds[-1] - 1 == l[0]:
2094 colors.append(
"white")
2097 bounds.append(bounds[-1])
2098 colors.append(
"white")
2099 cmap = mpl.colors.ListedColormap(colors)
2100 cmap.set_over(
'0.25')
2101 cmap.set_under(
'0.75')
2103 norm = mpl.colors.BoundaryNorm(bounds, cmap.N)
2104 cb2 = mpl.colorbar.ColorbarBase(ax, cmap=cmap,
2110 spacing=
'proportional',
2111 orientation=
'horizontal')
2120 mid = 1.0 / endres * float(l[0])
2125 l[2], xy=(mid, 1), xycoords=
'axes fraction',
2126 xytext=(mid + 0.025, float(npdb - 1) / 2.0 + 1.5), textcoords=
'axes fraction',
2127 arrowprops=dict(arrowstyle=
"->",
2128 connectionstyle=
"angle,angleA=0,angleB=90,rad=10"),
2130 extra_artists.append(t)
2133 title = ax.text(-0.005, 0.5, k, ha=
"right", va=
"center", rotation=90,
2136 extra_artists.append(title)
2138 labels = len(bounds) * [
" "]
2139 ax.set_xticklabels(labels)
2140 mid = 1.0 / endres * float(bounds[0])
2141 t = ax.annotate(bounds[0], xy=(mid, 0), xycoords=
'axes fraction',
2142 xytext=(mid - 0.01, -0.5), textcoords=
'axes fraction',)
2143 extra_artists.append(t)
2144 offsets = [0, -0.5, -1.0]
2146 for n
in range(1, len(bounds)):
2147 if bounds[n] == bounds[n - 1]:
2149 mid = 1.0 / endres * float(bounds[n])
2150 if (float(bounds[n]) - float(bounds[n - 1])) / max <= 0.01:
2152 offset = offsets[nclashes % 3]
2158 bounds[n], xy=(mid, 0), xycoords=
'axes fraction',
2159 xytext=(mid, -0.5 + offset), textcoords=
'axes fraction')
2162 bounds[n], xy=(mid, 0), xycoords=
'axes fraction',
2163 xytext=(mid, -0.5 + offset), textcoords=
'axes fraction', arrowprops=dict(arrowstyle=
"-"))
2164 extra_artists.append(t)
2166 cb2.add_lines(bounds, [
"black"] * len(bounds), [1] * len(bounds))
2170 k +
"structure.pdf",
2173 bbox_extra_artists=(extra_artists),
2174 bbox_inches=
'tight')
2177 def draw_coordinates_projection(self):
2178 import matplotlib.pyplot
as pp
2181 for name
in self.hier_geometry_pairs:
2182 for pt
in self.hier_geometry_pairs[name]:
2192 xpairs.append([x1, x2])
2193 ypairs.append([y1, y2])
2196 for xends, yends
in zip(xpairs, ypairs):
2201 pp.plot(xlist, ylist,
'b-', alpha=0.1)
2204 def get_prot_name_from_particle(self, particle):
2205 names = self.get_component_names()
2206 particle0 = particle
2208 while not name
in names:
2211 particle0 = h.get_particle()
2215 def get_particles_to_sample(self):
2225 if not self.rigidbodiesarefixed:
2226 for rb
in self.rigid_bodies:
2231 if rb
not in self.fixed_rigid_bodies:
2234 if not self.floppybodiesarefixed:
2235 for fb
in self.floppy_bodies:
2242 for srb
in self.super_rigid_bodies:
2245 rigid_bodies = list(srb[1])
2246 filtered_rigid_bodies = []
2247 for rb
in rigid_bodies:
2248 if rb
not in self.fixed_rigid_bodies:
2249 filtered_rigid_bodies.append(rb)
2250 srbtmp.append((srb[0], filtered_rigid_bodies))
2252 self.rigid_bodies = rbtmp
2253 self.floppy_bodies = fbtmp
2254 self.super_rigid_bodies = srbtmp
2256 ps[
"Rigid_Bodies_SimplifiedModel"] = (
2260 ps[
"Floppy_Bodies_SimplifiedModel"] = (
2263 ps[
"SR_Bodies_SimplifiedModel"] = (
2264 self.super_rigid_bodies,
2269 def set_output_level(self, level):
2270 self.output_level = level
2272 def _evaluate(self, deriv):
2273 """Evaluate the total score of all added restraints"""
2275 return r.evaluate(deriv)
2277 def get_output(self):
2281 output[
"SimplifiedModel_Total_Score_" +
2282 self.label] = str(self._evaluate(
False))
2283 output[
"SimplifiedModel_Linker_Score_" +
2284 self.label] = str(self.linker_restraints.unprotected_evaluate(
None))
2285 for name
in self.sortedsegments_cr_dict:
2286 partialscore = self.sortedsegments_cr_dict[name].evaluate(
False)
2287 score += partialscore
2289 "SimplifiedModel_Link_SortedSegments_" +
2294 partialscore = self.unmodeledregions_cr_dict[name].evaluate(
False)
2295 score += partialscore
2297 "SimplifiedModel_Link_UnmodeledRegions_" +
2302 for rb
in self.rigid_body_symmetries:
2304 output[name +
"_" +self.label]=str(rb.get_reference_frame().get_transformation_to())
2305 for name
in self.linker_restraints_dict:
2310 self.linker_restraints_dict[
2311 name].unprotected_evaluate(
2314 if len(self.reference_structures.keys()) != 0:
2315 rmsds = self.get_all_rmsds()
2318 "SimplifiedModel_" +
2321 self.label] = rmsds[
2324 if self.output_level ==
"high":
2327 output[
"Coordinates_" +
2328 p.get_name() +
"_" + self.label] = str(d)
2330 output[
"_TotalScore"] = str(score)
2333 def get_test_output(self):
2335 output = self.get_output()
2336 for n, p
in enumerate(self.get_particles_to_sample()):
2337 output[
"Particle_to_sample_" + str(n)] = str(p)
2339 output[
"Hierarchy_Dictionary"] = list(self.hier_dict.keys())
2340 output[
"Number_of_floppy_bodies"] = len(self.floppy_bodies)
2341 output[
"Number_of_rigid_bodies"] = len(self.rigid_bodies)
2342 output[
"Number_of_super_bodies"] = len(self.super_rigid_bodies)
2343 output[
"Selection_resolution_1"] = len(
2345 output[
"Selection_resolution_5"] = len(
2347 output[
"Selection_resolution_7"] = len(
2349 output[
"Selection_resolution_10"] = len(
2351 output[
"Selection_resolution_100"] = len(
2354 output[
"Selection_resolution=1"] = len(
2356 output[
"Selection_resolution=1,resid=10"] = len(
2358 for resolution
in self.hier_resolution:
2359 output[
"Hier_resolution_dictionary" +
2360 str(resolution)] = len(self.hier_resolution[resolution])
2361 for name
in self.hier_dict:
2363 "Selection_resolution=1,resid=10,name=" +
2371 "Selection_resolution=1,resid=10,name=" +
2373 ",ambiguous"] = len(
2378 name_is_ambiguous=
True,
2381 "Selection_resolution=1,resid=10,name=" +
2383 ",ambiguous"] = len(
2388 name_is_ambiguous=
True,
2391 "Selection_resolution=1,resrange=(10,20),name=" +
2400 "Selection_resolution=1,resrange=(10,20),name=" +
2402 ",ambiguous"] = len(
2407 name_is_ambiguous=
True,
2411 "Selection_resolution=10,resrange=(10,20),name=" +
2420 "Selection_resolution=10,resrange=(10,20),name=" +
2422 ",ambiguous"] = len(
2427 name_is_ambiguous=
True,
2431 "Selection_resolution=100,resrange=(10,20),name=" +
2440 "Selection_resolution=100,resrange=(10,20),name=" +
2442 ",ambiguous"] = len(
2447 name_is_ambiguous=
True,
2459 def __init__(self, *args, **kwargs):
2460 Representation.__init__(self, *args, **kwargs)
def link_components_to_rmf
Load coordinates in the current representation.
Select non water and non hydrogen atoms.
double get_volume_from_residue_type(ResidueType rt)
Return an estimate for the volume of a given residue.
A decorator to associate a particle with a part of a protein/DNA/RNA.
static Gaussian setup_particle(Model *m, ParticleIndex pi)
void show_molecular_hierarchy(Hierarchy h)
Print out the molecular hierarchy.
double get_drms(const Vector3DsOrXYZs0 &m1, const Vector3DsOrXYZs1 &m2)
Apply a SingletonFunction to a SingletonContainer to maintain an invariant.
atom::Hierarchies create_hierarchies(RMF::FileConstHandle fh, Model *m)
def shuffle_configuration
Shuffle configuration; used to restart the optimization.
A member of a rigid body, it has internal (local) coordinates.
static bool get_is_setup(const IMP::ParticleAdaptor &p)
Set of python classes to create a multi-state, multi-resolution IMP hierarchy.
static Atom setup_particle(Model *m, ParticleIndex pi, Atom other)
static Fragment setup_particle(Model *m, ParticleIndex pi)
Upper bound harmonic function (non-zero when feature > mean)
static Symmetric setup_particle(Model *m, ParticleIndex pi, Float symmetric)
def set_coordinates_from_rmf
Read and replace coordinates from an RMF file.
static XYZR setup_particle(Model *m, ParticleIndex pi)
A decorator for a particle which has bonds.
Select atoms which are selected by both selectors.
Rotation3D get_random_rotation_3d(const Rotation3D ¢er, double distance)
Pick a rotation at random near the provided one.
double get_mass(ResidueType c)
Get the mass from the residue type.
Color get_rgb_color(double f)
Return the color for f from the RGB color map.
double get_mass_from_number_of_residues(unsigned int num_aa)
Estimate the mass of a protein from the number of amino acids.
def set_rigid_bodies_as_fixed
Fix rigid bodies in their actual position.
def get_particles_from_selection_tuples
selection tuples must be [(r1,r2,"name1"),(r1,r2,"name2"),....
Add symmetric attribute to a particle.
double get_ball_radius_from_volume_3d(double volume)
Return the radius of a sphere with a given volume.
def setup_component_sequence_connectivity
Generate restraints between contiguous fragments in the hierarchy.
Add uncertainty to a particle.
A score on the distance between the surfaces of two spheres.
static Residue setup_particle(Model *m, ParticleIndex pi, ResidueType t, int index, int insertion_code)
static bool get_is_setup(const IMP::ParticleAdaptor &p)
static bool get_is_setup(const IMP::ParticleAdaptor &p)
static XYZ setup_particle(Model *m, ParticleIndex pi)
void read_pdb(TextInput input, int model, Hierarchy h)
def remove_floppy_bodies
Remove leaves of hierarchies from the floppy bodies list.
Vector3D get_random_vector_in(const Cylinder3D &c)
Generate a random vector in a cylinder with uniform density.
static Uncertainty setup_particle(Model *m, ParticleIndex pi, Float uncertainty)
def remove_floppy_bodies_from_component
Remove leaves of hierarchies from the floppy bodies list based on the component name.
Add resolution to a particle.
Bond create_bond(Bonded a, Bonded b, Bond o)
Connect the two wrapped particles by a custom bond.
Object used to hold a set of restraints.
static Molecule setup_particle(Model *m, ParticleIndex pi)
Rotation3D get_rotation_about_axis(const Vector3D &axis, double angle)
Generate a Rotation3D object from a rotation around an axis.
Select all CB ATOM records.
A Gaussian distribution in 3D.
static bool get_is_setup(Model *m, ParticleIndex pi)
static Hierarchy setup_particle(Model *m, ParticleIndex pi, ParticleIndexesAdaptor children=ParticleIndexesAdaptor())
Create a Hierarchy of level t by adding the needed attributes.
ParticleIndexPairs get_indexes(const ParticlePairsTemp &ps)
def set_rigid_bodies
Construct a rigid body from a list of subunits.
double get_volume_from_mass(double m, ProteinDensityReference ref=ALBER)
Estimate the volume of a protein from its mass.
The standard decorator for manipulating molecular structures.
Ints get_index(const ParticlesTemp &particles, const Subset &subset, const Subsets &excluded)
Store a list of ParticleIndexes.
def deprecated_method
Python decorator to mark a method as deprecated.
A decorator for a particle representing an atom.
void transform(XYZ a, const algebra::Transformation3D &tr)
Apply a transformation to the particle.
static Mass setup_particle(Model *m, ParticleIndex pi, Float mass)
static Bonded setup_particle(Model *m, ParticleIndex pi)
static bool get_is_setup(Model *m, ParticleIndex pi)
double get_rmsd(const Vector3DsOrXYZs0 &m1, const Vector3DsOrXYZs1 &m2)
def add_all_atom_densities
Decorates all specified particles as Gaussians directly.
def coarse_hierarchy
Generate all needed coarse grained layers.
def add_component_pdb
Add a component that has an associated 3D structure in a PDB file.
Basic utilities for handling cryo-electron microscopy 3D density maps.
void load_frame(RMF::FileConstHandle file, RMF::FrameID frame)
A decorator for a particle with x,y,z coordinates.
static bool get_is_setup(Model *m, ParticleIndex pi)
static Colored setup_particle(Model *m, ParticleIndex pi, Color color)
def set_chain_of_super_rigid_bodies
Make a chain of super rigid bodies from a list of hierarchies.
def add_component_sequence
Add the primary sequence for a single component.
Tools for clustering and cluster analysis.
static Resolution setup_particle(Model *m, ParticleIndex pi, Float resolution)
def create_components_from_rmf
still not working.
Transformation3D get_identity_transformation_3d()
Return a transformation that does not do anything.
Find all nearby pairs by testing all pairs.
Classes for writing output files and processing them.
Simple implementation of segments in 3D.
static bool get_is_setup(const IMP::ParticleAdaptor &p)
def deprecated_object
Python decorator to mark a class as deprecated.
A decorator for a residue.
Basic functionality that is expected to be used by a wide variety of IMP users.
General purpose algebraic and geometric methods that are expected to be used by a wide variety of IMP...
static bool get_is_setup(const IMP::ParticleAdaptor &p)
def add_component_density
Sets up a Gaussian Mixture Model for this component.
Sample using Monte Carlo.
The general base class for IMP exceptions.
Hierarchy create_simplified_along_backbone(Chain input, const IntRanges &residue_segments, bool keep_detailed=false)
IMP::core::RigidBody create_rigid_body(Hierarchy h)
static Reference setup_particle(Model *m, ParticleIndex pi, ParticleIndexAdaptor reference)
Rotation3D get_identity_rotation_3d()
Return a rotation that does not do anything.
void link_hierarchies(RMF::FileConstHandle fh, const atom::Hierarchies &hs)
Class to handle individual model particles.
Bond get_bond(Bonded a, Bonded b)
Get the bond between two particles.
std::string get_data_path(std::string file_name)
Return the full path to one of this module's data files.
Select atoms which are selected by either or both selectors.
def get_hierarchies_at_given_resolution
Get the hierarchies at the given resolution.
Store info for a chain of a protein.
Applies a PairScore to a Pair.
Select all CA ATOM records.
Python classes to represent, score, sample and analyze models.
static bool get_is_setup(Model *m, ParticleIndex pi)
double get_resolution(Model *m, ParticleIndex pi)
Estimate the resolution of the hierarchy as used by Representation.
A decorator for a rigid body.
Set up the representation of all proteins and nucleic acid macromolecules.
static bool get_is_setup(const IMP::ParticleAdaptor &p)
A dictionary-like wrapper for reading and storing sequence data.
Transformation3D get_transformation_aligning_first_to_second(Vector3Ds a, Vector3Ds b)
Output IMP model data in various file formats.
Functionality for loading, creating, manipulating and scoring atomic structures.
static Chain setup_particle(Model *m, ParticleIndex pi, std::string id)
Hierarchies get_leaves(const Selection &h)
def add_component_necklace
Generates a string of beads with given length.
An exception for an invalid value being passed to IMP.
def set_floppy_bodies_as_fixed
Fix floppy bodies in their actual position.
Select hierarchy particles identified by the biological name.
static RigidBody setup_particle(Model *m, ParticleIndex pi, ParticleIndexesAdaptor ps)
Select all ATOM and HETATM records with the given chain ids.
Support for the RMF file format for storing hierarchical molecular data and markup.
def set_rigid_body_from_hierarchies
Construct a rigid body from hierarchies (and optionally particles).
static bool get_is_setup(const IMP::ParticleAdaptor &p)
Transformation3D get_random_local_transformation(Vector3D origin, double max_translation=5., double max_angle_in_rad=0.26)
Get a local transformation.
def check_root
If the root hierarchy does not exist, construct it.
def add_component_beads
add beads to the representation
Inferential scoring building on methods developed as part of the Inferential Structure Determination ...
def add_component_hierarchy_clone
Make a copy of a hierarchy and append it to a component.
Harmonic function (symmetric about the mean)
A decorator for a particle with x,y,z coordinates and a radius.