MM-PBSA

(Note: These tutorials are meant to provide illustrative examples of how to use the AMBER software suite to carry out simulations that can be run on a simple workstation in a reasonable period of time. They do not necessarily provide the optimal choice of parameters or methods for the particular application area.)
Copyright Ross Walker 2006

MM-PBSA - SECTION 2

MM-PBSA

By Ross Walker & Thomas Steinbrecher

2) Run the production simulation and obtain an ensemble of snapshots.

The production phase of the simulation should be run using the same conditions as the final phase of equilibration to prevent an abrupt jump in the potential energy due to a change in simulation conditions.

We will run a total of 2 ns or production recording the coordinates every 10 ps. This should be sufficiently far apart that the structures are uncorrelated. Depending on your system you might obtain good results with snapshots taken closer together. As long as all the structures you obtain are uncorrelated the more snapshots you have the lower the statistical error of your results should be. Note for system such as the RAS-RAF complex we have here a simulation time of 2 ns is most likely too short to obtain a set of uncorrelated snapshots that adequately sample the equilibrium ensemble. A value of 20 ns or so would probably be more appropriate. However, this will suffice for the purposes of this tutorial.

Here is the input file:

prod.in

prod ras-raf
 &cntrl
  imin=0,irest=1,ntx=5,
  nstlim=250000,dt=0.002,
  ntc=2,ntf=2,
  cut=8.0, ntb=2, ntp=1, taup=2.0,
  ntpr=5000, ntwx=5000,
  ntt=3, gamma_ln=2.0,
  temp0=300.0, ig=-1,
 /

This should then be run 4 times to obtain 2 ns of simulation time. Since this is a simple periodic boundary PME simulation one can use PMEMD to do the simulation if required. This will typically offer better performance and scaling in parallel. Below is an example script I used on San Diego Supercomputer Center's Teragrid Cluster to run this job on 96 processors. The calculation took a total of 10 hours.

run.x

#SDSC Teragrid PBS Script
#PBS -j oe
#PBS -l nodes=48:ppn=2
#PBS -l walltime=12:00:00
#PBS -q dque
#PBS -V
#PBS -M name@email.com
#PBS -A account_no
#PBS -N run_pmemd_96

cd /gpfs/projects/prod/
mpirun -v -machinefile $PBS_NODEFILE -np 96 /usr/local/apps/amber9/exe/pmemd -O -i prod.in -o prod1.out \
-p ras-raf_solvated.prmtop -c equil.rst -r prod1.rst -x prod1.mdcrd
mpirun -v -machinefile $PBS_NODEFILE -np 96 /usr/local/apps/amber9/exe/pmemd -O -i prod.in -o prod2.out \
-p ras-raf_solvated.prmtop -c prod1.rst -r prod2.rst -x prod2.mdcrd
mpirun -v -machinefile $PBS_NODEFILE -np 96 /usr/local/apps/amber9/exe/pmemd -O -i prod.in -o prod3.out \
-p ras-raf_solvated.prmtop -c prod2.rst -r prod3.rst -x prod3.mdcrd
mpirun -v -machinefile $PBS_NODEFILE -np 96 /usr/local/apps/amber9/exe/pmemd -O -i prod.in -o prod4.out \
-p ras-raf_solvated.prmtop -c prod3.rst -r prod4.rst -x prod4.mdcrd
gzip -9 prod*.mdcrd

Here are the output files: prod.tar.gz (84.8 MB)

It is essential, for good results, that our system still be exploring equilibrium phase space during the production phase. We will check this in the same fashion as we did for the last equilibration step by plotting the density, temperature, total energy and backbone RMSD.

DENSITY	TEMPERATURE
TOTAL ENERGY	BACKBONE RMSD

Note the production RMSD does not look to be truly equilibrated while the other properties are essentially constant (note the small scales). Ideally we should probably run a much longer production run (ca. 20 ns). For the purposes of this tutorial, however, we will continue with what we have.

We can now proceed to section 3 where we will calculate the binding free energy. The first link will take you to the instructions for using (and installing) the Python script MMPBSA.py. The second link will take you to the instructions for using the Perl script mm_pbsa.pl.

CLICK HERE TO GO TO SECTION 3 FOR MMPBSA.py

CLICK HERE TO GO TO SECTION 3 FOR mm_pbsa.pl