Bio.PDB.ResidueDepth module¶
Calculation of residue depth using command line tool MSMS.
This module uses Michel Sanner’s MSMS program for the surface calculation. See: http://mgltools.scripps.edu/packages/MSMS
Residue depth is the average distance of the atoms of a residue from the solvent accessible surface.
Residue Depth:
from Bio.PDB.ResidueDepth import ResidueDepth
from Bio.PDB.PDBParser import PDBParser
parser = PDBParser()
structure = parser.get_structure("1a8o", "Tests/PDB/1A8O.pdb")
model = structure[0]
rd = ResidueDepth(model)
print(rd['A',(' ', 152, ' ')])
Direct MSMS interface, typical use:
from Bio.PDB.ResidueDepth import get_surface
surface = get_surface(model)
The surface is a Numeric array with all the surface vertices.
Distance to surface:
from Bio.PDB.ResidueDepth import min_dist
coord = (1.113, 35.393, 9.268)
dist = min_dist(coord, surface)
where coord is the coord of an atom within the volume bound by the surface (ie. atom depth).
To calculate the residue depth (average atom depth of the atoms in a residue):
from Bio.PDB.ResidueDepth import residue_depth
chain = model['A']
res152 = chain[152]
rd = residue_depth(res152, surface)
-
Bio.PDB.ResidueDepth.
get_surface
(model, PDB_TO_XYZR=None, MSMS='msms')¶ Represent molecular surface as a vertex list array.
Return a Numpy array that represents the vertex list of the molecular surface.
- Arguments:
PDB_TO_XYZR - deprecated, ignore this.
MSMS - msms executable (used as argument to subprocess.call)
-
Bio.PDB.ResidueDepth.
min_dist
(coord, surface)¶ Return minimum distance between coord and surface.
-
Bio.PDB.ResidueDepth.
residue_depth
(residue, surface)¶ Residue depth as average depth of all its atoms.
Return average distance to surface for all atoms in a residue, ie. the residue depth.
-
Bio.PDB.ResidueDepth.
ca_depth
(residue, surface)¶ Return CA depth.
-
class
Bio.PDB.ResidueDepth.
ResidueDepth
(model, pdb_file=None)¶ Bases:
Bio.PDB.AbstractPropertyMap.AbstractPropertyMap
Calculate residue and CA depth for all residues.
-
__init__
(self, model, pdb_file=None)¶ Initialize the class.
-