(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 2004
TUTORIAL 7 - SECTION 3
The Nudged Elastic Band Approach to Finding the Lowest
Energy
Pathway Between two States
By Ross Walker
3. Creating the NEB prmtop and inpcrd files
We now need to create a series of replicates (in this case 28 giving a total of 30 structures) that will stretch along our pathway and all be joined by springs. Ideally, so as to require the least amount of equilibration time, we should place these structures evenly along what we think is the best pathway. Lets for the moment though assume that we have no idea what the pathway looks like. In this situation what do we do? The simplest solution is to place a total of 15 structures on one end point and 15 on the other. We can then carry out simulated annealing on the system and the spring forces between the structures will cause then to stretch out along the pathway. Thus the NEB method should locate the pathway for us. This is the beauty of the simulated annealing approach to NEB implemented in AMBER 9. You don't necessarily need to provide an initial guess of the pathway. The situation where the distance between images (spring extension) is zero is explicitly handled in the code so that you do not obtain infinite spring forces on the first step.
In order to create our prmtop and inpcrd files we will use a tool that ships as part of the AMBER 9 installation: addles
This can do various manipulations of prmtop and inpcrd files. One such manipulation is to take a prmtop file for a single molecule and a series of inpcrd files and convert them into a LES (NEB) style prmtop file and inpcrd file representing all of the images. In our case we will use a single prmtop file for our molecule, since the topology is the same at both ends of the pathway and two inpcrd files representing the two pathway endpoints. These we will generate from the struct1.pdb and struct2.pdb files from the previous section.
The first stage is to create the prmtop file and the two inpcrd files. These will be regular prmtop and inpcrd files and we shall create them using tleap the command line version of Leap. This is the point at which we decide which force field/parameters etc that we will be using since the parameters in the single image prmtop file we create will be the ones used in the NEB simulation. For this tutorial we shall be using the FF03 force field. Make sure the AMBERHOME environment variable points to the amber9 installation and that $AMBERHOME/exe is in your path. Then fire up tleap:
tleap -s -f $AMBERHOME/dat/leap/cmd/leaprc.ff03
This will start up tleap with the ff03 force field. If we had any non-standard residues in our system now would be the time load any mol2 files, frcmod files etc. Next we enter the following commands into the tleap terminal to create the prmtop and inpcrd files:
str1=loadpdb struct1.pdb
str2=loadpdb struct2.pdb
saveamberparm str1 str1.prmtop str1.inpcrd
saveamberparm str2 str2.prmtop str2.inpcrd
quit
Input Files: struct1.pdb,
struct2.pdb
Output Files: str1.prmtop,
str1.inpcrd,
str2.prmtop, str2.inpcrd
The first thing that you should check here is that the two prmtop files are identical:
diff str1.prmtop str2.prmtop
This should return no differences with the exception of the line specifying the creation time. Since our two end points are the same molecule, just in different conformations, it is essential that the two prmtop files be identical. In fact I believe this is a requirement for NEB to work. So, since they are identical we can discard one of them. In this case I will keep str1.prmtop and use that from now on. At this point it is probably a good idea to make a quick visual check of the two end points using VMD. Load the str1.prmtop file and then load each of the two inpcrd files (selecting rst7 [vmd 1.8.3] or amber7 restart [vmd 1.8.4 or later] as the file type). You should then have a two frame view where you can flip between the two end points. You should check that this is what you expect. I always find such visual checks are always good as they can save you from wasting a lot of time (both human and computer) and effort trying to simulate something you didn't intend to. So, if everything looks good the next step is to use these files to an NEB compatible prmtop and incprd file.
For this we will use the addles program. You should refer to section 6.10.2 of the AMBER 9 manual for details on the options available for creating NEB prmtop and inpcrd files. We start by appending the str2.inpcrd file to the end of the str1.inpcrd file.
cat str2.inpcrd >> str1.inpcrd
If we were also planning on including some potential intermediate structures we would also append these to the str1.inpcrd file. This would be done in the order in which each of the structures would appear along the pathway. Next we need to create an input file for addles:
addles.in |
file rprm name=(str1.prmtop)
read file rcrd name=(str1.inpcrd) pack=2 read file wprm name=(neb.prmtop) wovr file wcrd name=(neb.inpcrd) wovr action ˜ use original mass omas pimd ˜ make 30 copies of atom 1 to 22 (the whole system) space numc=30 pick #prt 1 22 done *EOD |
The first line simply specifies the name of the prmtop file for a single image. The second line specifies the coordinate file and also the number of structures that are included in this coordinate file (pack=2). Addles will assign coordinates to images in appropriate sized blocks. For example, if the inpcrd file has 4 sets of coordinates and you have 20 images then image 1-5 will have coordinate set 1, image 6-10 have set 2, and so on. In this case we will have 30 images in total and so the first 15 will have the first structure from the coordinate file and the last 15 the second structure from the coordinate file.
The next two lines specify the name of the prmtop and inpcrd files we want created and that we want them to be written (wovr). We then have an action section which states what we want addles to do. Anything beginning with a tilda '~' is a comment so the use original mass line is simply a comment. The actual action is 'omas'. This tells addles that we want it to preserve the mass from the original prmtop file. This option is a function of the locally enhanced sampling method that addles was originally written for. For NEB simulations we always want it to preserve the mass from the prmtop file. The pimd flag should also always be present when creating NEB prmtop files. It tells addles not to generate an exclusion list that crosses all of the different images. If we didn't specify this option it would not cause problems but it could lead to us having a very very large prmtop file. The next line is a comment. The penultimate line states that that we want 30 copies to be created (numc=30) and that we want atoms 1 to 22 to be replicated, in this case all of the atoms in our system. For the moment when running NEB simulations you MUST replicate ALL atoms as the NEB algorithm at present can only act on the entire system. The final line acts as a marker for the end of the input file. So now we can run addles to generate the prmtop and inpcrd files we need.
$AMBERHOME/exe/addles <addles.in >addles.out
Output Files: addles.out, neb.prmtop, neb.inpcrd
The addles.out file is not required beyond this point but you should take a look at it to see if there are any errors reported.
Now we have the NEB prmtop and inpcrd files prepared we are ready to proceed to the next stage where we will run NEB MD and slowly heat the system up.