4.1.3. Using qhull to define regions for hopping analysis — hop.qhull
¶
Interface to some functions of the `qhull`_ (or rather the qconvex) program. `qhull`_ must be installed separately (see links).
The main functionality is to define a region in space within the convex hull of
a protein. The hull is typically defined by a selection of atoms and written as
a “density” file for use in hop
.
4.1.3.1. Example¶
In this example the convex hull of the C-alpha atoms is computed. Initially, points must be extracted from the structure to a file:
hop.qhull.points_from_selection(psf='protein.psf', pdb='protein.pdb', filename='ca_100%.dat')
and saved to file ca_100%.dat
.
This is usually too large and also entails regions of the hydration
shell outside of interal cavities. A relatively robust workaround for
roughly globular proteins is to shrink the convex hull, using the
scale
argument of hop.qhull.make_ca_points()
. Shrinking to
70% appears to be a good starting point:
hop.qhull.points_from_selection(psf='protein.psf', pdb='protein.pdb', filename='ca_70%.dat', scale=0.7)
The convex hull itself is generated from the datafile of the points:
Q70 = hop.qhull.ConvexHull('ca_70%.dat', workdir='cavity70%')
Another density grid b
(such as a real water density for the bulk) is
currently required to generate a pseudo density based on the convex
hull. The real density provides the grid on which the convex hull is
mapped:
b = hop.sitemap.Density(filename='bulk')
QD70 = Q70.Density(b)
(This maps out sites at the threshold level set in b
; change it
with the hop.sitemap.Density.map_sites()
method if required.)
- Insert a bulk density::
- QD70.site_insert_bulk(b)
-
class
hop.qhull.
ConvexHull
(coordinates, workdir=None, prefix=None)¶ The convex hull of a set of points.
The convex hull is calculated with the `qhull`_ program.
Compute convex hull and populate data structures.
Arguments: - coordinates: input suitable for qconvex
- workdir: store intermediate files in workdir (tmp dir by default)
- prefix: filename prefix for intermediate output files
-
Density
(density, fillvalue=None)¶ Create a Density object of the interior of the convex hall.
Uses another Density object density as a template for the grid.
Note
This is rather slow and should be optimized.
-
point_inside
(point)¶ Check if point [x,y,z] is inside the polyhedron defined by planes.
Iff for all i: plane[i]([x,y,z]) = n*[x,y,z] + p < 0 <==> [x,y,z] inside
(i.e. [x,y,z] is under all planes and the planes completely define the enclosed space
-
points_inside
(points)¶ Return bool array for all points:
True: inside False: outside
Arguments: - points = [[x1,y1,z1], …] or an iterator that supplies points
- planes: normal forms of planes
Returns: Array with truth values such as [True, False, True, …]
-
read_planes
()¶ Read planes from qconvex n file.
Numpy array [[n1,n2,n3,-p], …] for planes n*x = -p.
Planes are oriented and point outwards.
-
read_vertices
()¶ Read vertices from qconvex p file.
Numpy array of points [[x,y,z], …]
-
wd
(*args)¶ Return path in workdir.
-
class
hop.qhull.
VertexPDBWriter
(filename)¶ PDB writer that implements a subset of the PDB 3.2 standard. http://www.wwpdb.org/documentation/format32/v3.2.html
-
ATOM
(serial=None, name=None, altLoc=None, resName=None, chainID=None, resSeq=None, iCode=None, x=None, y=None, z=None, occupancy=1.0, tempFactor=0.0, element=None, charge=0)¶ Write ATOM record. http://www.wwpdb.org/documentation/format32/sect9.html Only some keword args are optional (altLoc, iCode, chainID), for some defaults are set.
All inputs are cut to the maximum allowed length. For integer numbers the highest-value digits are chopped (so that the serial and reSeq wrap); for strings the trailing characters are chopped.
Note: Floats are not checked and can potentially screw up the format.
-
REMARK
(*remark)¶ Write generic REMARK record (without number). http://www.wwpdb.org/documentation/format32/remarks1.html http://www.wwpdb.org/documentation/format32/remarks2.html
-
TITLE
(*title)¶ Write TITLE record. http://www.wwpdb.org/documentation/format32/sect2.html
-
write
(coordinates, name='CA', resname='VRT', resid=1)¶ Write coordinates as CA.
-
-
hop.qhull.
points_from_selection
(*args, **kwargs)¶ Create a list of points from selected atoms in a format suitable for
qhull
.points_from_selection(topology, structure, selection=”name CA”, filename=”points.dat”, scale=None)
Arguments: - psf: Charmm topology file
- pdb: coordinates
- selection: MDAnalysis select_atoms() selection string [C-alpha atoms]
- filename: name of the output file; used as input for
ConvexHull
- scale: scale points around the centre of geometry; values of 0.5 - 0.7 typically ensure that the convex hull is inside the protein; default is to not to scale, i.e. scale = 1.
-
hop.qhull.
write_coordinates
(filename, points, scale=None)¶ Write an array of points to a file suitable for qhull.