1 """@namespace IMP.pmi.topology
2 Set of Python classes to create a multi-state, multi-resolution IMP hierarchy.
3 * Start by creating a System with `model = IMP.Model(); s = IMP.pmi.topology.System(model)`. The System will store all the states.
4 * Then call System.create_state(). You can easily create a multistate system by calling this function multiple times.
5 * For each State, call State.create_molecule() to add a Molecule (a uniquely named polymer). This function returns the Molecule object which can be passed to various PMI functions.
6 * Some useful functions to help you set up your Molecules:
7 * Access the sequence residues with slicing (Molecule[a:b]) or functions like Molecule.get_atomic_residues() and Molecule.get_non_atomic_residues(). These functions all return Python sets for easy set arithmetic using & (and), | (or), - (difference)
8 * Molecule.add_structure() to add structural information from a PDB file.
9 * Molecule.add_representation() to create a representation unit - here you can choose bead resolutions as well as alternate representations like densities or ideal helices.
10 * Molecule.create_clone() lets you set up a molecule with identical representations, just a different chain ID. Use Molecule.create_copy() if you want a molecule with the same sequence but that allows custom representations.
11 * Once data has been added and representations chosen, call System.build() to create a canonical IMP hierarchy.
12 * Following hierarchy construction, setup rigid bodies, flexible beads, etc in IMP::pmi::dof.
13 * Check your representation with a nice printout: IMP::atom::show_with_representation()
14 See a [comprehensive example](https://integrativemodeling.org/nightly/doc/ref/pmi_2multiscale_8py-example.html) for using these classes.
16 Alternatively one can construct the entire topology and degrees of freedom via formatted text file with TopologyReader and IMP::pmi::macros::BuildSystem(). This is used in the [PMI tutorial](@ref rnapolii_stalk).
17 Note that this only allows a limited set of the full options available to PMI users (rigid bodies only, fixed resolutions).
20 from __future__
import print_function
28 from collections
import defaultdict, namedtuple
29 from .
import system_tools
30 from bisect
import bisect_left
31 from math
import pi,cos,sin
32 from operator
import itemgetter
35 def _build_ideal_helix(model, residues, coord_finder):
36 """Creates an ideal helix from the specified residue range
37 Residues MUST be contiguous.
38 This function actually adds them to the TempResidue hierarchy
43 for n, tempres
in enumerate(residues):
44 if tempres.get_has_structure():
45 raise ValueError(
"You tried to build ideal_helix for a residue "
46 "that already has structure: %s" % tempres)
47 if n > 0
and tempres.get_index() != prev_idx + 1:
48 raise ValueError(
"Passed non-contiguous segment to "
49 "build_ideal_helix for %s" % tempres.get_molecule())
53 rp.set_name(
"Residue_%i" % tempres.get_index())
61 x = 2.3 * cos(n * 2 * pi / 3.6)
62 y = 2.3 * sin(n * 2 * pi / 3.6)
63 z = 6.2 / 3.6 / 2 * n * 2 * pi / 3.6
71 tempres.set_structure(this_res)
72 created_hiers.append(this_res)
73 prev_idx = tempres.get_index()
74 coord_finder.add_residues(created_hiers)
76 class _SystemBase(object):
77 """The base class for System, State and Molecule
78 classes. It contains shared functions in common to these classes
81 def __init__(self,model=None):
92 def _create_hierarchy(self):
93 """create a new hierarchy"""
97 def _create_child(self,parent_hierarchy):
98 """create a new hierarchy, set it as child of the input
100 child_hierarchy=self._create_hierarchy()
101 parent_hierarchy.add_child(child_hierarchy)
102 return child_hierarchy
105 """Build the coordinates of the system.
106 Loop through stored(?) hierarchies and set up coordinates!"""
111 class System(_SystemBase):
112 """This class initializes the root node of the global IMP.atom.Hierarchy."""
113 def __init__(self,model=None,name="System"):
114 _SystemBase.__init__(self,model)
115 self._number_of_states = 0
116 self._protocol_output = []
121 self.hier=self._create_hierarchy()
122 self.hier.set_name(name)
123 self.hier._pmi2_system = weakref.ref(self)
125 def get_states(self):
129 """Makes and returns a new IMP.pmi.topology.State in this system"""
130 self._number_of_states+=1
131 state =
State(self,self._number_of_states-1)
132 self.states.append(state)
136 return self.hier.get_name()
139 """Returns the total number of states generated"""
140 return self._number_of_states
142 def get_hierarchy(self):
146 """Build all states"""
148 for state
in self.states:
149 state.build(**kwargs)
154 """Capture details of the modeling protocol.
155 @param p an instance of IMP.pmi.output.ProtocolOutput or a subclass.
157 self._protocol_output.append(p)
160 for state
in self.states:
161 state._add_protocol_output(p, self)
167 """Stores a list of Molecules all with the same State index.
168 Also stores number of copies of each Molecule for easy Selection.
170 def __init__(self,system,state_index):
171 """Define a new state
172 @param system the PMI System
173 @param state_index the index of the new state
174 @note It's expected that you will not use this constructor directly,
175 but rather create it with pmi::System::create_state()
177 self.model = system.get_hierarchy().get_model()
179 self.hier = self._create_child(system.get_hierarchy())
180 self.short_name = self.long_name =
"State_"+str(state_index)
181 self.hier.set_name(self.short_name)
182 self.molecules = IMP.pmi.tools.OrderedDict()
185 self._protocol_output = []
186 for p
in system._protocol_output:
187 self._add_protocol_output(p, system)
195 return self.system.__repr__()+
'.'+self.hier.get_name()
197 def _add_protocol_output(self, p, system):
198 state = p._add_state(self)
199 self._protocol_output.append((p, state))
200 state.model = system.model
201 state.prot = self.hier
204 """Return a dictionary where key is molecule name and value
205 are the list of all copies of that molecule in setup order"""
206 return self.molecules
209 """Access a molecule by name and copy number
210 @param name The molecule name used during setup
211 @param copy_num The copy number based on input order.
212 Default: 0. Set to 'all' to get all copies
214 if name
not in self.molecules:
215 raise KeyError(
"Could not find molname %s" % name)
216 if copy_num ==
'all':
217 return self.molecules[name]
219 return self.molecules[name][copy_num]
222 """Create a new Molecule within this State
223 @param name the name of the molecule (string);
224 it must not be already used
225 @param sequence sequence (string)
226 @param chain_id Chain ID to assign to this molecule
229 if name
in self.molecules:
230 raise ValueError(
'Cannot use a molecule name already used')
232 mol =
Molecule(self,name,sequence,chain_id,copy_num=0,is_nucleic=is_nucleic)
233 self.molecules[name] = [mol]
236 def get_hierarchy(self):
239 def get_number_of_copies(self,molname):
240 return len(self.molecules[molname])
242 def _register_copy(self,molecule):
243 molname = molecule.get_hierarchy().get_name()
244 self.molecules[molname].append(molecule)
247 """call build on all molecules (automatically makes clones)"""
249 for molname
in self.molecules:
250 for mol
in reversed(self.molecules[molname]):
257 _PDBElement = namedtuple(
'PDBElement', [
'offset',
'filename',
'chain_id'])
259 class _RepresentationHandler(object):
260 """Handle PMI representation and use it to populate that of any attached
261 ProtocolOutput objects"""
262 def __init__(self, name, pos, pdb_elements):
265 self.last_index =
None
266 self.last_pdb_index =
None
267 self.pdb_for_residue = {}
268 for residues, pdb
in pdb_elements:
270 self.pdb_for_residue[r.get_index()] = pdb
272 def _get_pdb(self, h):
273 """Return a PDBElement if the given hierarchy was read from a
277 return self.pdb_for_residue.get(rind,
None)
279 def __call__(self, res):
280 """Handle a single residue"""
281 if len(self.pos) == 0:
284 pi = h.get_particle_index()
286 if self.last_index
is None or pi != self.last_index:
287 pdb = self._get_pdb(h)
292 fragi = frag.get_particle_index()
294 if self.last_pdb_index
is not None \
295 and self.last_pdb_index == fragi:
297 self.last_pdb_index = fragi
298 indices = frag.get_residue_indexes()
299 for p, state
in self.pos:
300 p.add_pdb_element(state, self.name,
301 indices[0], indices[-1], pdb.offset,
302 pdb.filename, pdb.chain_id, frag)
305 indices = frag.get_residue_indexes()
306 for p, state
in self.pos:
307 p.add_bead_element(state, self.name,
308 indices[0], indices[-1], 1, h)
311 for p, state
in self.pos:
312 p.add_bead_element(state, self.name, resind, resind, 1, h)
314 raise TypeError(
"Unhandled hierarchy %s" % str(h))
320 """Stores a named protein chain.
321 This class is constructed from within the State class.
322 It wraps an IMP.atom.Molecule and IMP.atom.Copy.
323 Structure is read using this class.
324 Resolutions and copies can be registered, but are only created
325 when build() is called
328 def __init__(self,state,name,sequence,chain_id,copy_num,mol_to_clone=None,is_nucleic=None):
329 """The user should not call this directly; instead call State::create_molecule()
330 @param state The parent PMI State
331 @param name The name of the molecule (string)
332 @param sequence Sequence (string)
333 @param chain_id The chain of this molecule
334 @param copy_num Store the copy number
335 @param mol_to_clone The original molecule (for cloning ONLY)
336 @note It's expected that you will not use this constructor directly,
337 but rather create a Molecule with pmi::State::create_molecule()
340 self.model = state.get_hierarchy().get_model()
342 self.sequence = sequence
344 self.mol_to_clone = mol_to_clone
345 self.is_nucleic=is_nucleic
346 self.representations = []
347 self._pdb_elements = []
348 self._represented = IMP.pmi.tools.OrderedSet()
349 self.coord_finder = _FindCloseStructure()
350 self._ideal_helices = []
353 self.hier = self._create_child(self.state.get_hierarchy())
354 self.hier.set_name(name)
356 self._name_with_copy =
"%s.%d" % (name, copy_num)
359 self.chain.set_sequence(self.sequence)
362 for ns,s
in enumerate(sequence):
364 self.residues.append(r)
372 return self.state.__repr__()+
'.'+self.
get_name()+
'.'+ \
375 def __getitem__(self,val):
376 if isinstance(val,int):
377 return self.residues[val]
378 elif isinstance(val,str):
379 return self.residues[int(val)-1]
380 elif isinstance(val,slice):
381 return IMP.pmi.tools.OrderedSet(self.residues[val])
383 raise TypeError(
"Indexes must be int or str")
386 """Return the IMP Hierarchy corresponding to this Molecule"""
390 """Return this Molecule name"""
391 return self.hier.get_name()
394 """Return the State containing this Molecule"""
398 """Returns list of OrderedSets with requested ideal helices"""
399 return self._ideal_helices
402 """get residue range from a to b, inclusive.
403 Use integers to get 0-indexing, or strings to get PDB-indexing"""
404 if isinstance(a,int)
and isinstance(b,int)
and isinstance(stride,int):
405 return IMP.pmi.tools.OrderedSet(self.residues[a:b+1:stride])
406 elif isinstance(a,str)
and isinstance(b,str)
and isinstance(stride,int):
407 return IMP.pmi.tools.OrderedSet(self.residues[int(a)-1:int(b):stride])
409 raise TypeError(
"Range ends must be int or str. "
410 "Stride must be int.")
413 """ Return all modeled TempResidues as a set"""
414 all_res = IMP.pmi.tools.OrderedSet(self.residues)
418 """ Return set of TempResidues that have representation"""
419 return self._represented
422 """ Return a set of TempResidues that have associated structure coordinates """
423 atomic_res = IMP.pmi.tools.OrderedSet()
424 for res
in self.residues:
425 if res.get_has_structure():
430 """ Return a set of TempResidues that don't have associated structure coordinates """
431 non_atomic_res = IMP.pmi.tools.OrderedSet()
432 for res
in self.residues:
433 if not res.get_has_structure():
434 non_atomic_res.add(res)
435 return non_atomic_res
438 """Create a new Molecule with the same name and sequence but a higher copy number.
439 Returns the Molecule. No structure or representation will be copied!
440 @param chain_id Chain ID of the new molecule
443 copy_num=self.state.get_number_of_copies(self.
get_name()))
444 self.state._register_copy(mol)
448 """Create a Molecule clone (automatically builds same structure and representation)
449 @param chain_id If you want to set the chain ID of the copy to something
450 @note You cannot add structure or representations to a clone!
453 copy_num=self.state.get_number_of_copies(self.
get_name()),
455 self.state._register_copy(mol)
459 offset=0,model_num=
None,ca_only=
False,
461 """Read a structure and store the coordinates.
462 @eturn the atomic residues (as a set)
463 @param pdb_fn The file to read
464 @param chain_id Chain ID to read
465 @param res_range Add only a specific set of residues from the PDB file.
466 res_range[0] is the starting and res_range[1] is the
467 ending residue index.
468 @param offset Apply an offset to the residue indexes of the PDB
469 file. This number is added to the PDB sequence.
470 @param model_num Read multi-model PDB and return that model
471 @param ca_only Only read the CA positions from the PDB file
472 @param soft_check If True, it only warns if there are sequence
473 mismatches between the PDB and the Molecule (FASTA)
474 sequence, and uses the sequence from the PDB.
475 If False (Default), it raises an error when there
476 are sequence mismatches.
477 @note If you are adding structure without a FASTA file, set soft_check
480 if self.mol_to_clone
is not None:
481 raise ValueError(
'You cannot call add_structure() for a clone')
486 rhs = system_tools.get_structure(self.model,pdb_fn,chain_id,res_range,offset,ca_only=ca_only)
487 self.coord_finder.add_residues(rhs)
489 if len(self.residues)==0:
490 print(
"WARNING: Substituting PDB residue type with FASTA residue type. Potentially dangerous.")
493 self._pdb_elements.append((rhs,
494 _PDBElement(offset=offset, filename=pdb_fn, chain_id=chain_id)))
497 atomic_res = IMP.pmi.tools.OrderedSet()
498 for nrh,rh
in enumerate(rhs):
499 pdb_idx = rh.get_index()
500 raw_idx = pdb_idx - 1
503 while len(self.residues)<pdb_idx:
504 r =
TempResidue(self,
'A',len(self.residues)+1,len(self.residues))
505 self.residues.append(r)
508 internal_res = self.residues[raw_idx]
509 if len(self.sequence)<raw_idx:
511 internal_res.set_structure(rh,soft_check)
512 atomic_res.add(internal_res)
514 self.chain.set_sequence(self.sequence)
520 bead_extra_breaks=[],
521 bead_ca_centers=
True,
522 bead_default_coord=[0,0,0],
523 density_residues_per_component=
None,
525 density_force_compute=
False,
526 density_voxel_size=1.0,
527 setup_particles_as_densities=
False,
530 """Set the representation for some residues. Some options (beads, ideal helix)
531 operate along the backbone. Others (density options) are volumetric.
532 Some of these you can combine e.g., beads+densities or helix+densities
533 See @ref pmi_resolution
534 @param residues Set of PMI TempResidues for adding the representation.
535 Can use Molecule slicing to get these, e.g. mol[a:b]+mol[c:d]
536 If None, will select all residues for this Molecule.
537 @param resolutions Resolutions for beads representations.
538 If structured, will average along backbone, breaking at sequence breaks.
539 If unstructured, will just create beads.
540 Pass an integer or list of integers
541 @param bead_extra_breaks Additional breakpoints for splitting beads.
542 The value can be the 0-ordered position, after which it'll insert the break.
543 Alternatively pass PDB-style (1-ordered) indices as a string.
544 I.e., bead_extra_breaks=[5,25] is the same as ['6','26']
545 @param bead_ca_centers Set to True if you want the resolution=1 beads to be at CA centers
546 (otherwise will average atoms to get center). Defaults to True.
547 @param bead_default_coord Advanced feature. Normally beads are placed at the nearest structure.
548 If no structure provided (like an all bead molecule), the beads go here.
549 @param density_residues_per_component Create density (Gaussian Mixture Model)
550 for these residues. Must also supply density_prefix
551 @param density_prefix Prefix (assuming '.txt') to read components from or write to.
552 If exists, will read unless you set density_force_compute=True.
553 Will also write map (prefix+'.mrc').
554 Must also supply density_residues_per_component.
555 @param density_force_compute Set true to force overwrite density file.
556 @param density_voxel_size Advanced feature. Set larger if densities taking too long to rasterize.
557 Set to 0 if you don't want to create the MRC file
558 @param setup_particles_as_densities Set to True if you want each particle to be its own density.
559 Useful for all-atom models or flexible beads.
560 Mutually exclusive with density_ options
561 @param ideal_helix Create idealized helix structures for these residues at resolution 1.
562 Any other resolutions passed will be coarsened from there.
563 Resolution 0 will not work, you may have to use MODELLER to do that (for now).
564 @param color the color applied to the hierarchies generated.
565 Format options: tuple (r,g,b) with values 0 to 1
566 or float (from 0 to 1, a map from Blue to Green to Red)
567 or IMP.display.Color object
569 @note You cannot call add_representation multiple times for the
574 if self.mol_to_clone
is not None:
575 raise ValueError(
'You cannot call add_representation() for a clone.'
576 ' Maybe use a copy instead.')
580 res = IMP.pmi.tools.OrderedSet(self.residues)
582 res = IMP.pmi.tools.OrderedSet(self.residues)
584 res = IMP.pmi.tools.OrderedSet([residues])
585 elif hasattr(residues,
'__iter__'):
587 raise Exception(
'You passed an empty set to add_representation')
588 if type(residues)
is IMP.pmi.tools.OrderedSet
and type(next(iter(residues)))
is TempResidue:
590 elif type(residues)
is set
and type(next(iter(residues)))
is TempResidue:
591 res = IMP.pmi.tools.OrderedSet(residues)
592 elif type(residues)
is list
and type(residues[0])
is TempResidue:
593 res = IMP.pmi.tools.OrderedSet(residues)
595 raise Exception(
"You passed an iteratible of something other than TempResidue",res)
597 raise Exception(
"add_representation: you must pass a set of residues or nothing(=all residues)")
600 ov = res & self._represented
603 self._represented|=res
606 if not hasattr(resolutions,
'__iter__'):
607 if type(resolutions)
is int:
608 resolutions = [resolutions]
610 raise Exception(
"you tried to pass resolutions that are not int or list-of-int")
611 if len(resolutions)>1
and not ideal_helix:
613 if not r.get_has_structure():
614 raise Exception(
'You are creating multiple resolutions for '
615 'unstructured regions. This will have unexpected results.')
618 if density_residues_per_component
or density_prefix:
619 if not density_residues_per_component
and density_prefix:
620 raise Exception(
'If requesting density, must provide '
621 'density_residues_per_component AND density_prefix')
622 if density_residues_per_component
and setup_particles_as_densities:
623 raise Exception(
'Cannot create both volumetric density '
624 '(density_residues_per_component) AND '
625 'individual densities (setup_particles_as_densities) '
626 'in the same representation')
627 if len(resolutions)>1
and setup_particles_as_densities:
628 raise Exception(
'You have multiple bead resolutions but are attempting to '
629 'set them all up as individual Densities. '
630 'This could have unexpected results.')
635 raise Exception(
"For ideal helices, cannot build resolution 0: "
636 "you have to do that in MODELLER")
637 if 1
not in resolutions:
638 resolutions = [1] + list(resolutions)
639 self._ideal_helices.append(res)
643 if r.get_molecule()!=self:
644 raise Exception(
'You are adding residues from a different molecule to',self.__repr__())
648 for b
in bead_extra_breaks:
650 breaks.append(int(b)-1)
654 self.representations.append(_Representation(res,
659 density_residues_per_component,
661 density_force_compute,
663 setup_particles_as_densities,
667 def _all_protocol_output(self):
668 return self.state._protocol_output
671 """Create all parts of the IMP hierarchy
672 including Atoms, Residues, and Fragments/Representations and,
674 Will only build requested representations.
675 @note Any residues assigned a resolution must have an IMP.atom.Residue
676 hierarchy containing at least a CAlpha. For missing residues,
677 these can be constructed from the PDB file.
681 name=self.hier.get_name()
682 for po, state
in self._all_protocol_output():
683 po.create_component(state, name,
True,
684 asym_name=self._name_with_copy)
685 po.add_component_sequence(state, name, self.sequence,
686 asym_name=self._name_with_copy)
688 if self.mol_to_clone
is not None:
689 for nr,r
in enumerate(self.mol_to_clone.residues):
690 if r.get_has_structure():
691 clone = IMP.atom.create_clone(r.get_hierarchy())
692 self.residues[nr].set_structure(
694 for old_rep
in self.mol_to_clone.representations:
695 new_res = IMP.pmi.tools.OrderedSet()
696 for r
in old_rep.residues:
697 new_res.add(self.residues[r.get_internal_index()])
698 self._represented.add(self.residues[r.get_internal_index()])
699 new_rep = _Representation(new_res,
700 old_rep.bead_resolutions,
701 old_rep.bead_extra_breaks,
702 old_rep.bead_ca_centers,
703 old_rep.bead_default_coord,
704 old_rep.density_residues_per_component,
705 old_rep.density_prefix,
707 old_rep.density_voxel_size,
708 old_rep.setup_particles_as_densities,
711 self.representations.append(new_rep)
712 self.coord_finder = self.mol_to_clone.coord_finder
715 no_rep = [r
for r
in self.residues
if r
not in self._represented]
717 print(
'WARNING: Residues without representation in molecule',
718 self.
get_name(),
':',system_tools.resnums2str(no_rep))
721 for rep
in self.representations:
723 _build_ideal_helix(self.model,rep.residues,self.coord_finder)
727 for rep
in self.representations:
728 built_reps += system_tools.build_representation(
729 self, rep, self.coord_finder)
733 for br
in built_reps:
734 self.hier.add_child(br)
738 rephandler = _RepresentationHandler(self._name_with_copy,
739 list(self._all_protocol_output()),
741 for res
in self.residues:
742 idx = res.get_index()
747 residue_index=res.get_index(),
748 resolution=1).get_selected_particles()
761 self._represented = IMP.pmi.tools.OrderedSet([a
for a
in self._represented])
766 """Helpful utility for getting particles at all resolutions from this molecule.
767 Can optionally pass a set of residue indexes"""
769 raise Exception(
"Cannot get all resolutions until you build the Molecule")
770 if residue_indexes
is None:
771 residue_indexes = [r.get_index()
for r
in self.
get_residues()]
772 ps = IMP.pmi.tools.select_at_all_resolutions(self.
get_hierarchy(),
773 residue_indexes=residue_indexes)
778 class _Representation(object):
779 """Private class just to store a representation request"""
786 density_residues_per_component,
788 density_force_compute,
790 setup_particles_as_densities,
793 self.residues = residues
794 self.bead_resolutions = bead_resolutions
795 self.bead_extra_breaks = bead_extra_breaks
796 self.bead_ca_centers = bead_ca_centers
797 self.bead_default_coord = bead_default_coord
798 self.density_residues_per_component = density_residues_per_component
799 self.density_prefix = density_prefix
800 self.density_force_compute = density_force_compute
801 self.density_voxel_size = density_voxel_size
802 self.setup_particles_as_densities = setup_particles_as_densities
803 self.ideal_helix = ideal_helix
806 class _FindCloseStructure(object):
807 """Utility to get the nearest observed coordinate"""
810 def add_residues(self,residues):
817 self.coords.append([idx,xyz])
820 self.coords.append([idx,xyz])
822 raise(
"_FindCloseStructure: wrong selection")
824 self.coords.sort(key=itemgetter(0))
825 def find_nearest_coord(self,query):
828 keys = [r[0]
for r
in self.coords]
829 pos = bisect_left(keys,query)
832 elif pos == len(self.coords):
833 ret = self.coords[-1]
835 before = self.coords[pos - 1]
836 after = self.coords[pos]
837 if after[0] - query < query - before[0]:
844 """A dictionary-like wrapper for reading and storing sequence data"""
845 def __init__(self,fasta_fn,name_map=None):
846 """read a fasta file and extract all the requested sequences
847 @param fasta_fn sequence file
848 @param name_map dictionary mapping the fasta name to final stored name
850 self.sequences = IMP.pmi.tools.OrderedDict()
851 self.read_sequences(fasta_fn,name_map)
853 return len(self.sequences)
854 def __contains__(self,x):
855 return x
in self.sequences
856 def __getitem__(self,key):
859 allseqs = list(self.sequences.keys())
860 return self.sequences[allseqs[key]]
862 raise Exception(
"You tried to access sequence num",key,
"but there's only",len(self.sequences.keys()))
864 return self.sequences[key]
866 return self.sequences.__iter__()
869 for s
in self.sequences:
870 ret +=
'%s\t%s\n'%(s,self.sequences[s])
872 def read_sequences(self,fasta_fn,name_map=None):
875 with open(fasta_fn,
'r') as fh:
876 for (num, line)
in enumerate(fh):
877 if line.startswith(
'>'):
879 self.sequences[code] = seq.strip(
'*')
880 code = line.rstrip()[1:]
881 if name_map
is not None:
883 code = name_map[code]
892 "Found FASTA sequence before first header at line %d: %s" % (num + 1, line))
895 self.sequences[code] = seq.strip(
'*')
899 """Data_structure for reading and storing sequence data from pdb"""
900 def __init__(self,model,pdb_fn,name_map=None):
901 """read a pdb file and returns all sequences for each contiguous fragment
903 @param name_map dictionary mapping the pdb chain id to final stored name
915 self.sequences = IMP.pmi.tools.OrderedDict()
916 self.read_sequences(pdb_fn,name_map)
923 def read_sequences(self,pdb_fn,name_map):
924 read_file = IMP.atom.read_pdb
925 if pdb_fn.endswith(
'.cif'):
926 read_file = IMP.atom.read_mmcif
928 cs=IMP.atom.get_by_type(t,IMP.atom.CHAIN_TYPE)
936 print(
"Chain ID %s not in name_map, skipping" % id)
938 rs=IMP.atom.get_by_type(c,IMP.atom.RESIDUE_TYPE)
945 isprotein=dr.get_is_protein()
946 isrna=dr.get_is_rna()
947 isdna=dr.get_is_dna()
951 rids_olc_dict[rid]=olc
953 if dr.get_residue_type() == IMP.atom.DADE: olc=
"A"
954 if dr.get_residue_type() == IMP.atom.DURA: olc=
"U"
955 if dr.get_residue_type() == IMP.atom.DCYT: olc=
"C"
956 if dr.get_residue_type() == IMP.atom.DGUA: olc=
"G"
957 if dr.get_residue_type() == IMP.atom.DTHY: olc=
"T"
959 rids_olc_dict[rid]=olc
961 if dr.get_residue_type() == IMP.atom.ADE: olc=
"A"
962 if dr.get_residue_type() == IMP.atom.URA: olc=
"U"
963 if dr.get_residue_type() == IMP.atom.CYT: olc=
"C"
964 if dr.get_residue_type() == IMP.atom.GUA: olc=
"G"
965 if dr.get_residue_type() == IMP.atom.THY: olc=
"T"
967 rids_olc_dict[rid]=olc
968 group_rids=self.group_indexes(rids)
969 contiguous_sequences=IMP.pmi.tools.OrderedDict()
970 for group
in group_rids:
972 for i
in range(group[0],group[1]+1):
973 sequence_fragment+=rids_olc_dict[i]
974 contiguous_sequences[group]=sequence_fragment
975 self.sequences[id]=contiguous_sequences
981 def group_indexes(self,indexes):
982 from itertools
import groupby
984 for k, g
in groupby(enumerate(indexes),
lambda x:x[0]-x[1]):
985 group = [x[1]
for x
in g]
986 ranges.append((group[0], group[-1]))
991 '''This function computes and prints the alignment between the
992 fasta file and the pdb sequence, computes the offsets for each contiguous
994 @param fasta_sequences IMP.pmi.topology.Sequences object
995 @param pdb_sequences IMP.pmi.topology.PDBSequences object
996 @param show boolean default False, if True prints the alignments.
997 The input objects should be generated using map_name dictionaries such that fasta_id
998 and pdb_chain_id are mapping to the same protein name. It needs BioPython.
999 Returns a dictionary of offsets, organized by peptide range (group):
1000 example: offsets={"ProtA":{(1,10):1,(20,30):10}}'''
1001 from Bio
import pairwise2
1002 from Bio.pairwise2
import format_alignment
1003 from Bio.SubsMat
import MatrixInfo
as matlist
1004 matrix = matlist.blosum62
1006 raise Exception(
"Fasta sequences not type IMP.pmi.topology.Sequences")
1008 raise Exception(
"pdb sequences not type IMP.pmi.topology.PDBSequences")
1009 offsets=IMP.pmi.tools.OrderedDict()
1010 for name
in fasta_sequences.sequences:
1012 seq_fasta=fasta_sequences.sequences[name]
1013 if name
not in pdb_sequences.sequences:
1014 print(
"Fasta id %s not in pdb names, aligning against every pdb chain" % name)
1015 pdbnames=pdb_sequences.sequences.keys()
1018 for pdbname
in pdbnames:
1019 for group
in pdb_sequences.sequences[pdbname]:
1020 if group[1]-group[0]+1<7:
continue
1021 seq_frag_pdb=pdb_sequences.sequences[pdbname][group]
1023 print(
"########################")
1025 print(
"protein name",pdbname)
1026 print(
"fasta id", name)
1027 print(
"pdb fragment",group)
1028 align=pairwise2.align.localms(seq_fasta, seq_frag_pdb, 2, -1, -.5, -.1)[0]
1030 offset=a[3]+1-group[0]
1032 print(
"alignment sequence start-end",(a[3]+1,a[4]+1))
1033 print(
"offset from pdb to fasta index",offset)
1034 print(format_alignment(*a))
1035 if name
not in offsets:
1037 if group
not in offsets[pdbname]:
1038 offsets[pdbname][group]=offset
1040 if group
not in offsets[pdbname]:
1041 offsets[pdbname][group]=offset
1050 """Temporarily stores residue information, even without structure available."""
1052 def __init__(self,molecule,code,index,internal_index,is_nucleic=None):
1053 """setup a TempResidue
1054 @param molecule PMI Molecule to which this residue belongs
1055 @param code one-letter residue type code
1056 @param index PDB index
1057 @param internal_index The number in the sequence
1060 self.molecule = molecule
1061 self.rtype = IMP.pmi.tools.get_residue_type_from_one_letter_code(code,is_nucleic)
1062 self.pdb_index = index
1063 self.internal_index = internal_index
1067 self._structured =
False
1072 return str(self.state_index)+
"_"+self.molecule.get_name()+
"_"+str(self.copy_index)+
"_"+self.get_code()+str(self.get_index())
1074 return self.__str__()
1077 return (self.state_index, self.molecule, self.copy_index, self.rtype, self.pdb_index, self.internal_index)
1078 def __eq__(self,other):
1079 return type(other)==type(self)
and self.__key() == other.__key()
1081 return hash(self.__key())
1083 return self.pdb_index
1084 def get_internal_index(self):
1085 return self.internal_index
1088 def get_residue_type(self):
1090 def get_hierarchy(self):
1092 def get_molecule(self):
1093 return self.molecule
1094 def get_has_structure(self):
1095 return self._structured
1096 def set_structure(self,res,soft_check=False):
1097 if res.get_residue_type()!=self.get_residue_type():
1100 print(
'WARNING: Inconsistency between FASTA sequence and PDB sequence. FASTA type',\
1101 self.get_index(),self.hier.get_residue_type(),
1102 'and PDB type',res.get_residue_type())
1103 self.hier.set_residue_type((self.get_residue_type()))
1104 self.rtype = self.get_residue_type()
1106 raise Exception(
'ERROR: PDB residue index',self.get_index(),
'is',
1108 'and sequence residue is',self.get_code())
1110 for a
in res.get_children():
1111 self.hier.add_child(a)
1113 a.get_particle().set_name(
'Atom %s of residue %i'%(atype.__str__().strip(
'"'),
1114 self.hier.get_index()))
1115 self._structured =
True
1118 """Automatically setup Sytem and Degrees of Freedom with a formatted text file.
1119 The file is read in and each part of the topology is stored as a
1120 ComponentTopology object for input into IMP::pmi::macros::BuildSystem.
1121 The topology file should be in a simple pipe-delimited format:
1123 |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|flags|
1124 |Rpb1 |blue |1WCM.fasta|1WCM:A|1WCM.pdb|A|1,1140 |0|10|0|1|1,3|1||
1125 |Rpb1 |blue |1WCM.fasta|1WCM:A|1WCM.pdb|A|1141,1274|0|10|0|2|1,3|1||
1126 |Rpb1 |blue |1WCM.fasta|1WCM:A|1WCM.pdb|A|1275,END |0|10|0|3|1,3|1||
1127 |Rpb2 |red |1WCM.fasta|1WCM:B|1WCM.pdb|B|all |0|10|0|4|2,3|2||
1128 |Rpb2.1 |green |1WCM.fasta|1WCM:B|1WCM.pdb|B|all |0|10|0|4|2,3|2||
1132 These are the fields you can enter:
1133 - `component_name`: Name of the component (chain). Serves as the parent
1134 hierarchy for this structure.
1135 - `color`: The color used in the output RMF file. Uses chimera names or R,G,B values
1136 - `fasta_fn`: Name of FASTA file containing this component.
1137 - `fasta_id`: String found in FASTA sequence header line. The sequence read
1138 from the file is assumed to be a protein sequence. If it should instead
1139 be treated as RNA or DNA, add an ',RNA' or ',DNA' suffix. For example,
1140 a `fasta_id` of 'myseq,RNA' will read the sequence 'myseq' from the
1141 FASTA file and treat it as RNA.
1142 - `pdb_fn`: Name of PDB or mmCIF file with coordinates (if available).
1143 If left empty, will set up as BEADS (you can also specify "BEADS")
1144 Can also write "IDEAL_HELIX".
1145 - `chain`: Chain ID of this domain in the PDB file.
1146 - `residue_range`: Comma delimited pair defining range.
1147 Can leave empty or use 'all' for entire sequence from PDB file.
1148 The second item in the pair can be END to select the last residue in the
1150 - `pdb_offset`: Offset to sync PDB residue numbering with FASTA numbering.
1151 - `bead_size`: The size (in residues) of beads used to model areas not
1152 covered by PDB coordinates. These will be automatically built.
1153 - `em_residues`: The number of Gaussians used to model the electron
1154 density of this domain. Set to zero if no EM fitting will be done.
1155 The GMM files will be written to <gmm_dir>/<component_name>_<em_res>.txt (and .mrc)
1156 - `rigid_body`: Leave empty if this object is not in a rigid body.
1157 Otherwise, this is a number corresponding to the rigid body containing
1158 this object. The number itself is just used for grouping things.
1159 - `super_rigid_body`: Like a rigid_body, except things are only occasionally rigid
1160 - `chain_of_super_rigid_bodies` For a polymer, create SRBs from groups.
1161 - `flags` additional flags for advanced options
1162 @note All filenames are relative to the paths specified in the constructor.
1171 @param topology_file Pipe-delimited file specifying the topology
1172 @param pdb_dir Relative path to the pdb directory
1173 @param fasta_dir Relative path to the fasta directory
1174 @param gmm_dir Relative path to the GMM directory
1176 self.topology_file = topology_file
1177 self.molecules = IMP.pmi.tools.OrderedDict()
1178 self.pdb_dir = pdb_dir
1179 self.fasta_dir = fasta_dir
1180 self.gmm_dir = gmm_dir
1181 self._components = self.
read(topology_file)
1185 "Use 'get_components()' instead of 'component_list'.")
1186 def __get_component_list(self):
return self._components
1187 component_list = property(__get_component_list)
1189 def write_topology_file(self,outfile):
1190 with open(outfile,
"w")
as f:
1191 f.write(
"|molecule_name|color|fasta_fn|fasta_id|pdb_fn|chain|"
1192 "residue_range|pdb_offset|bead_size|em_residues_per_gaussian|"
1193 "rigid_body|super_rigid_body|chain_of_super_rigid_bodies|\n")
1194 for c
in self._components:
1195 output = c.get_str()+
'\n'
1200 """ Return list of ComponentTopologies for selected components
1201 @param topology_list List of indices to return"""
1202 if topology_list ==
"all":
1203 topologies = self._components
1206 for i
in topology_list:
1207 topologies.append(self._components[i])
1210 def get_molecules(self):
1211 return self.molecules
1213 def read(self, topology_file, append=False):
1214 """Read system components from topology file. append=False will erase
1215 current topology and overwrite with new
1218 is_directories =
False
1223 with open(topology_file)
as infile:
1225 if line.lstrip()==
"" or line[0]==
"#":
1227 elif line.split(
'|')[1].strip()
in (
"molecule_name"):
1229 is_directories =
False
1232 elif line.split(
'|')[1] ==
"component_name":
1235 "Old-style topology format (using "
1236 "|component_name|) is deprecated. Please switch to "
1237 "the new-style format (using |molecule_name|)\n")
1239 is_directories =
False
1241 elif line.split(
'|')[1] ==
"directories":
1243 "Setting directories in the topology file "
1244 "is deprecated. Please do so through the "
1245 "TopologyReader constructor. Note that new-style "
1246 "paths are relative to the current working "
1247 "directory, not the topology file.\n")
1248 is_directories =
True
1249 elif is_directories:
1250 fields = line.split(
'|')
1251 setattr(self, fields[1],
1254 new_component = self._parse_line(line, linenum, old_format)
1255 self._components.append(new_component)
1257 return self._components
1259 def _parse_line(self, component_line, linenum, old_format):
1260 """Parse a line of topology values and matches them to their key.
1261 Checks each value for correct syntax
1262 Returns a list of Component objects
1266 values = [s.strip()
for s
in component_line.split(
'|')]
1271 c.molname = values[1]
1273 c._domain_name = values[2]
1276 names = values[1].split(
'.')
1278 c.molname = names[0]
1281 c.molname = names[0]
1282 c.copyname = names[1]
1284 c.molname = names[0]
1285 c.copyname = names[1]
1286 errors.append(
"Molecule name should be <molecule.copyID>")
1287 errors.append(
"For component %s line %d " % (c.molname,linenum))
1288 c._domain_name = c.molname +
'.' + c.copyname
1289 colorfields = values[2].split(
',')
1290 if len(colorfields)==3:
1291 c.color = [float(x)
for x
in colorfields]
1292 if any([x>1
for x
in c.color]):
1293 c.color=[x/255
for x
in c.color]
1296 c._orig_fasta_file = values[3]
1297 c.fasta_file = values[3]
1298 fasta_field = values[4].split(
",")
1299 c.fasta_id = fasta_field[0]
1301 if len(fasta_field) > 1:
1302 c.fasta_flag = fasta_field[1]
1303 c._orig_pdb_input = values[5]
1304 pdb_input = values[5]
1305 tmp_chain = values[6]
1308 bead_size = values[9]
1311 rbs = srbs = csrbs =
''
1317 if c.molname
not in self.molecules:
1318 self.molecules[c.molname] = _TempMolecule(c)
1321 c._orig_fasta_file = self.molecules[c.molname].orig_component.fasta_file
1322 c.fasta_id = self.molecules[c.molname].orig_component.fasta_id
1323 self.molecules[c.molname].add_component(c,c.copyname)
1326 c.fasta_file = os.path.join(self.fasta_dir,c._orig_fasta_file)
1328 errors.append(
"PDB must have BEADS, IDEAL_HELIX, or filename")
1329 errors.append(
"For component %s line %d is not correct"
1330 "|%s| was given." % (c.molname,linenum,pdb_input))
1331 elif pdb_input
in (
"IDEAL_HELIX",
"BEADS"):
1332 c.pdb_file = pdb_input
1334 c.pdb_file = os.path.join(self.pdb_dir,pdb_input)
1337 if len(tmp_chain)==1
or len(tmp_chain)==2:
1340 errors.append(
"PDB Chain identifier must be one or two characters.")
1341 errors.append(
"For component %s line %d is not correct"
1342 "|%s| was given." % (c.molname,linenum,tmp_chain))
1346 if rr.strip()==
'all' or str(rr)==
"":
1347 c.residue_range =
None
1348 elif len(rr.split(
','))==2
and self._is_int(rr.split(
',')[0])
and (self._is_int(rr.split(
',')[1])
or rr.split(
',')[1] ==
'END'):
1350 c.residue_range = (int(rr.split(
',')[0]), rr.split(
',')[1])
1351 if c.residue_range[1] !=
'END':
1352 c.residue_range = (c.residue_range[0], int(c.residue_range[1]))
1354 if old_format
and c.residue_range[1] == -1:
1355 c.residue_range = (c.residue_range[0],
'END')
1357 errors.append(
"Residue Range format for component %s line %d is not correct" % (c.molname, linenum))
1358 errors.append(
"Correct syntax is two comma separated integers: |start_res, end_res|. end_res can also be END to select the last residue in the chain. |%s| was given." % rr)
1359 errors.append(
"To select all residues, indicate |\"all\"|")
1362 if self._is_int(offset):
1363 c.pdb_offset=int(offset)
1364 elif len(offset)==0:
1367 errors.append(
"PDB Offset format for component %s line %d is not correct" % (c.molname, linenum))
1368 errors.append(
"The value must be a single integer. |%s| was given." % offset)
1371 if self._is_int(bead_size):
1372 c.bead_size=int(bead_size)
1373 elif len(bead_size)==0:
1376 errors.append(
"Bead Size format for component %s line %d is not correct" % (c.molname, linenum))
1377 errors.append(
"The value must be a single integer. |%s| was given." % bead_size)
1380 if self._is_int(emg):
1382 c.density_prefix = os.path.join(self.gmm_dir,c.get_unique_name())
1383 c.gmm_file = c.density_prefix+
'.txt'
1384 c.mrc_file = c.density_prefix+
'.gmm'
1386 c.em_residues_per_gaussian=int(emg)
1388 c.em_residues_per_gaussian = 0
1390 c.em_residues_per_gaussian = 0
1392 errors.append(
"em_residues_per_gaussian format for component "
1393 "%s line %d is not correct" % (c.molname, linenum))
1394 errors.append(
"The value must be a single integer. |%s| was given." % emg)
1398 if not self._is_int(rbs):
1399 errors.append(
"rigid bodies format for component "
1400 "%s line %d is not correct" % (c.molname, linenum))
1401 errors.append(
"Each RB must be a single integer, or empty. "
1402 "|%s| was given." % rbs)
1403 c.rigid_body = int(rbs)
1407 srbs = srbs.split(
',')
1409 if not self._is_int(i):
1410 errors.append(
"super rigid bodies format for component "
1411 "%s line %d is not correct" % (c.molname, linenum))
1412 errors.append(
"Each SRB must be a single integer. |%s| was given." % srbs)
1413 c.super_rigid_bodies = srbs
1417 if not self._is_int(csrbs):
1418 errors.append(
"em_residues_per_gaussian format for component "
1419 "%s line %d is not correct" % (c.molname, linenum))
1420 errors.append(
"Each CSRB must be a single integer. |%s| was given." % csrbs)
1421 c.chain_of_super_rigid_bodies = csrbs
1425 raise ValueError(
"Fix Topology File syntax errors and rerun: " \
1426 +
"\n".join(errors))
1432 """Change the GMM dir"""
1433 self.gmm_dir = gmm_dir
1434 for c
in self._components:
1435 c.gmm_file = os.path.join(self.gmm_dir,c.get_unique_name()+
".txt")
1436 c.mrc_file = os.path.join(self.gmm_dir,c.get_unique_name()+
".mrc")
1437 print(
'new gmm',c.gmm_file)
1440 """Change the PDB dir"""
1441 self.pdb_dir = pdb_dir
1442 for c
in self._components:
1443 if not c._orig_pdb_input
in (
"",
"None",
"IDEAL_HELIX",
"BEADS"):
1444 c.pdb_file = os.path.join(self.pdb_dir,c._orig_pdb_input)
1447 """Change the FASTA dir"""
1448 self.fasta_dir = fasta_dir
1449 for c
in self._components:
1450 c.fasta_file = os.path.join(self.fasta_dir,c._orig_fasta_file)
1452 def _is_int(self, s):
1456 return float(s).is_integer()
1461 """Return list of lists of rigid bodies (as domain name)"""
1462 rbl = defaultdict(list)
1463 for c
in self._components:
1465 rbl[c.rigid_body].append(c.get_unique_name())
1469 """Return list of lists of super rigid bodies (as domain name)"""
1470 rbl = defaultdict(list)
1471 for c
in self._components:
1472 for rbnum
in c.super_rigid_bodies:
1473 rbl[rbnum].append(c.get_unique_name())
1477 """Return list of lists of chains of super rigid bodies (as domain name)"""
1478 rbl = defaultdict(list)
1479 for c
in self._components:
1480 for rbnum
in c.chain_of_super_rigid_bodies:
1481 rbl[rbnum].append(c.get_unique_name())
1484 class _TempMolecule(object):
1485 """Store the Components and any requests for copies"""
1487 self.molname = init_c.molname
1489 self.domains = IMP.pmi.tools.OrderedDefaultDict(list)
1490 self.add_component(init_c,init_c.copyname)
1491 self.orig_copyname = init_c.copyname
1492 self.orig_component = self.domains[init_c.copyname][0]
1493 def add_component(self,component,copy_id):
1494 self.domains[copy_id].append(component)
1495 component.domainnum = len(self.domains[copy_id])-1
1497 return ','.join(
'%s:%i'%(k,len(self.domains[k]))
for k
in self.domains)
1499 class _Component(object):
1500 """Stores the components required to build a standard IMP hierarchy
1501 using IMP.pmi.BuildModel()
1505 self.copyname =
None
1507 self.fasta_file =
None
1508 self._orig_fasta_file =
None
1509 self.fasta_id =
None
1510 self.fasta_flag =
None
1511 self.pdb_file =
None
1512 self._orig_pdb_input =
None
1514 self.residue_range =
None
1517 self.em_residues_per_gaussian = 0
1520 self.density_prefix =
''
1522 self.rigid_body =
None
1523 self.super_rigid_bodies = []
1524 self.chain_of_super_rigid_bodies = []
1527 return ",".join(
"%s" % x
for x
in l)
1530 return self.get_str()
1533 return "%s.%s.%i"%(self.molname,self.copyname,self.domainnum)
1536 res_range = self.residue_range
1537 if self.residue_range
is None:
1540 if self.copyname!=
'':
1541 name +=
'.'+self.copyname
1542 if self.chain
is None:
1547 if isinstance(color, list):
1548 color=
','.join([str(x)
for x
in color])
1549 a=
'|'+
'|'.join([name,color,self._orig_fasta_file,self.fasta_id,
1550 self._orig_pdb_input,chain,self._l2s(list(res_range)),
1551 str(self.pdb_offset),str(self.bead_size),
1552 str(self.em_residues_per_gaussian),
1553 str(self.rigid_body)
if self.rigid_body
else '',
1554 self._l2s(self.super_rigid_bodies),
1555 self._l2s(self.chain_of_super_rigid_bodies)])+
'|'
1560 def __get_name(self):
return self.molname
1561 name = property(__get_name)
1565 "Use 'get_unique_name()' instead of 'domain_name'.")
1566 def __get_domain_name(self):
return self._domain_name
1567 domain_name = property(__get_domain_name)
1571 '''Extends the functionality of IMP.atom.Molecule'''
1573 def __init__(self,hierarchy):
1574 IMP.atom.Molecule.__init__(self,hierarchy)
1583 def get_extended_name(self):
1584 return self.get_name()+
"."+\
1585 str(self.get_copy_index())+\
1586 "."+str(self.get_state_index())
1588 def get_sequence(self):
1591 def get_residue_indexes(self):
1594 def get_residue_segments(self):
1595 return IMP.pmi.tools.Segments(self.get_residue_indexes())
1601 s=
'PMIMoleculeHierarchy '
1604 s+=
" "+
"State "+str(self.get_state_index())
1605 s+=
" "+
"N residues "+str(len(self.get_sequence()))
def __init__
setup a TempResidue
def build
call build on all molecules (automatically makes clones)
def set_pdb_dir
Change the PDB dir.
static bool get_is_setup(const IMP::ParticleAdaptor &p)
Hierarchy get_parent() const
Get the parent particle.
def get_atomic_residues
Return a set of TempResidues that have associated structure coordinates.
A decorator to associate a particle with a part of a protein/DNA/RNA.
Extends the functionality of IMP.atom.Molecule.
def get_residues
Return all modeled TempResidues as a set.
std::string get_unique_name(std::string templ)
Return a unique name produced from the string.
static bool get_is_setup(const IMP::ParticleAdaptor &p)
static Atom setup_particle(Model *m, ParticleIndex pi, Atom other)
def build
Build all states.
def __init__
read a fasta file and extract all the requested sequences
static XYZR setup_particle(Model *m, ParticleIndex pi)
def __init__
read a pdb file and returns all sequences for each contiguous fragment
def fasta_pdb_alignments
This function computes and prints the alignment between the fasta file and the pdb sequence...
def get_ideal_helices
Returns list of OrderedSets with requested ideal helices.
def get_chains_of_super_rigid_bodies
Return list of lists of chains of super rigid bodies (as domain name)
def __init__
The user should not call this directly; instead call State::create_molecule()
void handle_use_deprecated(std::string message)
def set_gmm_dir
Change the GMM dir.
def residue_range
get residue range from a to b, inclusive.
def get_molecule
Access a molecule by name and copy number.
def __init__
Define a new state.
def add_representation
Set the representation for some residues.
static State setup_particle(Model *m, ParticleIndex pi, unsigned int index)
def create_molecule
Create a new Molecule within this State.
def build
Create all parts of the IMP hierarchy including Atoms, Residues, and Fragments/Representations and...
static Residue setup_particle(Model *m, ParticleIndex pi, ResidueType t, int index, int insertion_code)
char get_one_letter_code(ResidueType c)
Get the 1-letter amino acid code from the residue type.
def get_non_atomic_residues
Return a set of TempResidues that don't have associated structure coordinates.
def get_name
Return this Molecule name.
def get_hierarchy
Return the IMP Hierarchy corresponding to this Molecule.
def get_components
Return list of ComponentTopologies for selected components.
Class for storing model, its restraints, constraints, and particles.
Stores a named protein chain.
static bool get_is_setup(Model *m, ParticleIndex pi)
A decorator for keeping track of copies of a molecule.
Select all non-alternative ATOM records.
static Hierarchy setup_particle(Model *m, ParticleIndex pi, ParticleIndexesAdaptor children=ParticleIndexesAdaptor())
Create a Hierarchy of level t by adding the needed attributes.
def set_fasta_dir
Change the FASTA dir.
The standard decorator for manipulating molecular structures.
Ints get_index(const ParticlesTemp &particles, const Subset &subset, const Subsets &excluded)
Data_structure for reading and storing sequence data from pdb.
def deprecated_method
Python decorator to mark a method as deprecated.
A decorator for a particle representing an atom.
std::string get_relative_path(std::string base, std::string relative)
Return a path to a file relative to another file.
A decorator for a particle with x,y,z coordinates.
def create_clone
Create a Molecule clone (automatically builds same structure and representation)
def add_structure
Read a structure and store the coordinates.
int get_state_index(Hierarchy h)
Walk up the hierarchy to find the current state.
def add_protocol_output
Capture details of the modeling protocol.
def get_molecules
Return a dictionary where key is molecule name and value are the list of all copies of that molecule ...
static Copy setup_particle(Model *m, ParticleIndex pi, Int number)
def read
Read system components from topology file.
def get_state
Return the State containing this Molecule.
A decorator for a residue.
General purpose algebraic and geometric methods that are expected to be used by a wide variety of IMP...
Automatically setup Sytem and Degrees of Freedom with a formatted text file.
The general base class for IMP exceptions.
def get_rigid_bodies
Return list of lists of rigid bodies (as domain name)
Associate an integer "state" index with a hierarchy node.
Residue get_residue(Atom d, bool nothrow=false)
Return the Residue containing this atom.
Class to handle individual particles of a Model object.
Stores a list of Molecules all with the same State index.
def get_represented
Return set of TempResidues that have representation.
Store info for a chain of a protein.
int get_copy_index(Hierarchy h)
Walk up the hierarchy to find the current copy index.
Python classes to represent, score, sample and analyze models.
A dictionary-like wrapper for reading and storing sequence data.
def create_copy
Create a new Molecule with the same name and sequence but a higher copy number.
Functionality for loading, creating, manipulating and scoring atomic structures.
std::string get_chain_id(Hierarchy h)
Walk up the hierarchy to determine the chain id.
def get_particles_at_all_resolutions
Helpful utility for getting particles at all resolutions from this molecule.
static Chain setup_particle(Model *m, ParticleIndex pi, std::string id)
A decorator for a molecule.
Select hierarchy particles identified by the biological name.
def get_number_of_states
Returns the total number of states generated.
def get_super_rigid_bodies
Return list of lists of super rigid bodies (as domain name)
Temporarily stores residue information, even without structure available.
def create_state
Makes and returns a new IMP.pmi.topology.State in this system.