PUFF: Pulsed Unconstrained Fluctating Forces

PUFF is a method for generating large forced motions. PUFF generates large conformational changes by defining target momentums (and rotational momentums) for groups of atoms. This drives conformational change if and only if the target momentum is sufficient to overcome internal energy barriers. The result is that critical forces can be derived from the trajectories.

The canonical test case for large forced motions is Atomic Force Microscopy (AFM) experiments that mechanically unfold the the 27th domain of titin. The AFM experiments directly measure the forces required to unfold the protein. PUFF has been used to unfold titin, and most importantly, PUFF can reproduce the magnitude of forces required to unfold the protein

As well, PUFF allows a straightforward analysis of intermediate states. PUFF generates trajectories of unfolding, but only for target momentum values above a certain value. Below this value, the system behaves intermittently, or not at all. Trajectories pulled at these intermittent velocities display interesting oscillations in the hydrogen bonding:

How PUFF works

The general schematic of PUFF is this:

pick a bunch of atoms
assign a target momentum (or a target rotational momentum)
divide the simulation into a bunch of short constant energy simulation (pulses)
at the beginning of each pulse, set the momentum of the atoms to the target values

This should be compared with the most used method in the literature, Steered Molecular Dynamics:

pick a bunch of atoms
define a path for a bunch of ghost points, usually defined in terms of a constant velocity, or constant rotational velocity
attach harmonic springs from the atoms to the ghost points
simulate in MD, letting the harmonic springs generate the forces to make the atoms follow the ghost points

SMD has produced wonderful trajectories where microscopic atomic-level details have been observed, such as in the unfolding of titin. However, the problem is that the forces are typically of an order of magnitude larger than that measured by the AFM experiments. This has been found irrespective of implicit or explicit solvent representations, and should be attributed to the use of harmonic constraints. In contrast, PUFF can generate better pulses. The downside with using PUFF is the overhead with running short pulses of MD, and then, subsequently, stringing them together.

Download and Installation

The PUFF package is a series of PYTHON scripts that runs AMBER and GROMACS (Oct, 2010). It has been tested on a MacOSX and LINUX environment:

DOWNLOAD

You can run the scripts directly in its directory, or put the scripts in the PATH of your shell resource file. It's also a good idea to add the location of the scripts to the PYTHONPATH environment variable so that you can import the library into your own PYTHON scripts.

Pre-requisites. The PUFF scripts relies on quite a few other packages. For basic functionality, you will need to install the following before anything will work:

AMBER a standard molecular-dynamics package
PYTHON the interpreter for the python scripts
NUMPY provides robust numerical arrays for analysis

AMBER costs unfortunately $400 for academics. PUFF requires AMBER on the command-line: make sure the AMBER programs are in the system PATH.

PYTHON is available on most systems, but the right version is required as the RIP scripts are written for PYTHON 2.5 (I had erroneously said 2.6, which is wrong! Wrong! WRONG!). They do not work with PYTHON 2.6 and PYTHON 3.0 as the external modules have not been translated to PYTHON 2.6 or higher.

For the PYTHON library NUMPY, the best way is to install the python package manager PIP, then at the command line:

sudo pip install numpy

Optional Packages. Although the flexibility analysis generates results as text files, there is an option that generates viewable images in HTML format from the text files. To do that, you will also need to install the following:

PYMOL for display of protein structures
MATPLOTLIB to generate 2D graphs
PIL the Python Imaging Library
MAKO a template library to make the HTML file

In particular the scripts, SHOWPDB.PY and PYTRAJ.PY, are wrappers for PYMOL to generate images and movies of the simulations.

Once PYMOL is installed, you must edit SHOWPDB.PY and PYTRAJ.PY by finding the variable 'pymol' in the scripts (should be at the top of the file), and setting 'pymol' to the full pathname of the PYMOL executable, which on my machine is:

pymol = "/Applications/MacPyMOL.app/Contents/MacOS/MacPyMOL"

For the PYTHON libraries MATPLOTLIB, PIL and MAKO, the best way is to install with PIP:

sudo pip install matplotlib pil mako

Before you run the scripts, don't forget to

source ~/.cshrc  [or]  source ~/.bash_profile

Example

An example is given of the Protein G beta-hairpin structure. The ends to be pulled are the first 2 residues at one end, and the last 2 residues at the other end.

In the 'example' directory, two files are given: 1) hairpin.pdb, a PDB file, where all hydrogen atoms have been stripped, and residue names renamed for AMBER, as well, a TER tag must be added to the C-terminus or any breaks in the chain; and 2) hairpin.puff-config config file, which should contain the following:

{
   'force_field'        : 'AMBER',
   'top_dir'            : 'hairpin'
   'in_pdb'             : 'hairpin.pdb',
   'n_equilibrate_step' : 1000,
   'temp_back'          : 300,
   'sim_dir'            : 'hairpin/v0.40_t50',
   'n_step'             : 50000,
   'n_step_per_pulse'   : 100,
   'pulse_strategy'     : 'PushApartByVel([0, 1], [14, 15], 0.400000)',
}

The parameters in the hairpin.puff-config file are:

'force_field' tells the library what MD package to use. Currently implemented: 'AMBER' and 'GROMACS'. However at the moment, 'GROMACS' requires the previous installation use of the GROMACS port of the AMBER force-field.
'top_dir' the directory where the minimization, prequilibration trajectories are stored
'in_pdb' the PDB of the structure you want to pull
'n_equilibrate_step' the number of simulation steps (in fs) of the pre-equilibration step
'temp_back' is the background temperature is used in the prequilibration step
'sim_dir' is the directory where the trajectory for the pulling is stored
'n_step' the number of steps (in fs) for the pulling
'n_step_per_pulse' the number of steps between the application of the force (also in fs)
'pulse_strategy' is the Python routine that implements the forces 'PushApartByVel([0, 1], [14, 15], 0.400000)'
1. the first parameter is the list of residues at one end of the pulling
2. the second parameter is the list of residues at the other end
3. the target velocity of the pulling

To run, go to the examples directory and run

python ../scripts/puff.py .

where the script PUFF.PY will attempt to run any file with the extension .puff-config.

The results are stored in the 'sim_dir', which in this case, is 'hairpin/v0.40_t50'. Apart from the trajectory itself, there is a file called 'md.acc' which records the parameters of the pulling during the simulation. Each row corresponds to the beginning of each pulse. There are six values for each row, corresponding to each pulse. They are

the change in velocity to be applied between the two ends to be pulled
the pre-pulse velocity between the two ends
the distance between the two ends
the target velocity for the two ends
the work done in the pulling
the combined mass of the two end groups in Daltons
a different calculation of the work

Alternatively, you can extract the parameters yourself from the 'md.trj' coordinate files, and the 'md.vel' for the velocities.

General Schema of the Library

PUFF is implemented as a PYTHON wrapper around the AMBER and GROMACS molecular dynamics package. Portions of the wrapper library also work with NAMD.

The PUFF scripts rests on PDBSTRUCT, a light-weight library that reads straight-forward PDB files and returns a PYTHON object that allows fine-grained manipulation of a protein structure. The PDBSTRUCT objects can be written back to PDB files or MD restart files.

The PUFF scripts provides a abstraction layer to the MD simulations, where in the future, other different packages may be plugged in. At the moment, only AMBER has been extensively tested, where AMBER was chosen because it supports implicit solvent force-fields. The scripts pre-processes PDB files (removing alternate conformations, adding terminii groups, defining disulfide bonds, etc.) so that they can be automatically fed into the AMBER topology preparation program TLEAP. The PUFF scripts then run short 100 fs constant-energy MD trajectories. Between each run, velocities are directly changed, which effectively implements RIP in 100 fs pulses.

More Info

Reference: "An Improved Strategy for Generating Forces in Steered Molecular Dynamics: The Mechanical Unfolding of Titin, e2lip3 and Ubiquitin" by Bosco K. Ho and David A. Agard, PLoS ONE (2010) 5(9): e13068.

This work was funded by HHMI in the laboratory of David Agard.

I have released the code under the GPL2 license. Everybody is free to do whatever you want with it, but if you release a modified version of it, it too must be free to use and modify.

If you have questions or bug reports, email: boscoh@gmail.com.