You are reading the documentation for version 3.0 of ProMod3.
You may also want to read the documentation for:
1.3
2.0
2.1
3.1
3.2
Generate
|
Parameters: |
|
---|
GetDisulfidBridges
(all_pos, res_indices)¶Returns: | Pairs of indices (i,j), where res_indices[i] and res_indices[j] are assumed to have a disulfid bridge (CYS-SG pairs within 2.5 A). |
---|---|
Return type: | |
Parameters: |
|
SetupSystem
(all_pos, res_indices, loop_length, is_n_ter, is_c_ter, disulfid_bridges)¶SetupSystem
(all_pos, res_indices, loop_start_indices, loop_lengths, is_n_ter, is_c_ter, disulfid_bridges)Setup a new system for the loop / surrounding defined by all_pos and res_indices. Positions are taken from all_pos[res_indices[i]] and each of them must have all heavy atoms defined. Residue all_pos[i] is assumed to be peptide bound to all_pos[i-1] unless res. i is labelled as N-terminal. If all_pos[i-1] has an unset C (could happen for gaps while modelling), res. i will still be treated as if peptide-bound and the hydrogen attached to N is built assuming a helical structure as an approximation. Similarly, we do not generate OXT atoms for residues followed by a gap unless they are specifically labelled as C-terminal.
Each loop is defined by an entry in loop_start_indices and loop_lengths. For loop i_loop, res_indices[loop_start_indices[i_loop]] is the N-stem and res_indices[loop_start_indices[i_loop] + loop_lengths[i_loop] - 1] is the C-stem of the loop. The loop indices are expected to be contiguous in res_indices between the stems. The stems of the loop are kept fixed (N, CA, CB for N-stem / CA, CB, C, O for C-stem) unless they are terminal residues. Overlapping loops are merged and 0-length loops are removed.
If loop_length is given, a single loop with start index 0 and the given loop length is defined.
If a system was setup previously, it is overwritten completely.
If possible, this uses the “CPU” platform in OpenMM using the env. variable
PM3_OPENMM_CPU_THREADS
to define the number of desired threads (1 thread
is used if variable is undefined).
Parameters: |
|
---|---|
Raises: |
|
UpdatePositions
(all_pos, res_indices)¶Updates the positions in the system. Even though SetupSystem()
is
already very fast, this can speed up resetting a simulation. The data must
be consistent with the data used in the last SetupSystem()
call.
Parameters: |
|
---|---|
Raises: |
|
ExtractLoopPositions
(loop_pos)¶ExtractLoopPositions
(out_pos, res_indices)Extracts loop positions from the current simulation. If the simulation was
run (via GetSimulation()
), we internally have new positions for the
residues corresponding to all_pos[res_indices[i]] passed in
SetupSystem()
. This function then extracts these positions back into a
compatible output storage.
If loop_pos is passed, only the loops are stored. The loop residues are
stored contiguously and the loops are stored in the order given by
GetLoopStartIndices()
/ GetLoopLengths()
. The loop_pos storage
must have at least GetNumLoopResidues()
residues and must have
matching amino acid types with respect to the loop residues passed as
all_pos of SetupSystem()
.
If out_pos and res_indices are passed, residues must have matching amino
acid types for out_pos[res_indices[i]] and all_pos[res_indices[i]] of
SetupSystem()
. Only loop residues with i in the ranges defined by
GetLoopStartIndices()
/ GetLoopLengths()
are touched.
Parameters: |
|
---|---|
Raises: |
|
GetSimulation
()¶Returns: | Simulation object setup by SetupSystem() . Use this to run
MM simulations. |
---|---|
Return type: | ost.mol.mm.Simulation |
GetNumLoopResidues
()¶Returns: | Number of loop residues (incl. stems) of current simulation setup. |
---|---|
Return type: | int |
GetLoopStartIndices
()¶Returns: | Start indices of loops (see SetupSystem() ). |
---|---|
Return type: | list of int |
GetLoopLengths
()¶Returns: | Lengths of loops (see SetupSystem() ). |
---|---|
Return type: | list of int |
GetForcefieldAminoAcids
()¶Returns: | Forcefield-specific amino acid type for each residue of last
SetupSystem() call. |
---|---|
Return type: | list of ForcefieldAminoAcid |
GetIndexing
()¶The atoms of residue i are stored from idx[i] to idx[i+1]-1, where
idx is the list returned by this function. The atoms themselves are
ordered according to the indexing defined by ForcefieldLookup
.
Returns: | Indexing to positions vector used by the simulation object. The last item of the list contains the number of atoms in the system. |
---|---|
Return type: | list of int |
GetCpuPlatformSupport
()¶Returns: | True, if we will use OpenMM’s “CPU” platform (enabled by default if platform is available). False, if we use “Reference” platform. |
---|---|
Return type: | bool |
SetCpuPlatformSupport
(cpu_platform_support)¶Override “CPU” platform support setting. Useful to force the use of the OpenMM’s “Reference” platform for testing (by default we use “CPU” if it is available).
Parameters: | cpu_platform_support (bool ) – True, if “CPU” platform desired (ignored if
platform not available). False, if “Reference”
platform desired. |
---|
The ForcefieldLookup
class and its helpers define a fast way to extract
FF specific data for amino acids in a protein. We distinguish amino acid types
(and a few variants thereof) which may all be N- and/or C-terminal.
promod3.loop.
ForcefieldLookup
¶This class provides all functionality to generate ost.mol.mm.Simulation
objects. Specifically, we can:
GetNumAtoms()
(note that only OXT indexing depends on whether a
residue is terminal)ForcefieldConnectivity
), which include all
possible bonds, angles, dihedrals, impropers and LJ pairsThere is functionality to adapt the lookup and store it as needed or you can
load a predefined one with GetDefault()
or LoadCHARMM()
.
The atom indexing and GetAA()
are independent of the loaded file.
GetDefault
()¶Returns: | A default singleton instance (shared throughout the Python process) of this class with all data defined. Using this instance has the advantage that the object is only loaded once! |
---|---|
Return type: | ForcefieldLookup |
SetDefault
(new_default)¶Sets default singleton instance for all future GetDefault()
calls.
Parameters: | new_default (ForcefieldLookup ) – Lookup object to use as the new default. |
---|
LoadCHARMM
()¶Returns: | Predefined lookup object extracted from a CHARMM forcefield (loaded from disk) |
---|---|
Return type: | ForcefieldLookup |
Load
(filename)¶LoadPortable
(filename)¶Loads raw binary file generated with Save()
(optimized for fast
reading) / portable file generated with SavePortable()
(slower but
less machine-dependent).
Parameters: | filename (str ) – Path to the file from which to load. |
---|---|
Returns: | The loaded lookup object |
Return type: | ForcefieldLookup |
Raises: | RuntimeError if file cannot be opened or if
file cannot be parsed (see here for details). |
Save
(filename)¶SavePortable
(filename)¶Saves a raw / portable binary representation. Use portable files for distribution and convert locally to raw files. See here for details.
Parameters: | filename (str ) – Path to the file where it will be saved. |
---|---|
Raises: | RuntimeError if file cannot be opened. |
GetAA
(ff_aa)¶Returns: | Amino acid type for given ff_aa |
---|---|
Return type: | ost.conop.AminoAcid |
Parameters: | ff_aa (ForcefieldAminoAcid ) – Forcefield-specific amino acid type |
GetNumAtoms
(ff_aa, is_nter, is_cter)¶Returns: | Number of atoms for given input. |
---|---|
Return type: | |
Parameters: |
|
GetHeavyIndex
(ff_aa, atom_idx)¶GetHeavyIndex
(ff_aa, atom_name)Returns: | Internal index for given heavy atom in [0, |
---|---|
Return type: | |
Parameters: |
|
GetHydrogenIndex
(ff_aa, atom_idx)¶GetHydrogenIndex
(ff_aa, atom_name)Returns: | Internal index for given hydrogen atom in [0, |
---|---|
Return type: | |
Parameters: |
|
GetOXTIndex
(ff_aa, is_nter)¶Returns: | Internal index of OXT atom for C-terminal residue |
---|---|
Return type: | |
Parameters: |
|
GetFudgeLJ
()¶Returns: | Dampening factor for LJ 1,4 interactions (see
ost.mol.mm.Topology.SetFudgeLJ() ) |
---|---|
Return type: | float |
GetFudgeQQ
()¶Returns: | Dampening factor for electrostatic 1,4 interactions (see
ost.mol.mm.Topology.SetFudgeQQ() ) |
---|---|
Return type: | float |
GetMasses
(ff_aa, is_nter, is_cter)¶Returns: | Mass for each atom (see |
---|---|
Return type: |
|
Parameters: |
|
GetCharges
(ff_aa, is_nter, is_cter)¶Returns: | Charge for each atom (see |
---|---|
Return type: |
|
Parameters: |
|
GetSigmas
(ff_aa, is_nter, is_cter)¶Returns: | Sigma in nm for each atom
(see |
---|---|
Return type: |
|
Parameters: |
|
GetEpsilons
(ff_aa, is_nter, is_cter)¶Returns: | Epsilon in kJ/mol for each atom
(see |
---|---|
Return type: |
|
Parameters: |
|
GetInternalConnectivity
(ff_aa, is_nter, is_cter)¶Returns: | Internal connectivity of a residue |
---|---|
Return type: | |
Parameters: |
|
GetPeptideBoundConnectivity
(ff_aa_one, ff_aa_two, is_nter, is_cter)¶Returns: | All connectivity which include peptide bond between two residues
(additional to |
---|---|
Return type: | |
Parameters: |
|
GetDisulfidConnectivity
()¶Returns: | All connectivity which include disulfid bridge between two cysteins
(additional to GetInternalConnectivity() ) |
---|---|
Return type: | ForcefieldConnectivity |
SetFudgeLJ
(fudge)¶Set value for future GetFudgeLJ()
calls to fudge.
SetFudgeQQ
(fudge)¶Set value for future GetFudgeQQ()
calls to fudge.
SetMasses
(ff_aa, is_nter, is_cter, masses)¶Set value for future GetMasses()
calls to masses.
SetCharges
(ff_aa, is_nter, is_cter, charges)¶Set value for future GetCharges()
calls to charges.
SetSigmas
(ff_aa, is_nter, is_cter, sigmas)¶Set value for future GetSigmas()
calls to sigmas.
SetEpsilons
(ff_aa, is_nter, is_cter, epsilons)¶Set value for future GetEpsilons()
calls to epsilons.
SetInternalConnectivity
(ff_aa, is_nter, is_cter, connectivity)¶Set value for future GetInternalConnectivity()
calls to
connectivity.
SetPeptideBoundConnectivity
(ff_aa_one, ff_aa_two, is_nter, is_cter, connectivity)¶Set value for future GetPeptideBoundConnectivity()
calls to
connectivity.
SetDisulfidConnectivity
(connectivity)¶Set value for future GetDisulfidConnectivity()
calls to
connectivity.
promod3.loop.
ForcefieldAminoAcid
¶Enumerates the amino acid types for forcefields. The first 20 values
correspond to the 20 values of ost.conop.AminoAcid
. Additionally,
there are values for disulfid bridges (FF_CYS2), d-protonated histidine
(FF_HISD, default for ost.conop.HIS is FF_HISE) and FF_XXX for unknown
types. The full list of values is:
FF_ALA, FF_ARG, FF_ASN, FF_ASP, FF_GLN, FF_GLU, FF_LYS, FF_SER, FF_CYS, FF_MET, FF_TRP, FF_TYR, FF_THR, FF_VAL, FF_ILE, FF_LEU, FF_GLY, FF_PRO FF_HISE, FF_PHE, FF_CYS2, FF_HISD, FF_XXX
promod3.loop.
ForcefieldConnectivity
¶Contains lists of bonds, angles, dihedrals, impropers and LJ pairs (exclusions
are the combination of all bonds and 1,3 pairs of angles and are not stored
separately). Each type of connectivity has it’s own class (see below) storing
indices and parameters to be used for methods of ost.mol.mm.Topology
.
The indexing of atoms for internal connectivities is in [0, N-1], where N
= ForcefieldLookup.GetNumAtoms()
. For connectivities of pairs of
residues, atoms of the first residue are in [0, N1-1] and atoms of the
second one are in [N1, N1+N2-1], where N1 and N2 are the number of atoms
of the two residues. For disulfid bridges, N1 = N2 = GetNumAtoms(FF_CYS2,
False, False).
harmonic_bonds
¶List of harmonic bonds.
Type: | list of ForcefieldBondInfo |
---|
harmonic_angles
¶List of harmonic angles.
Type: | list of ForcefieldHarmonicAngleInfo |
---|
urey_bradley_angles
¶List of Urey-Bradley angles.
Type: | list of ForcefieldUreyBradleyAngleInfo |
---|
periodic_dihedrals
¶List of periodic dihedrals.
Type: | list of ForcefieldPeriodicDihedralInfo |
---|
periodic_impropers
¶List of periodic impropers.
Type: | list of ForcefieldPeriodicDihedralInfo |
---|
harmonic_impropers
¶List of harmonic impropers.
Type: | list of ForcefieldHarmonicImproperInfo |
---|
lj_pairs
¶List of LJ pairs.
Type: | list of ForcefieldLJPairInfo |
---|
promod3.loop.
ForcefieldBondInfo
¶Define harmonic bond (see ost.mol.mm.Topology.AddHarmonicBond()
)
promod3.loop.
ForcefieldHarmonicAngleInfo
¶Define harmonic angle (see ost.mol.mm.Topology.AddHarmonicAngle()
)
promod3.loop.
ForcefieldUreyBradleyAngleInfo
¶Define Urey-Bradley angle
(see ost.mol.mm.Topology.AddUreyBradleyAngle()
)
promod3.loop.
ForcefieldPeriodicDihedralInfo
¶Define periodic dihedral or improper (see
ost.mol.mm.Topology.AddPeriodicDihedral()
and
ost.mol.mm.Topology.AddPeriodicImproper()
)
promod3.loop.
ForcefieldHarmonicImproperInfo
¶Define harmonic improper (see ost.mol.mm.Topology.AddHarmonicImproper()
)
Enter search terms or a module, class or function name.
loop
- Loop Handling
ost.mol.mm
systems