Author: | Lukas Stelzl, Oliver Beckstein |
---|---|
Year: | 2011-2012 |
Copyright: | GNU Public License v2 |
With the help of this module, HOLE can be run on frames in a trajectory. Data can be combined and analyzed. HOLE [Smart1993] [Smart1996] must be installed separately.
References
[Smart1993] | O.S. Smart, J.M. Goodfellow and B.A. Wallace. The Pore Dimensions of Gramicidin A. Biophysical Journal 65:2455-2460, 1993. |
[Smart1996] | O.S. Smart, J.G. Neduvelil, X. Wang, B.A. Wallace, and M.S.P. Sansom. HOLE: A program for the analysis of the pore dimensions of ion channel structural models. J.Mol.Graph., 14:354–360, 1996. URL http://hole.biop.ox.ac.uk/hole. |
Gramicidin A (gA) channel:
from MDAnalysis.analysis.hole import HOLE, HOLEtraj
from MDAnalysis.tests.datafiles import PDB_HOLE
H = HOLE(PDB_HOLE, executable="~/hole2/bin/hole") # set path to your hole binary
H.run()
H.collect()
H.plot(linewidth=3, color="black", label=False)
Analyzing a trajectory:
u = MDAnalysis.Universe(psf, trajectory)
H = HOLEtraj(u, ...)
H.run()
H.plot3D()
The profiles are available as the attribute HOLEtraj.profiles (H.profiles in the example) and are indexed by frame number but can also be indexed by an arbitrary order parameter as shown in the next example.
In order to classify the HOLE profiles the RMSD to a reference structure is calculated for each trajectory frame (e.g. using the MDAnalysis.analysis.rms.RMSD analysis class). Then the HOLE profiles can be ordered by the RMSD, which acts as an order parameter.
import MDAnalysis.analysis.hole
import MDAnalysis.analysis.rms
MDAnalysis.start_logging()
ref = MDAnalysis.Universe("ref.psf", "ref.pdb") # reference structure
u = MDAnalysis.Universe("ref.psf", "traj.xtc") # trajectory
# calculate RMSD
R = MDAnalysis.analysis.rms.RMSD(u, reference=ref, select="protein", mass_weighted=True)
R.run()
# HOLE analysis with order parameters
H = MDAnalysis.analysis.hole.HOLEtraj(u, cvect=[0,0,1], orderparameters=R.rmsd[:,2])
H.run()
The HOLEtraj.profiles dictionary will have the order parameter as key for each frame. The plot functions will automatically sort the profiles by ascending order parameter. To access the individual profiles one can simply iterate over the sorted profiles (see HOLEtraj.sorted_profiles_iter())
for q, profile in H:
print "orderparameter = %g" % q
print "min(R) = %g" % profile.radius.min()
A profile is stored as a numpy.recarray:
frame rxncoord radius
The HOLE.profiles or HOLEtraj.profiles dictionary holds one profile for each key. By default the keys are the frame numbers but HOLEtraj can take the optional orderparameters keyword argument and load an arbitrary order parameter for each frame. In that case, the key becomes the orderparameter.
Note
The profiles dict is not ordered and hence one typically needs to manually order the keys first. Furthermore, duplicate keys are not possible: In the case of duplicate orderparameters, the last one read will be stored in the dict.
Run HOLE on a single frame or a DCD trajectory.
Only a subset of all HOLE control parameters is supported and can be set with keyword arguments.
HOLE (as a FORTRAN77 program) has a number of limitations when it comes to filename lengths (must be shorter than the empirically found HOLE.HOLE_MAX_LENGTH). This class tries to work around them by creating temporary symlinks to files when needed but this can still fail when permissions are not correctly set on the current directory.
Running HOLE with the HOLE class is a 3-step process:
- set up the class with all desired parameters
- run HOLE with HOLE.run()
- collect the data from the output file with HOLE.collect()
The class also provides some simple plotting functions of the collected data such as HOLE.plot() or HOLE.plot3D().
New in version 0.7.7.
Set up parameters to run HOLE on PDB filename.
Arguments: | filename
|
---|---|
Keywords: |
logfile
sphpdb
step
cpoint
cvect
shorto
ignore_residues
sphpdb
executable
|
After running HOLE.collect(), this dict contains all the HOLE profiles, indexed by the frame number. If only a single frame was analyzed then this will be HOLE.profiles[0]. Note that the order is random; one needs to sort the keys first.
Maximum number of characters in a filename (limitation of HOLE)
Return filename suitable for HOLE.
HOLE is limited to filenames <= HOLE.HOLE_MAX_LENGTH. This method
returns filename if HOLE can process it
returns a relative path (see os.path.relpath()) if that shortens the path sufficiently
creates a symlink to filename (os.symlink()) in a safe temporary directory and returns the path of the symlink. The temporary directory and the symlink are stored in HOLE.tempfiles and HOLE.tempdirs and deleted when the HOLE instance is deleted or garbage collected.
By default the temporary directory is created inside the current directory in order to keep that path name short. This can be changed with the tmpdir keyword (e.g. one can use “/tmp”).
Parse the output from a HOLE run into numpy recarrays.
Can deal with outputs containing multiple frames. Output format:
frame rxncoord radius
The method saves the result as HOLE.profiles, a dictionary indexed by the frame number. Each entry is a numpy.recarray.
If the keyword outdir is supplied (e.g. ”.”) then each profile is saved to a gzipped data file.
Keywords: |
|
---|
Process HOLE output to create a smooth pore surface suitable for VMD.
Takes the sphpdb file and feeds it to sph_process and sos_triangle as described under Visualization of HOLE results.
Load the output file filename into VMD by issuing in the tcl console
source hole.vmd
The level of detail is determined by HOLE.dotden (which can be overriden by keyword dotden).
The surface will be colored so that parts that are inaccessible to water (pore radius < 1.15 Å) are red, water accessible parts (1.15 Å > pore radius < 2.30 Å) are green and wide areas (pore radius > 2.30 Å are blue).
List of residues that are ignore by default. Can be changed with the ignore_residues keyword.
Return the minimum radius over all profiles as a function of q
Plot profiles in a 1D graph.
Lines are coloured according to the colour map cmap.
Keywords: |
|
---|
All other kwargs are passed to matplotlib.pyplot.plot().
Stacked graph of profiles.
Lines are coloured according to the colour map cmap.
Keywords: |
|
---|
Based on Stack Overflow post 3d plots using matplotlib.
Save profiles as a Python pickle file filename.
Load profiles dictionary with
import cPickle
profiles = cPickle.load(open(filename))
Return an iterator over profiles sorted by frame/order parameter q.
The iterator produces tuples (q, profile).
Analyze all frames in a trajectory.
The HOLE class provides a direct interface to HOLE. HOLE itself has limited support for analysing trajectories but cannot deal with all the trajectory formats understood by MDAnalysis. This class can take any universe and feed it to HOLE. By default it sequentially creates a PDB for each frame and runs HOLE on the frame.
Set up the class.
Arguments: |
|
---|
After running HOLE.collect(), this dict contains all the HOLE profiles, indexed by the frame number or the order parameter (if orderparameters was supplied). Note that the order is random; one needs to sort the keys first.
Guess a point inside the pore.
This method simply uses the center of geometry of the protein selection as a guess. selection is “protein” by default.
Return the minimum radius over all profiles as a function of q
Plot profiles in a 1D graph.
Lines are coloured according to the colour map cmap.
Keywords: |
|
---|
All other kwargs are passed to matplotlib.pyplot.plot().
Stacked graph of profiles.
Lines are coloured according to the colour map cmap.
Keywords: |
|
---|
Based on Stack Overflow post 3d plots using matplotlib.
Run HOLE on the whole trajectory and collect profiles.
Keyword arguments start, stop, and step can be used to only analyse part of the trajectory.
Save profiles as a Python pickle file filename.
Load profiles dictionary with
import cPickle
profiles = cPickle.load(open(filename))
Return an iterator over profiles sorted by frame/order parameter q.
The iterator produces tuples (q, profile).
Built-in HOLE radii (based on simple.rad from the HOLE distribution): van der Waals radii are AMBER united atom from Weiner et al. (1984), JACS, vol 106 pp765-768. Simple - Only use one value for each element C O H etc. Added radii for K+, NA+, CL- (Pauling hydration radius from Hille 2002). The data file can be written with the convenience function write_simplerad2().
Write the built-in radii in SIMPLE2_RAD to filename.
Does nothing if filename already exists.
Return sequence v as a string of numbers with spaces as separators.
In the special case of None, the empty string “” is returned.
Raised when an external application failed.
The error code is specific for the application.
New in version 0.7.7.
x.__init__(...) initializes x; see help(type(x)) for signature