Biskit :: Dock :: Complex :: Complex :: Class Complex
[hide private]
[frames] | no frames]

Class Complex

source code

Manage receptor and ligand structure, store additional information, calculate some properties.

Instance Methods [hide private]
  __init__(self, rec_model=None, lig_model=None, ligMatrix=None, info={})
str version(self)
Version of Dock.Complex
  __getstate__(self)
Called before pickling.
dict __getitem__(self, key)
Returns C.info[ key ]
  __setitem__(self, key, v)
  __delitem__(self, key)
  __contains__(self, key)
  __setstate__(self, state)
called for unpickling the object.
  __defaults(self)
backwards compatibility to earlier pickled models
[str] keys(self)
Keys of info dictionary.
1|0 has_key(self, k)
Check if info dictionary has key.
[any] values(self, keys=[], default=None)
Use:
any OR None get(self, key, default=None)
Use:
PCRModel rec(self)
Return PCRModel object.
PCRModel lig(self, force=0, cache=1)
Return ligand structure transformed.
PCRModel model(self)
Model with both receptor and ligand.
array, vector ligMatrix(self)
Return transformation matrix to transform original Ligand PCRModel into docked conformation.
  setLigMatrix(self, m)
Set a ligand translation/rotation matrix.
dict getInfo(self)
Return contents of the info dictionary.
  writeLigand(self, fname, ter=0)
Write transformed ligand to pdb file.
  writeReceptor(self, fname, ter=0)
Write receptor to pdb without TER statement (unless removeTer!=0).
  writeComplex(self, fname, ter=0)
Write single pdb containing both receptor and ligand separated by END but not TER (unless ter!=0).
  slim(self)
Remove coordinates and atoms of ligand and receptor from memory, if they can be restored from file, compress contact matrix.
float contactsOverlap(self, ref, cutoff=None)
Fraction of overlapping residue-residue contacts between this and reference complex.
int contactsShared(self, reference, cutoff=None)
Number of equal residue-residue contacts in this and reference complex.
int contactsDiff(self, ref, cutoff=None)
Number of different residue-residue contacts in this and reference complex.
float fractionNativeContacts(self, ref, cutoff=None)
Fraction of native residue-residue contacts.
(float, float) fractionNativeSurface(self, cont, contRef)
fraction of atoms/residues that are involved in any contacts in both complexes.
float rmsLig(self, ref)
Rms of ligand from reference ligand after superimposing the two receptors.
float rmsInterface(self, ref, cutoff=4.5, fit=1)
Rmsd between this and reference interface.
[(str.str)..] contactResPairs(self, cm=None)
Get list of residue names paired up in contacts.
dict contactResDistribution(self, cm=None)
Count occurrence of residues in protein-protein interface.
dict contactTypeDistribution(self, cm=None)
Count occurrence of residue types aromatic (a), charged(c), polar(p), else/unpolar(u).
dict resPairCounts(self, cm=None)
Count how often (all possible) res-res combinations occur in the given contacts.
dict OR None loadResContacts(self)
Uncompress residue contact matrix if necessary.
array resContacts(self, cutoff=4.5, maskRec=None, maskLig=None, refComplex=None, force=0, cache=1, cache_pw=0)
Matrix of all residue - residue contacts between receptor and ligand.
2D array __alignMatrixDimension(self, cm, thisSeq, castSeq, axis=0)
Correct one dimension of contactMatrix by inserting and deleting columns, so that it can be later compared to contact matrices based on slightly different sequences.
array __unmaskedMatrix(self, contacts, rec_mask, lig_mask)
Map contacts between selected rec and lig atoms back to all atoms matrix.
array __atomContacts(self, cutoff, rec_mask, lig_mask, cache)
Intermolecular distances below cutoff after applying the two masks.
array atomContacts(self, cutoff=4.5, rec_mask=None, lig_mask=None, cache=0, map_back=1)
Find all inter-molecular atom-atom contacts between rec and lig
array __resContacts(self, cutoff, maskRec=None, maskLig=None, cache=0)
Find all inter-molecule residue-residue contacts between receptor and ligand.
array __atom2residueMatrix(self, m)
Reduce binary matrix of n x k atoms to binary matrix of i x j residues.
[1|0], [1|0], [1|0], [1|0] equalAtoms(self, ref)
Compare two complexes' atom content of receptor and ligand.
[int], [int], [int], [int] compareAtoms(self, ref)
Compare atom content and oder of rec and lig of 2 complexes.
Complex take(self, rec_pos, lig_pos)
Get copy of this complex with given atoms of rec and lig.
Complex compress(self, rec_mask, lig_mask)
Compress complex using a rec and lig mask.
float contPairScore(self, cutoff=6.0)
Score interaction surface residue pairs.
(float, float) or dict interfaceArea(self, profiles=0, log=StdLog(), verbose=1)
Calculate the difference between the surface area of the complex vs.
(array, array, array) prosa2003Energy(self)
Calculate Prosa2003 total energies for the complex.
dict foldXEnergy(self, force=1)
Calculate E_rec and/or E_lig only if not given:
float conservationScore(self, cons_type='cons_ent', ranNr=150, log=StdLog(), verbose=1)
Score of conserved residue pairs in the interaction surface.
array rtTuple2matrix(self, r, t)
Put rotation and translation matrix into single 4x4 matrix.
array extractLigandMatrix(self, lig)
Find transformation matrix for rigid body-transformed ligand.
array, array __findTransformation(self, x, y)
Match two arrays by rotation and translation.
array __pairwiseDistances(self, u, v)
pairwise distance between 2 3-D numpy arrays of atom coordinates.
(array, array) __extractLigandMatrix(self, fcomplex)
Compare structure from hex complex with original ligand pdb and store transformation matrix of ligand in self.ligandMatrix.

Method Details [hide private]

__init__(self, rec_model=None, lig_model=None, ligMatrix=None, info={})
(Constructor)

source code 
Parameters:
  • rec_model (PDBModel OR PCRModel) - model of original receptor conformation
  • lig_model (PDBModel OR PCRModel) - model of original ligand conformation
  • ligMatrix (matrix) - Numeric array 4 by 4, ligand transformation matrix
  • info (dict) - {'eshape':-123.3, 'rms':12.2 }

version(self)

source code 

Version of Dock.Complex
Returns: str
version of class

__getstate__(self)

source code 

Called before pickling.

__getitem__(self, key)
(Indexing operator)

source code 
Returns: dict
C.info[ key ]

__setitem__(self, key, v)
(Index assignment operator)

source code 

__delitem__(self, key)
(Index deletion operator)

source code 

__contains__(self, key)
(In operator)

source code 

__setstate__(self, state)

source code 

called for unpickling the object.

__defaults(self)

source code 

backwards compatibility to earlier pickled models

keys(self)

source code 

Keys of info dictionary.
Returns: [str]
list of keys

has_key(self, k)

source code 

Check if info dictionary has key.
Parameters:
  • k (str) - key to test
Returns: 1|0
dict has key

values(self, keys=[], default=None)

source code 

Use:
  values( keys=None, default=None ) -> list of info dict values
Parameters:
  • keys ([str]) - give only values of these keys
  • default (any) - only used with keys!=[], default value for missing keys
Returns: [any]
all values (default) or values of requested keys (ordered)

get(self, key, default=None)

source code 

Use:
  C.get(k[,default]) -> C.info[k] if C.info.has_key(k), else
                        default defaults to None.
Parameters:
  • key (str) - key to call
  • default (any) - default value if dictionary doesn't have key (default: None)
Returns: any OR None
dictionary value of key OR else default

rec(self)

source code 

Return PCRModel object.
Returns: PCRModel
receptor model

lig(self, force=0, cache=1)

source code 

Return ligand structure transformed. Use cached object, if available.
Parameters:
  • force (1|0) - force new transformation (default: 0)
  • cache (1|0) - cache transformed model (default: 1)
Returns: PCRModel
ligand model

model(self)

source code 

Model with both receptor and ligand.
Returns: PCRModel
single PDBModel with first rec then lig

ligMatrix(self)

source code 

Return transformation matrix to transform original Ligand PCRModel into docked conformation.
Returns: array, vector
tuple with rotation and transformation matrix

setLigMatrix(self, m)

source code 

Set a ligand translation/rotation matrix.
Parameters:
  • m (array) - translation array 4x4 with rotation

getInfo(self)

source code 

Return contents of the info dictionary.
Returns: dict
info dictionary

writeLigand(self, fname, ter=0)

source code 

Write transformed ligand to pdb file. Remove TER record for xplor usage (unless removeTer!=0).
Parameters:
  • fname (str) - output filename
  • ter (0|1|2) - option for how to treat the TER statement:
                 - 0, don't write any TER statements (default)
                 - 1, restore original TER statements
                 - 2, put TER between all detected chains
    

writeReceptor(self, fname, ter=0)

source code 

Write receptor to pdb without TER statement (unless removeTer!=0).
Parameters:
  • fname (str) - output file name
  • ter (0|1|2) - option for how to treat the TER statement:
                 - 0, don't write any TER statements (default)
                 - 1, restore original TER statements
                 - 2, put TER between all detected chains
    

writeComplex(self, fname, ter=0)

source code 

Write single pdb containing both receptor and ligand separated by END but not TER (unless ter!=0).
Parameters:
  • fname (str) - filename of pdb
  • ter - option for how to treat the TER statement:
                 - 0, don't write any TER statements (default)
                 - 1, restore original TER statements
                 - 2, put TER between all detected chains
    

slim(self)

source code 

Remove coordinates and atoms of ligand and receptor from memory, if they can be restored from file, compress contact matrix.

Note: CALLED BEFORE PICKLING

contactsOverlap(self, ref, cutoff=None)

source code 

Fraction of overlapping residue-residue contacts between this and reference complex.
Parameters:
  • ref (Complex) - reference complex
  • cutoff (float) - maximal atom-atom distance, None .. previous setting
Returns: float
fraction of contacts shared between this and ref (normalized to number of all contacts)

contactsShared(self, reference, cutoff=None)

source code 

Number of equal residue-residue contacts in this and reference complex.
Parameters:
  • reference (Complex) - reference complex
  • cutoff (float) - cutoff for atom-atom contact to be counted
Returns: int
the number or residue-residue contacts that are common to both this and reference:
 abs( N.sum( N.sum( contactMatrix_a - contactMatrix_b )))

contactsDiff(self, ref, cutoff=None)

source code 

Number of different residue-residue contacts in this and reference complex.
Parameters:
  • ref (Complex) - to compare this one with
  • cutoff (float) - maximal atom-atom distance, None .. previous setting
Returns: int
number of contacts different in this and refererence complex.

fractionNativeContacts(self, ref, cutoff=None)

source code 

Fraction of native residue-residue contacts.
Parameters:
  • ref (Complex) - native complex
  • cutoff (float) - maximal atom-atom distance, None .. previous setting
Returns: float
fraction of native contacts

fractionNativeSurface(self, cont, contRef)

source code 

fraction of atoms/residues that are involved in any contacts in both complexes.
Parameters:
  • cont (matrix) - contact matrix
  • contRef (matrix) - reference contact matrix
Returns: (float, float)
(fractRec, fractLig), fraction of atoms/residues that are involved in any contacts in both complexes

rmsLig(self, ref)

source code 

Rms of ligand from reference ligand after superimposing the two receptors.
Parameters:
  • ref (Complex) - reference complex, must have identical atoms
Returns: float
ligand rmsd

Note: not implemented

rmsInterface(self, ref, cutoff=4.5, fit=1)

source code 

Rmsd between this and reference interface. The interface is defined as any residue that has an atom which is within the distance given by |cutoff| from its partner.
Parameters:
  • ref (Complex) - reference complex
  • cutoff (float) - atom distance cutoff for interface residue definition (default: 4.5)
  • fit (1|0) - least-squares fit before calculating the rms (default: 1)
Returns: float
interface rmad

contactResPairs(self, cm=None)

source code 

Get list of residue names paired up in contacts.
Parameters:
  • cm (matrix) - pre-calculated contact matrix (default: None)
Returns: [(str.str)..]
list of tuples [('N','G'), ('P','C')..]

contactResDistribution(self, cm=None)

source code 

Count occurrence of residues in protein-protein interface.
Parameters:
  • cm (matrix) - pre-calculated contact matrix (default: None)
Returns: dict
{'A':3, 'C':1, .. } (20 standard amino acids)

contactTypeDistribution(self, cm=None)

source code 

Count occurrence of residue types aromatic (a), charged(c), polar(p), else/unpolar(u).
Parameters:
  • cm (matrix) - pre-calculated contact matrix (default: None)
Returns: dict
{'a':4, 'c':5, 'p':12, 'u':4}

resPairCounts(self, cm=None)

source code 

Count how often (all possible) res-res combinations occur in the given contacts.
Parameters:
  • cm (matrix) - pre-calculated contact matrix (default: None)
Returns: dict
{'AA':3,'AC':0, .. 'WW':2, 'WY':0 }

loadResContacts(self)

source code 

Uncompress residue contact matrix if necessary.
Returns: dict OR None
dict with contact matrix and parameters OR None

resContacts(self, cutoff=4.5, maskRec=None, maskLig=None, refComplex=None, force=0, cache=1, cache_pw=0)

source code 

Matrix of all residue - residue contacts between receptor and ligand. Result is cached.
Parameters:
  • cutoff (float) - float/int, cutoff in \AA for atom-atom contact to be counted ( default 4.5; if None, last one used or 4.5)
  • maskRec ([1|0]) - atom mask, receptor atoms to consider
  • maskLig ([1|0]) - atom mask, ligand atoms to consider
  • force (0|1) - re-calculate even if cached matrix is available (default: 0)
  • cache (0|1) - cache result (default: 1)
  • cache_pw (0|1) - cache pairwise atom distances (default:0)
  • refComplex (Complex) - if <> None, 'reshape' contact matrix, so that residue positions match residue positions in matrices from castComplex. Useful if calling complex is a crystallized reference structure with slight sequence variations from the docked receptor and ligand.
Returns: array
residue contact matrix, 2-D array(residues_receptor x residues_ligand) of 0 or 1

__alignMatrixDimension(self, cm, thisSeq, castSeq, axis=0)

source code 

Correct one dimension of contactMatrix by inserting and deleting columns, so that it can be later compared to contact matrices based on slightly different sequences.
Parameters:
  • cm (array) - contact matrix, 2D matrix of residue contacts recceptor x ligand sequence
  • thisSeq (string) - AA sequence of this dimension of the contactMatrix
  • castSeq (string) - AA sequence of this dimension in the other contact
  • axis (1|0) - which dimension to adapt (0=receptor, 1=ligand)
Returns: 2D array
contact matrix with residue contacts compatible to refSeq.

__unmaskedMatrix(self, contacts, rec_mask, lig_mask)

source code 

Map contacts between selected rec and lig atoms back to all atoms matrix.
Parameters:
  • contacts (array) - contact matrix, array sum_rec_mask x sum_lig_mask
  • rec_mask ([1|0]) - atom mask
  • lig_mask ([1|0]) - atom mask
Returns: array
atom contact matrix, array N_atoms_rec x N_atoms_lig

__atomContacts(self, cutoff, rec_mask, lig_mask, cache)

source code 

Intermolecular distances below cutoff after applying the two masks.
Parameters:
  • cutoff (float) - cutoff for atom-atom contact in \AA
  • rec_mask ([1|0]) - atom mask
  • lig_mask ([1|0]) - atom mask
  • cache (1|0) - cache pairwise atom distance matrix
Returns: array
atom contact matrix, array sum_rec_mask x sum_lig_mask

atomContacts(self, cutoff=4.5, rec_mask=None, lig_mask=None, cache=0, map_back=1)

source code 

Find all inter-molecular atom-atom contacts between rec and lig
Parameters:
  • cutoff (float) - cutoff for atom - atom contact in \AA
  • rec_mask ([1|0]) - atom mask (default: all heavy)
  • lig_mask ([1|0]) - atom mask (default: all heavy)
  • cache (1|0) - cache pairwise atom distance matrix (default: 0)
  • map_back (1|0) - map masked matrix back to matrix for all atoms
Returns: array
atom contact matrix, Numpy array N(atoms_lig) x N(atoms_rec)

__resContacts(self, cutoff, maskRec=None, maskLig=None, cache=0)

source code 

Find all inter-molecule residue-residue contacts between receptor and ligand. A contact between A and B is set if any heavy atom of A is within |cutoff| A of B.
Parameters:
  • cutoff (float) - distance cutoff in \AA
  • maskRec ([1|0]) - atom mask (default: all heavy)
  • maskLig ([1|0]) - atom mask (default: all heavy)
  • cache (1|0) - cache pairwise atom distance matrix to pw_dist (default:0)
Returns: array
residue contact matrix, 2-D Numpy N.array(residues_receptor x residues_ligand) where 1-contact, 0-no contact

__atom2residueMatrix(self, m)

source code 

Reduce binary matrix of n x k atoms to binary matrix of i x j residues.
Parameters:
  • m (array) - atom contact matrix, array n x k with 1(contact) or 0(no contact)
Returns: array
residue contact matrix, 2-D numpy array(residues_receptor x residues_ligand)

equalAtoms(self, ref)

source code 

Compare two complexes' atom content of receptor and ligand. Apply to SORTED models without HETATOMS. Coordinates are not checked.
Parameters:
  • ref (Complex) - reference complex
Returns: [1|0], [1|0], [1|0], [1|0]
(mask_rec, mask_lig, mask_rec_ref, mask_lig_ref), 4 atom masks for all equal atoms

Note: in some rare cases m1.equalAtoms( m2 ) gives a different result than m2.equalAtoms( m1 ). This is due to the used SequenceMatcher class.

compareAtoms(self, ref)

source code 

Compare atom content and oder of rec and lig of 2 complexes. Get list of all atom positions of rec and lig in both complexes that have a matching atom (by residue and atom name) in the other complex. Indices for this complex are returned in the right oder to match the atom order of ref. I.e., in contrast to equalAtoms, castAtoms can equalize two complexes where atoms are ordered differently within the residues.
Parameters:
  • ref (Complex) - reference complex
Returns: [int], [int], [int], [int]
(ind_rec, ind_lig, ind_rec_ref, ind_lig_ref), 4 lists of indices

take(self, rec_pos, lig_pos)

source code 

Get copy of this complex with given atoms of rec and lig.
Parameters:
  • rec_pos ([int]) - receptor indices to take
  • lig_pos ([int]) - ligand indices to take
Returns: Complex
new complex

compress(self, rec_mask, lig_mask)

source code 

Compress complex using a rec and lig mask.
Parameters:
  • rec_mask ([1|0]) - atom mask
  • lig_mask ([1|0]) - atom mask
Returns: Complex
compressed complex

contPairScore(self, cutoff=6.0)

source code 

Score interaction surface residue pairs. Info on Scoring matrix see Biskit.molUtils
Parameters:
  • cutoff (float) - CB-CB distance cutoff for defining a contact (default: 6.0)
Returns: float
score

interfaceArea(self, profiles=0, log=StdLog(), verbose=1)

source code 

Calculate the difference between the surface area of the complex vs. its free components in square angstrom.
Parameters:
  • profiles (1|0 (default: 0)) - option to return the lig, rec and com profiles rather than the value (both absolute and relative values are returned)
  • log (Biskit.LogFile) - log file [STDOUT]
  • verbose (bool | int) - give progress report [1]
Returns: (float, float) or dict
AS area, MS area OR a dictionary of lig, rec and com profiles

prosa2003Energy(self)

source code 

Calculate Prosa2003 total energies for the complex.
Returns: (array, array, array)
Prosa energy profiles, N.array(pair energy, surface energy, combined energy )

foldXEnergy(self, force=1)

source code 

Calculate E_rec and/or E_lig only if not given:
 E_com - (E_rec + E_lig).
Parameters:
  • force (1|0) - force calc of E_rec and E_lig
Returns: dict
dict with the different fold-X energy terms

conservationScore(self, cons_type='cons_ent', ranNr=150, log=StdLog(), verbose=1)

source code 

Score of conserved residue pairs in the interaction surface. Optionally, normalized by radom surface contacts.
Parameters:
  • cons_type (str) - precalculated conservation profile name, see Biskit.PDBDope.
  • ranNr (int) - number of random matricies to use (default: 150)
  • log (Biskit.LogFile) - log file [STDOUT]
  • verbose (bool | int) - give progress report [1]
Returns: float
conservation score

rtTuple2matrix(self, r, t)

source code 

Put rotation and translation matrix into single 4x4 matrix.
Parameters:
  • r (array) - rotation matric, array 3x3 of float
  • t (vector) - translation vector, array 1x3 of float
Returns: array
rotation/translation matrix, array 4x4 of float

extractLigandMatrix(self, lig)

source code 

Find transformation matrix for rigid body-transformed ligand.
Parameters:
  • lig (PDBModel) - ligand model
Returns: array
rotation/translation matrix, array 4x4 of float

__findTransformation(self, x, y)

source code 

Match two arrays by rotation and translation. Returns the rotation matrix and the translation vector. Back transformation: for atom i new coordinates will be:
   y_new[i] = N.dot(r, y[i]) + t
for all atoms in one step:
   y_new = N.dot(y, N.transpose(r)) + t
Parameters:
  • x (array) - coordinates
  • y (array) - coordinates
Returns: array, array
rotation matrix, translation vector

Author: Michael Habeck

__pairwiseDistances(self, u, v)

source code 

pairwise distance between 2 3-D numpy arrays of atom coordinates.
Parameters:
  • u (array) - coordinates
  • v (array) - coordinates
Returns: array
Numpy array len(u) x len(v)

Author: Wolfgang Rieping.

__extractLigandMatrix(self, fcomplex)

source code 

Compare structure from hex complex with original ligand pdb and store transformation matrix of ligand in self.ligandMatrix.
Parameters:
  • fcomplex (complec) - pdb file with hex complex
Returns: (array, array)
rotation matrix and translation matrix as tuple