Analysis plugins

Analysis plugins consist of a book-keeping class derived from gromacs.analysis.core.Plugin and a “worker” class (a child of gromacs.analysis.core.Worker), which contains the actual analysis code.

Plugins

Plugin modules are named like the plugin class but the filename is all lower case. All plugin classes are available in the gromacs.analysis.plugins name space.

analysis.plugins – Plugin Modules

Classes for gromacs.analysis.core.Simulation that provide code to analyze trajectory data.

New analysis plugins should follow the API sketched out in gromacs.analysis.core; see an example for use there.

List of plugins

Right now the number of plugins is limited. Feel free to contribute your own by sending it to the package author. You will be acknowledged in the list below.

Plugins for analysis.
plugin author description
CysAccessibility [1] estimate accessibility of Cys residues by water
HelixBundle [1] g_bundle analysis of helices
Distances [1] time series of distances
MinDistances [1] time series of shortest distances
COM [1] time series of centres of mass
Dihedrals [1] analysis of dihedral angles
RMSF [1] calculate root mean square fluctuations
RMSD [1] calculate root mean square distance
Energy [1] terms from the energy file
Plugins for trajectory manipulation and status queries.
plugin author description
Trajectories [1] write xy-fitted trajectories
StripWater [1] remove solvent (and optionally fit to reference)
Ls [1] simple ls (for testing)

Footnotes

[1](1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12) Oliver Beckstein <oliver.beckstein@bioch.ox.ac.uk>

Plugin classes

class gromacs.analysis.plugins.CysAccessibility(name=None, simulation=None, **kwargs)

CysAccessibility plugin.

For each frame of a trajectory, the shortest distance of all water oxygens to all cysteine sulphur atoms is computed. For computational efficiency, only distances smaller than a cutoff are taken into account. A histogram of the distances shows how close water molecules can get to cysteines. The closest approach distance should be indicative of the reactivity of the SH group with crosslinking agents.

class CysAccessibility(cysteines[, cys_cutoff[, name[, simulation]]])
Arguments:
name : string

plugin name (used to access it)

simulation : instance

The gromacs.analysis.Simulation instance that owns the plugin.

cysteines : list

list of all resids (eg from the sequence) that are used as labels or in the form ‘Cys<resid>’. (required)

cys_cutoff : number

cutoff in nm for the minimum S-OW distance [1.0]

Note that all Cys residues in the protein are analyzed. Therefore, the list of cysteine labels must contain as many entries as there are cysteines in the protein. There are no sanity checks.

Registers the plugin with the simulation class.

Specific keyword arguments are listed below, all other kwargs are passed through.

Arguments:
name : string

Name of the plugin. Should differ for different instances. Defaults to the class name.

simulation : Simulation instance

The Simulation instance that owns this plugin instance. Can be None but then the register() method has to be called manually with a simulation instance later.

kwargs

All other keyword arguments are passed to the Worker.

CysAccessibility.worker_class
alias of _CysAccessibility
class gromacs.analysis.plugins.HelixBundle(name=None, simulation=None, **kwargs)

HelixBundle plugin.

gromacs.g_bundle() helix analysis

class HelixBundle([helixtable, offset, with_kinks[, name[, simulation]]])
Arguments:
helixtable

reST table with columns “name”, “top”, “bottom”, “kink”; see gromacs.analysis.plugins.helixbundle for details

offset

add the offset to the residue numbers in helixtable [0]

helixndx

provide a index file with the appropriate groups instead of the table; also requires na

na

number of helices

with_kinks

take kinks into account [True]

name

plugin name [HelixBundle]

simulation

The gromacs.analysis.Simulation instance that owns the plugin. [None]

Registers the plugin with the simulation class.

Specific keyword arguments are listed below, all other kwargs are passed through.

Arguments:
name : string

Name of the plugin. Should differ for different instances. Defaults to the class name.

simulation : Simulation instance

The Simulation instance that owns this plugin instance. Can be None but then the register() method has to be called manually with a simulation instance later.

kwargs

All other keyword arguments are passed to the Worker.

HelixBundle.worker_class
alias of _HelixBundle
class gromacs.analysis.plugins.Distances(name=None, simulation=None, **kwargs)

Distances plugin.

The distance between the center of mass of two index groups are calculated for each time step and written to files.

class Distances(groups, ndx[, cutoff[, name[, simulation]]])
Arguments:
name : string

plugin name (used to access it)

simulation : instance

The gromacs.analysis.Simulation instance that owns the plugin.

groups : list of index group names

The first entry is the primary group, the second is the *secondary group.

ndx : index filename or list

All index files that contain the listed groups.

cutoff : float

A contact is recorded if the distance is <cutoff [0.6 nm]

Example:

Generate index files with the groups of interest, for instance with gromacs.cbook.IndexBuilder:

from gromacs.cbook import IndexBuilder
A_grp, A_ndx = IndexBuilder(tpr, ['@a 62549 & r NA'], names=['Na1_ion'], offset=-9, 
                            out_ndx='Na1.ndx', name_all="Na1").combine()
B = IndexBuilder(tpr, ['S312:OG','T313:OG1','A38:O','I41:O','A309:O'], offset=-9, 
                      out_ndx='Na1_site.ndx', name_all="Na1_site")
B_grp, B_ndx = B.combine()                            
all_ndx_files = [A_ndx, B_ndx]

To calculate the distance between “Na1” and the “Na1_site”, create an instance with the appropriate parameters and add them to a gromacs.analysis.Simulation instance:

dist_Na1_site = Distances(name='Dsite', groups=['Na1', 'Na1_site'], ndx=all_ndx_files)
S.add_plugin(dist_Na1_site)

Registers the plugin with the simulation class.

Specific keyword arguments are listed below, all other kwargs are passed through.

Arguments:
name : string

Name of the plugin. Should differ for different instances. Defaults to the class name.

simulation : Simulation instance

The Simulation instance that owns this plugin instance. Can be None but then the register() method has to be called manually with a simulation instance later.

kwargs

All other keyword arguments are passed to the Worker.

Distances.worker_class
alias of _Distances
class gromacs.analysis.plugins.MinDistances(name=None, simulation=None, **kwargs)

MinDistances plugin.

The minimum distances between the members of at least two index groups and the number of contacts are calculated for each time step and written to files.

class Distances(groups, ndx[, cutoff[, name[, simulation]]])
Arguments:
name : string

plugin name (used to access it)

simulation : instance

The gromacs.analysis.Simulation instance that owns the plugin.

groups : list of index group names

The first entry is the primary group. All other entries are secondary groups and the plugin calculates the minimum distance between members of the primary group and the members of each secondary group.

ndx : index filename or list

All index files that contain the listed groups.

cutoff : float

A contact is recorded if the distance is <cutoff [0.6 nm]

Example:

Generate index files with the groups of interest, for instance with gromacs.cbook.IndexBuilder:

from gromacs.cbook import IndexBuilder
A_grp, A_ndx = IndexBuilder(tpr, ['@a 62549 & r NA'], names=['Na1_ion'], offset=-9, 
                            out_ndx='Na1.ndx', name_all="Na1").combine()
B = IndexBuilder(tpr, ['S312:OG','T313:OG1','A38:O','I41:O','A309:O'], offset=-9, 
                      out_ndx='Na1_site.ndx', name_all="Na1_site")
B_grp, B_ndx = B.combine()                            
all_ndx_files = [A_ndx, B_ndx]

To calculate the distance between “Na1” and the “Na1_site”, create an instance with the appropriate parameters and add them to a gromacs.analysis.Simulation instance:

dist_Na1_site = Distances(name='Dsite', groups=['Na1', 'Na1_site'], ndx=all_ndx_files)
S.add_plugin(dist_Na1_site)

To calculate the individual distances:

dist_Na1_res = Distances(name='Dres', groups=['Na1']+B.names, ndx=all_ndx_files)
S.add_plugin(dist_Na1_res)

(Keeping the second IndexBuilder instance B allows us to directly use all groups without typing them, B.names = ['A309_O', 'S312_OG', 'I41_O', 'T313_OG1', 'A38_O'].)

Registers the plugin with the simulation class.

Specific keyword arguments are listed below, all other kwargs are passed through.

Arguments:
name : string

Name of the plugin. Should differ for different instances. Defaults to the class name.

simulation : Simulation instance

The Simulation instance that owns this plugin instance. Can be None but then the register() method has to be called manually with a simulation instance later.

kwargs

All other keyword arguments are passed to the Worker.

MinDistances.worker_class
alias of _MinDistances
class gromacs.analysis.plugins.COM(name=None, simulation=None, **kwargs)

COM plugin.

Calculate the centre of mass (COM) of various index groups.

class COM(group_names[, ndx[, offset[, name[, simulation]]]])

Registers the plugin with the simulation class.

Specific keyword arguments are listed below, all other kwargs are passed through.

Arguments:
name : string

Name of the plugin. Should differ for different instances. Defaults to the class name.

simulation : Simulation instance

The Simulation instance that owns this plugin instance. Can be None but then the register() method has to be called manually with a simulation instance later.

kwargs

All other keyword arguments are passed to the Worker.

COM.worker_class
alias of _COM
class gromacs.analysis.plugins.Dihedrals(name=None, simulation=None, **kwargs)

Dihedrals plugin.

class Dihedrals(dihedrals[, labels[, name[, simulation]]])
Keywords:
dihedrals

list of tuples; each tuple contains atom indices that define the dihedral.

labels

optional list of labels for the dihedrals. Must have as many entries as dihedrals.

name : string

plugin name (used to access it)

simulation : instance

The gromacs.analysis.Simulation instance that owns the plugin.

Registers the plugin with the simulation class.

Specific keyword arguments are listed below, all other kwargs are passed through.

Arguments:
name : string

Name of the plugin. Should differ for different instances. Defaults to the class name.

simulation : Simulation instance

The Simulation instance that owns this plugin instance. Can be None but then the register() method has to be called manually with a simulation instance later.

kwargs

All other keyword arguments are passed to the Worker.

Dihedrals.worker_class
alias of _Dihedrals
class gromacs.analysis.plugins.RMSF(name=None, simulation=None, **kwargs)

RMSF plugin.

Compute the root mean square fluctuations (RMSF) of the C-alpha atoms. The trajectory is always fitted to the reference structure in the tpr file.

class RMSF([name[, simulation]])
Arguments:
name : string

plugin name (used to access it)

simulation : instance

The gromacs.analysis.Simulation instance that owns the plugin.

Registers the plugin with the simulation class.

Specific keyword arguments are listed below, all other kwargs are passed through.

Arguments:
name : string

Name of the plugin. Should differ for different instances. Defaults to the class name.

simulation : Simulation instance

The Simulation instance that owns this plugin instance. Can be None but then the register() method has to be called manually with a simulation instance later.

kwargs

All other keyword arguments are passed to the Worker.

RMSF.worker_class
alias of _RMSF
class gromacs.analysis.plugins.RMSD(name=None, simulation=None, **kwargs)

RMSD plugin.

Calculation of the root mean square distance (RMSD) of a protein structure over the course of a MD simulation.

The trajectory is always fitted to the reference structure in the tpr file.

class RMSD([name[, simulation]])
Arguments:
name : string

plugin name (used to access it)

simulation : instance

The gromacs.analysis.Simulation instance that owns the plugin.

Registers the plugin with the simulation class.

Specific keyword arguments are listed below, all other kwargs are passed through.

Arguments:
name : string

Name of the plugin. Should differ for different instances. Defaults to the class name.

simulation : Simulation instance

The Simulation instance that owns this plugin instance. Can be None but then the register() method has to be called manually with a simulation instance later.

kwargs

All other keyword arguments are passed to the Worker.

RMSD.worker_class
alias of _RMSD
class gromacs.analysis.plugins.Energy(name=None, simulation=None, **kwargs)

Energy plugin.

Analysis of terms in the Gromacs energy (edr) file.

class Energy([name[, simulation]])
Arguments:
name : string

plugin name (used to access it)

simulation : instance

The gromacs.analysis.Simulation instance that owns the plugin.

Registers the plugin with the simulation class.

Specific keyword arguments are listed below, all other kwargs are passed through.

Arguments:
name : string

Name of the plugin. Should differ for different instances. Defaults to the class name.

simulation : Simulation instance

The Simulation instance that owns this plugin instance. Can be None but then the register() method has to be called manually with a simulation instance later.

kwargs

All other keyword arguments are passed to the Worker.

Energy.worker_class
alias of _Energy
class gromacs.analysis.plugins.Trajectories(name=None, simulation=None, **kwargs)

Trajectories plugin.

Write new xy-fitted trajectories (see gromacs.cbook.trj_fitandcenter()),

class Trajectories([name[, simulation]])
Arguments:
name : string

plugin name (used to access it)

simulation : instance

The gromacs.analysis.Simulation instance that owns the plugin.

Registers the plugin with the simulation class.

Specific keyword arguments are listed below, all other kwargs are passed through.

Arguments:
name : string

Name of the plugin. Should differ for different instances. Defaults to the class name.

simulation : Simulation instance

The Simulation instance that owns this plugin instance. Can be None but then the register() method has to be called manually with a simulation instance later.

kwargs

All other keyword arguments are passed to the Worker.

Trajectories.worker_class
alias of _Trajectories
class gromacs.analysis.plugins.StripWater(name=None, simulation=None, **kwargs)

StripWater plugin.

Write a new trajectory which has the water index group removed.

class StripWater([selection[, name[, simulation]]])
Arguments:
selection

optional selection for the water instead of “SOL”

name : string

plugin name (used to access it)

simulation : instance

The gromacs.analysis.Simulation instance that owns the plugin.

Registers the plugin with the simulation class.

Specific keyword arguments are listed below, all other kwargs are passed through.

Arguments:
name : string

Name of the plugin. Should differ for different instances. Defaults to the class name.

simulation : Simulation instance

The Simulation instance that owns this plugin instance. Can be None but then the register() method has to be called manually with a simulation instance later.

kwargs

All other keyword arguments are passed to the Worker.

StripWater.worker_class
alias of _StripWater
class gromacs.analysis.plugins.Ls(name=None, simulation=None, **kwargs)

ls plugin.

This simply lists the files on disk. It is useful for testing the plugin architecture.

class Ls([name[, simulation]])
Arguments:
name : string

plugin name (used to access it)

simulation : instance

The gromacs.analysis.Simulation instance that owns the plugin.

Registers the plugin with the simulation class.

Specific keyword arguments are listed below, all other kwargs are passed through.

Arguments:
name : string

Name of the plugin. Should differ for different instances. Defaults to the class name.

simulation : Simulation instance

The Simulation instance that owns this plugin instance. Can be None but then the register() method has to be called manually with a simulation instance later.

kwargs

All other keyword arguments are passed to the Worker.

Ls.worker_class
alias of _Ls

Developer notes

In principle all that needs to be done to automatically load plugins is to add their name to __plugins__. See the source code for further comments and how the auto loading of plugins is done.

gromacs.analysis.plugins.__plugins__
All available plugin names are listed here. Because this is used to automatically set up imports a module file must be named like the plugin class it contains but in all lower case. For example, the Distances plugin class is contained in the module distances (the file plugins/distances.py).
gromacs.analysis.plugins.__plugin_classes__
Gives access to all available plugin classes (or use the module __dict__)

Helper modules

Various helper classes that can be used by plugins. Because they are not necessarily restricted to a single plugin they have been moved into separate modules for code reuse.

analysis.plugins.dist — Helper Class for g_dist

dist contains helper classes for other analysis plugins that want to make use of the Gromacs command g_dist.

Overview

The task we are solving is to analyze output from

g_dist -f md.xtc -s md.tpr -n cys_ow.ndx -dist 1.0 | bzip2 -vc > mindist_C60_OW_1nm.dat.bz2 

and produce a histogram of minimum contact distances. This should provide an estimate for water accessibility of the atom (here: SG of Cys60).

File format

g_dist with the -dist CUTOFF option writes to stdout the identity of all atoms within the cutoff distance and the distance itself:

Selected 22: 'CYSH_CYSH_60_&_SG'
Selected 25: 'OW'
....
t: 184  6682 SOL 35993 OW  0.955138 (nm)
t: 184  10028 SOL 46031 OW  0.803889 (nm)
t: 185  6682 SOL 35993 OW  0.879949 (nm)
t: 185  10028 SOL 46031 OW  0.738299 (nm)
t: 186  6682 SOL 35993 OW  0.897016 (nm)
t: 186  10028 SOL 46031 OW  0.788268 (nm)
t: 187  6682 SOL 35993 OW  0.997688 (nm)
....

Classes

class gromacs.analysis.plugins.dist.Mindist(datasource, cutoff=None)

The Mindist class allows analysis of the output from g_dist -dist CUTOFF.

Output is read from a file or stream. The raw data is transformed into a true ‘mindist’ time series (available in the Mindist.distances attribute): for each frame only the shortest distance is stored (whereas g_dist provides all distances below the cutoff).

Todo:
  • Save analysis to pickle or data files.
  • Export data as simple data files for plotting in other programs.

Note

gromacs.tools.G_mindist is apparently providing exactly the service that is required: a timeseries of the minimum distance between two groups. Feel free to use that tool instead :-).

Read mindist data from file or stream.

Arguments:
datasource

a filename (plain, gzip, bzip2) or file object

cutoff

the -dist CUTOFF that was provided to g_dist; if supplied we work around a bug in g_dist (at least in Gromacs 4.0.2) in which sometimes numbers >> CUTOFF are printed.

histogram(nbins=None, lo=None, hi=None, midpoints=False, normed=True)

Returns a distribution or histogram of the minimum distances:

If no values for the bin edges are given then they are set to 0.1 below and 0.1 above the minimum and maximum values seen in the data.

If the number of bins is not provided then it is set so that on average 100 counts come to a bin. Set nbins manually if the histogram only contains a single bin (and then get more data)!

Keywords:
nbins : int

number of bins

lo : float

lower edge of histogram

hi : float

upper edge of histogram

midpoints : boolean

False: return edges. True: return midpoints

normed : boolean

True: return probability distribution. False: histogram

hist
Histogram of the minimum distances.
dist
Distribution of the minimum distances.
edges
Edges of the histogram of the minimum distances.
midpoints
Midpoints of the histogram of the minimum distances.
plot(**kwargs)

Plot histograms with matplotlib’s plot() function:

plot(**histogramargs, **plotargs)

Arguments for both Mindist.histogram() and pylab.plot() can be provided (qv).

class gromacs.analysis.plugins.dist.GdistData(stream)

Object that represents the standard output of g_dist -dist CUTOFF.

Initialize from a stream (e.g. a pipe) and the iterate over the instance to get the data line by line. Each line consists of a tuple

(frame, distance)

Initialize with an open stream to the data (eg stdin or file).

Arguments:
stream

open stream (file or pipe or really any iterator providing the data from g_dist); the stream is not closed automatically when the iterator completes.

__iter__()
Iterator that filters the input stream and returns (frame, distance) tuples.

gromacs.analysis.plugins.gridmatmd — Lipid bilayer analysis helper

This helper module contains code to drive GridMAT-MD.pl, available from the GridMAT-MD home page and written by WJ Allen et al [Allen2009] . The GromacsWrapper distribution comes with version 1.0.2 of GridMAT-MD.pl and includes a small patch so that it can accept filenames on the command line.

References

[Allen2009]W. J. Allen, J. A. Lemkul, and D. R. Bevan. (2009) “GridMAT-MD: A Grid-based Membrane Analysis Tool for Use With Molecular Dynamics.” J. Comput. Chem. 30 (12): 1952-1958.

Module contents

class gromacs.analysis.plugins.gridmatmd.GridMatMD(config, filenames)

Analysis of lipid bilayers with GridMAT-MD.

It requires a configuration file and a list of structure files (gro or pdb) as input. See the documentation (pdf) for the format of the config file. Note that the bilayer keyword will be ignored in the config file.

Set up GridMAT-MD analysis.

Arguments:
config : filename

input file for GridMAT-MD (see docs)

filenames : list or glob-pattern

list of gro or pdb files, or a glob pattern that creates such a list

imshow(name, **kwargs)
Display array name with pylab.imshow.
run()
Run analysis on all files and average results.
run_frame(frame)

Run GridMAT-MD on a single frame and return results.

Arguments:frame is a filename (gro or pdb)
Returns:a dict of GridMapData objects; the keys are “top”, “bottom”, “average”
save(name)
Save object as pickle.
class gromacs.analysis.plugins.gridmatmd.GridMatData(filename, shape=None, delta=None)

Represent GridMatMD data file.

The loaded array data is accessible as a numpy array in GridMatData.array and bins and midpoints as GridMatData.bins and GridMatData.midpoints respectively.

Load the data into a numpy array.

The filename is an output file from GridMAT-MD. shape and delta are optional. The shape of the array is parsed from the filename if not provided. The spacing is set to (1,1) if not provided.

Arguments:
filename

2D grid as written by GridMAT-MD

shape

Shape tuple (NX, NY) of the array in filename.

delta

Tuple of bin sizes of grid (DX, DY).

imshow(**kwargs)
Display data as a 2D image using pylab.imshow().
parse_filename(filename)
Get dimensions from filename
class gromacs.analysis.plugins.gridmatmd.Grid2D(data, bins)

Represents a 2D array with bin sizes.

Addition and subtraction of grids is defined for the arrays and the bins. Multiplication and division with scalars is also defined. Each operation returns a new Grid2D object.

(Actually, it should work for arrays of any dimension, not just 2D.)

Initialize the Grid2D instance.

Arguments:
data

array data, e.g. a list of array

bins

tuple of lists of bin edges, one for each dimension

imshow(**kwargs)
Display data as a 2D image using pylab.imshow().
__add__(other)
Add arrays and bins (really only makes sense when averaging).
__sub__(other)
Subtract other from self (also subtracts bins... which is odd but consistent).
__mul__(x)
Multiply arrays (and bins) by a scalar x.
__div__(x)
Divide arrays (and bins) by a scalar x.