IMP Reference Guide
2.6.1
The Integrative Modeling Platform
|
Set up the representation of all proteins and nucleic acid macromolecules. More...
Inherits object.
Set up the representation of all proteins and nucleic acid macromolecules.
Create the molecular hierarchies, representation, sequence connectivity for the various involved proteins and nucleic acid macromolecules:
Create a protein, DNA or RNA, represent it as a set of connected balls of appropriate radii and number of residues, PDB at given resolution(s), or ideal helices.
How to use the SimplifiedModel class (typical use):
see test/test_hierarchy_contruction.py
examples:
1) Create a chain of helices and flexible parts
c_1_119 =self.add_component_necklace("prot1",1,119,20) c_120_131 =self.add_component_ideal_helix("prot1",resolutions=[1,10],resrange=(120,131)) c_132_138 =self.add_component_beads("prot1",[(132,138)]) c_139_156 =self.add_component_ideal_helix("prot1",resolutions=[1,10],resrange=(139,156)) c_157_174 =self.add_component_beads("prot1",[(157,174)]) c_175_182 =self.add_component_ideal_helix("prot1",resolutions=[1,10],resrange=(175,182)) c_183_194 =self.add_component_beads("prot1",[(183,194)]) c_195_216 =self.add_component_ideal_helix("prot1",resolutions=[1,10],resrange=(195,216)) c_217_250 =self.add_component_beads("prot1",[(217,250)])
self.set_rigid_body_from_hierarchies(c_120_131) self.set_rigid_body_from_hierarchies(c_139_156) self.set_rigid_body_from_hierarchies(c_175_182) self.set_rigid_body_from_hierarchies(c_195_216)
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, c_217_250]
self.set_chain_of_super_rigid_bodies(clist,2,3)
self.set_super_rigid_bodies(["prot1"])
Definition at line 45 of file pmi/representation.py.
Public Member Functions | |
def | __init__ |
Constructor. More... | |
def | add_all_atom_densities |
Decorates all specified particles as Gaussians directly. More... | |
def | add_component_beads |
add beads to the representation More... | |
def | add_component_density |
Sets up a Gaussian Mixture Model for this component. More... | |
def | add_component_hierarchy_clone |
Make a copy of a hierarchy and append it to a component. More... | |
def | add_component_necklace |
Generates a string of beads with given length. More... | |
def | add_component_pdb |
Add a component that has an associated 3D structure in a PDB file. More... | |
def | add_component_sequence |
Add the primary sequence for a single component. More... | |
def | check_root |
If the root hierarchy does not exist, construct it. More... | |
def | coarse_hierarchy |
Generate all needed coarse grained layers. More... | |
def | create_components_from_rmf |
still not working. More... | |
def | get_hierarchies_at_given_resolution |
Get the hierarchies at the given resolution. More... | |
def | get_particles_from_selection_tuples |
selection tuples must be [(r1,r2,"name1"),(r1,r2,"name2"),.... More... | |
def | link_components_to_rmf |
Load coordinates in the current representation. More... | |
def | remove_floppy_bodies |
Remove leaves of hierarchies from the floppy bodies list. More... | |
def | remove_floppy_bodies_from_component |
Remove leaves of hierarchies from the floppy bodies list based on the component name. More... | |
def | set_chain_of_super_rigid_bodies |
Make a chain of super rigid bodies from a list of hierarchies. More... | |
def | set_coordinates_from_rmf |
Read and replace coordinates from an RMF file. More... | |
def | set_floppy_bodies_as_fixed |
Fix floppy bodies in their actual position. More... | |
def | set_rigid_bodies |
Construct a rigid body from a list of subunits. More... | |
def | set_rigid_bodies_as_fixed |
Fix rigid bodies in their actual position. More... | |
def | set_rigid_body_from_hierarchies |
Construct a rigid body from hierarchies (and optionally particles). More... | |
def | setup_component_sequence_connectivity |
Generate restraints between contiguous fragments in the hierarchy. More... | |
def | shuffle_configuration |
Shuffle configuration; used to restart the optimization. More... | |
def IMP.pmi.representation.Representation.__init__ | ( | self, | |
m, | |||
upperharmonic = True , |
|||
disorderedlength = True |
|||
) |
Constructor.
m | the model |
upperharmonic | This flag uses either harmonic (False) or upperharmonic (True) in the intra-pair connectivity restraint. |
disorderedlength | This flag uses either disordered length calculated for random coil peptides (True) or zero surface-to-surface distance between beads (False) as optimal distance for the sequence connectivity restraint. |
Definition at line 70 of file pmi/representation.py.
def IMP.pmi.representation.Representation.add_all_atom_densities | ( | self, | |
name, | |||
hierarchies = None , |
|||
selection_tuples = None , |
|||
particles = None , |
|||
resolution = 0 , |
|||
output_map = None , |
|||
voxel_size = 1.0 |
|||
) |
Decorates all specified particles as Gaussians directly.
name | component name |
hierarchies | set up GMM for some hierarchies |
selection_tuples | (list of tuples) example (first_residue,last_residue,component_name) |
particles | set up GMM for particles directly |
resolution | usual PMI resolution for selecting particles from the hierarchies |
Definition at line 723 of file pmi/representation.py.
def IMP.pmi.representation.Representation.add_component_beads | ( | self, | |
name, | |||
ds, | |||
colors = None , |
|||
incoord = None |
|||
) |
add beads to the representation
name | the component name |
ds | a list of tuples corresponding to the residue ranges of the beads |
colors | a list of colors associated to the beads |
incoord | the coordinate tuple corresponding to the position of the beads |
Definition at line 478 of file pmi/representation.py.
def IMP.pmi.representation.Representation.add_component_density | ( | self, | |
name, | |||
hierarchies = None , |
|||
selection_tuples = None , |
|||
particles = None , |
|||
resolution = 0.0 , |
|||
num_components = 10 , |
|||
inputfile = None , |
|||
outputfile = None , |
|||
outputmap = None , |
|||
kernel_type = None , |
|||
covariance_type = 'full' , |
|||
voxel_size = 1.0 , |
|||
out_hier_name = '' , |
|||
sampled_points = 1000000 , |
|||
num_iter = 100 , |
|||
simulation_res = 1.0 , |
|||
multiply_by_total_mass = True , |
|||
transform = None , |
|||
intermediate_map_fn = None , |
|||
density_ps_to_copy = None , |
|||
use_precomputed_gaussians = False |
|||
) |
Sets up a Gaussian Mixture Model for this component.
Can specify input GMM file or it will be computed.
name | component name |
hierarchies | set up GMM for some hierarchies |
selection_tuples | (list of tuples) example (first_residue,last_residue,component_name) |
particles | set up GMM for particles directly |
resolution | usual PMI resolution for selecting particles from the hierarchies |
inputfile | read the GMM from this file |
outputfile | fit and write the GMM to this file (text) |
outputmap | after fitting, create GMM density file (mrc) |
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 |
covariance_type | for fitting the GMM. options are 'full', 'diagonal' and 'spherical' |
voxel_size | for creating the intermediate density map and output map. lower number increases accuracy but also rasterizing time grows |
out_hier_name | name of the output density hierarchy |
sampled_points | number of points to sample. more will increase accuracy and fitting time |
num_iter | num GMM iterations. more will increase accuracy and fitting time |
multiply_by_total_mass | multiply the weights of the GMM by this value (only works on creation!) |
transform | for input file only, apply a transformation (eg for multiple copies same GMM) |
intermediate_map_fn | for debugging, this will write the intermediate (simulated) map |
density_ps_to_copy | in case you already created the appropriate GMM (eg, for beads) |
use_precomputed_gaussians | Set this flag and pass fragments - will use roughly spherical Gaussian setup |
Definition at line 596 of file pmi/representation.py.
def IMP.pmi.representation.Representation.add_component_hierarchy_clone | ( | self, | |
name, | |||
hierarchy | |||
) |
Make a copy of a hierarchy and append it to a component.
Definition at line 783 of file pmi/representation.py.
def IMP.pmi.representation.Representation.add_component_necklace | ( | self, | |
name, | |||
begin, | |||
end, | |||
length, | |||
color = None , |
|||
incoord = None |
|||
) |
Generates a string of beads with given length.
Definition at line 572 of file pmi/representation.py.
def IMP.pmi.representation.Representation.add_component_pdb | ( | self, | |
name, | |||
pdbname, | |||
chain, | |||
resolutions, | |||
color = None , |
|||
resrange = None , |
|||
offset = 0 , |
|||
cacenters = True , |
|||
show = False , |
|||
isnucleicacid = False , |
|||
readnonwateratoms = False , |
|||
read_ca_cb_only = False |
|||
) |
Add a component that has an associated 3D structure in a PDB file.
Reads the PDB, and constructs the fragments corresponding to contiguous sequence stretches.
name | (string) the name of the component |
pdbname | (string) the name of the PDB file |
chain | (string or integer) can be either a string (eg, "A") or an integer (eg, 0 or 1) in case you want to get the corresponding chain number in the PDB. |
resolutions | (integers) a list of integers that corresponds to the resolutions that have to be generated |
color | (float from 0 to 1) the color applied to the hierarchies generated |
resrange | (tuple of integers): the residue range to extract from the PDB. It is a tuple (beg,end). If not specified, all residues belonging to the specified chain are read. |
offset | (integer) specifies the residue index offset to be applied when reading the PDB (the FASTA sequence is assumed to start from residue 1, so use this parameter if the PDB file does not start at residue 1) |
cacenters | (boolean) if True generates resolution=1 beads centered on C-alpha atoms. |
show | (boolean) print out the molecular hierarchy at the end. |
isnucleicacid | (boolean) use True if you're reading a PDB with nucleic acids. |
readnonwateratoms | (boolean) if True fixes some pathological PDB. |
read_ca_cb_only | (boolean) if True, only reads CA/CB |
Definition at line 288 of file pmi/representation.py.
def IMP.pmi.representation.Representation.add_component_sequence | ( | self, | |
name, | |||
filename, | |||
id = None , |
|||
offs = None , |
|||
format = 'FASTA' |
|||
) |
Add the primary sequence for a single component.
name | Human-readable name of the component |
filename | Name of the FASTA file |
id | Identifier of the sequence in the FASTA file header (if not provided, use name instead) |
Definition at line 172 of file pmi/representation.py.
def IMP.pmi.representation.Representation.check_root | ( | self, | |
name, | |||
protein_h, | |||
resolution | |||
) |
If the root hierarchy does not exist, construct it.
Definition at line 1043 of file pmi/representation.py.
def IMP.pmi.representation.Representation.coarse_hierarchy | ( | self, | |
name, | |||
start, | |||
end, | |||
resolutions, | |||
isnucleicacid, | |||
input_hierarchy, | |||
protein_h, | |||
type, | |||
color | |||
) |
Generate all needed coarse grained layers.
name | name of the protein |
resolutions | list of resolutions |
protein_h | root hierarchy |
input_hierarchy | hierarchy to coarse grain |
type | a string, typically "pdb" or "helix" |
Definition at line 1054 of file pmi/representation.py.
def IMP.pmi.representation.Representation.create_components_from_rmf | ( | self, | |
rmfname, | |||
frameindex | |||
) |
still not working.
create the representation (i.e. hierarchies) from the rmf file. it will be stored in self.prot, which will be overwritten. load the coordinates from the rmf file at frameindex.
Definition at line 1671 of file pmi/representation.py.
def IMP.pmi.representation.Representation.get_hierarchies_at_given_resolution | ( | self, | |
resolution | |||
) |
Get the hierarchies at the given resolution.
The map between resolution and hierarchies is cached to accelerate the selection algorithm. The cache is invalidated when the representation was changed by any add_component_xxx.
Definition at line 1194 of file pmi/representation.py.
def IMP.pmi.representation.Representation.get_particles_from_selection_tuples | ( | self, | |
selection_tuples, | |||
resolution = None |
|||
) |
selection tuples must be [(r1,r2,"name1"),(r1,r2,"name2"),....
]
Definition at line 1918 of file pmi/representation.py.
def IMP.pmi.representation.Representation.link_components_to_rmf | ( | self, | |
rmfname, | |||
frameindex | |||
) |
Load coordinates in the current representation.
This should be done only if the hierarchy self.prot is identical to the one as stored in the rmf i.e. all components were added.
Definition at line 1657 of file pmi/representation.py.
def IMP.pmi.representation.Representation.remove_floppy_bodies | ( | self, | |
hierarchies | |||
) |
Remove leaves of hierarchies from the floppy bodies list.
Given a list of hierarchies, get the leaves and remove the corresponding particles from the floppy bodies list. We need this function because sometimes we want to constrain the floppy bodies in a rigid body - for instance when you want to associate a bead with a density particle.
Definition at line 1879 of file pmi/representation.py.
def IMP.pmi.representation.Representation.remove_floppy_bodies_from_component | ( | self, | |
component_name | |||
) |
Remove leaves of hierarchies from the floppy bodies list based on the component name.
Definition at line 1864 of file pmi/representation.py.
def IMP.pmi.representation.Representation.set_chain_of_super_rigid_bodies | ( | self, | |
hiers, | |||
min_length = None , |
|||
max_length = None , |
|||
axis = None |
|||
) |
Make a chain of super rigid bodies from a list of hierarchies.
Takes a linear list of hierarchies (they are supposed to be sequence-contiguous) and produces a chain of super rigid bodies with given length range, specified by min_length and max_length.
Definition at line 1816 of file pmi/representation.py.
def IMP.pmi.representation.Representation.set_coordinates_from_rmf | ( | self, | |
component_name, | |||
rmf_file_name, | |||
rmf_frame_number, | |||
rmf_component_name = None , |
|||
check_number_particles = True , |
|||
representation_name_to_rmf_name_map = None , |
|||
state_number = 0 , |
|||
save_file = False |
|||
) |
Read and replace coordinates from an RMF file.
Replace the coordinates of particles with the same name. It assumes that the RMF and the representation have the particles in the same order.
component_name | Component name |
rmf_component_name | Name of the component in the RMF file (if not specified, use component_name ) |
representation_name_to_rmf_name_map | a dictionary that map the original rmf particle name to the recipient particle component name |
save_file | save a file with the names of particles of the component |
Definition at line 935 of file pmi/representation.py.
def IMP.pmi.representation.Representation.set_floppy_bodies_as_fixed | ( | self, | |
floppybodiesarefixed = True |
|||
) |
Fix floppy bodies in their actual position.
The get_particles_to_sample() function will return just the rigid bodies (if they are not fixed).
Definition at line 1961 of file pmi/representation.py.
def IMP.pmi.representation.Representation.set_rigid_bodies | ( | self, | |
subunits | |||
) |
Construct a rigid body from a list of subunits.
Each subunit is a tuple that identifies the residue ranges and the component name (as used in create_component()).
subunits: [(name_1,(first_residue_1,last_residue_1)),(name_2,(first_residue_2,last_residue_2)),.....] or [name_1,name_2,(name_3,(first_residue_3,last_residue_3)),.....]
example: ["prot1","prot2",("prot3",(1,10))]
sometimes, we know about structure of an interaction and here we make such PPIs rigid
Definition at line 1730 of file pmi/representation.py.
def IMP.pmi.representation.Representation.set_rigid_bodies_as_fixed | ( | self, | |
rigidbodiesarefixed = True |
|||
) |
Fix rigid bodies in their actual position.
The get_particles_to_sample() function will return just the floppy bodies (if they are not fixed).
Definition at line 1953 of file pmi/representation.py.
def IMP.pmi.representation.Representation.set_rigid_body_from_hierarchies | ( | self, | |
hiers, | |||
particles = None |
|||
) |
Construct a rigid body from hierarchies (and optionally particles).
hiers | list of hierarchies to make rigid |
particles | list of particles to add to the rigid body |
Definition at line 1693 of file pmi/representation.py.
def IMP.pmi.representation.Representation.setup_component_sequence_connectivity | ( | self, | |
name, | |||
resolution = 10 , |
|||
scale = 1.0 |
|||
) |
Generate restraints between contiguous fragments in the hierarchy.
The linkers are generated at resolution 10 by default.
Definition at line 1427 of file pmi/representation.py.
def IMP.pmi.representation.Representation.shuffle_configuration | ( | self, | |
max_translation = 300.0 , |
|||
max_rotation = 2.0 * pi , |
|||
avoidcollision = True , |
|||
cutoff = 10.0 , |
|||
niterations = 100 , |
|||
bounding_box = None , |
|||
excluded_rigid_bodies = None , |
|||
hierarchies_excluded_from_collision = None |
|||
) |
Shuffle configuration; used to restart the optimization.
The configuration of the system is initialized by placing each rigid body and each bead randomly in a box with a side of max_translation
angstroms, and far enough from each other to prevent any steric clashes. The rigid bodies are also randomly rotated.
avoidcollision | check if the particle/rigid body was placed close to another particle; uses the optional arguments cutoff and niterations |
bounding_box | defined by ((x1,y1,z1),(x2,y2,z2)) |
Definition at line 1218 of file pmi/representation.py.