Amber Workshop - tutorial 3

(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

We have our original FBP1_noH_min.pdb and FBP2_noH_min.pdb files from above - these suffice for the protein reference state without the extra proton, but we now need to calculate the energy needed to add a proton to this structure. In order to create a structure that has a protonated version of glutamic acid we need to use the residue name GLH.

In xleap if you want to see the currently available units you can use the list command. E.g

Here we can see a residue called GLH which is the protonated version of GLU. Also defined is ASP and ASH. Let's quickly take a look at the difference between GLU and GLH. We can do this by editing the two units:

We should have two model windows up showing the two units. You should be able to see that GLH has a protonated hydroxyl oxygen.

If we select Unit -> Calculate Net Charge we will also see that the de-protonated GLU residue has a charge of -1.0 while the protonated residue has a charge of zero.

So, by changing the residue name in the pdb we can control how many hydrogens xleap will add. So let's create prmtop and inpcrd files for our minimised FBP1 and FBP2 with residues 7, 10 and 15 protonated in turn.

First off we will protonate residue 7. Let's create a copy of our minimised pdb files:

ATOM     69  C   SER     6      -4.542  -8.449  -5.931
ATOM     70  O   SER     6      -4.564  -7.264  -5.599
ATOM     71  N   GLU     7      -3.473  -8.991  -6.525
ATOM     72  H   GLU     7      -3.580  -9.974  -6.761
ATOM     73  CA  GLU     7      -2.188  -8.352  -6.892
ATOM     74  HA  GLU     7      -1.575  -9.134  -7.340
ATOM     75  CB  GLU     7      -2.380  -7.273  -7.977
ATOM     76 2HB  GLU     7      -3.083  -6.518  -7.625
ATOM     77 3HB  GLU     7      -1.418  -6.790  -8.149
ATOM     78  CG  GLU     7      -2.867  -7.854  -9.314
ATOM     79 2HG  GLU     7      -2.122  -8.568  -9.673
ATOM     80 3HG  GLU     7      -3.798  -8.400  -9.142
ATOM     81  CD  GLU     7      -3.106  -6.785 -10.403
ATOM     82  OE1 GLU     7      -2.607  -5.639 -10.295
ATOM     83  OE2 GLU     7      -3.798  -7.101 -11.403
ATOM     84  C   GLU     7      -1.321  -7.810  -5.738
ATOM     85  O   GLU     7      -0.100  -7.730  -5.893
ATOM     86  N   TRP     8      -1.892  -7.504  -4.573
ATOM     87  H   TRP     8      -2.902  -7.549  -4.528

We change all GLU's for residue 7 here to GLH.

cp FBP1_noH_min.pdb FBP1_10H_min.pdb
cp FBP2_noH_min.pdb FBP2_10H_min.pdb
cp FBP1_noH_min.pdb FBP1_15H_min.pdb
cp FBP2_noH_min.pdb FBP2_15H_min.pdb

Edit each pdb in turn and replace the relevant residue with GLH or ASH respetively.

$AMBERHOME/exe/xleap -s -f $AMBERHOME/dat/leap/cmd/leaprc.ff99
FBP1_7H=loadpdb FBP1_7H_min.pdb
FBP1_10H=loadpdb FBP1_10H_min.pdb
FBP1_15H=loadpdb FBP1_15H_min.pdb
FBP2_7H=loadpdb FBP2_7H_min.pdb
FBP2_10H=loadpdb FBP2_10H_min.pdb
FBP2_15H=loadpdb FBP2_15H_min.pdb

In each case leap should, as we expect, report that it added a single missing hydrogen. Now save the prmtop and inpcrd files:

saveamberparm FBP1_7H FBP1_7H.prmtop FBP1_7H.inpcrd
saveamberparm FBP1_10H FBP1_10H.prmtop FBP1_10H.inpcrd
saveamberparm FBP1_15H FBP1_15H.prmtop FBP1_15H.inpcrd
saveamberparm FBP2_7H FBP2_7H.prmtop FBP2_7H.inpcrd
saveamberparm FBP2_10H FBP2_10H.prmtop FBP2_10H.inpcrd
saveamberparm FBP2_15H FBP2_15H.prmtop FBP2_15H.inpcrd

At this point we should really minimise the position of the extra proton that xleap added in each of these structures. We could do this by running 20 steps of minimisation using the IBELLY option to keep the rest of the protein fixed. However, since we are running everything manually here instead of as part of an automated script we will skip this step to save time.

We now need to find the energy of each of these protonated structures. This is easily done by doing exactly the same minimisation procedure as we did above but specifying only a single step:

This gives us the energy of the structure in each inpcrd file without doing any minimisation.

We will now run this on our 6 protonated structures we created above. This will allow us to fill in boxes C-H of our table. These calculations are so fast (as only 1 step is needed) that we can just have the results printed to standard output (the terminal) instead of creating an output file for each one. We can then just simply note down the final energy in each case. To write the output to the terminal we use the keyword 'stdout' in place of the specification of an output file when running sander.

$AMBERHOME/exe/sander -O -i min_1step.in -o stdout -p FBP1_7H.prmtop -c FBP1_7H.inpcrd

   NSTEP       ENERGY          RMS            GMAX         NAME    NUMBER
      1      -1.5081E+03     3.9648E+00     9.4896E+01     CD         81

 BOND    =       28.0200  ANGLE   =       62.8789  DIHED      =      304.8228
 VDWAALS =     -253.3751  EEL     =    -2559.3930  EGB        =     -892.9541
 1-4 VDW =      118.3223  1-4 EEL =     1683.5319  RESTRAINT  =        0.0000

$AMBERHOME/exe/sander -O -i min_1step.in -o stdout -p FBP1_10H.prmtop -c FBP1_10H.inpcrd
$AMBERHOME/exe/sander -O -i min_1step.in -o stdout -p FBP1_15H.prmtop -c FBP1_15H.inpcrd
$AMBERHOME/exe/sander -O -i min_1step.in -o stdout -p FBP2_7H.prmtop -c FBP2_7H.inpcrd
$AMBERHOME/exe/sander -O -i min_1step.in -o stdout -p FBP2_10H.prmtop -c FBP2_10H.inpcrd
$AMBERHOME/exe/sander -O -i min_1step.in -o stdout -p FBP2_15H.prmtop -c FBP2_15H.inpcrd

Structure 1	PROTEIN NO H (A) -1529.5 KCal/mol	PROTEIN GLU 7 H (C) -1508.1 KCal/mol	PROTEIN GLU 10 H (E) -1506.8 KCal/mol	PROTEIN ASP 15 H (G) -1507.5 KCal/mol
Structure 2	PROTEIN NO H (B) -1516.4 KCal/mol	PROTEIN GLU 7 H (D) -1493.3 KCal.mol	PROTEIN GLU 10 H (F) -1494.6 KCal/mol	PROTEIN ASP 15 H (H) -1494.1 KCal/mol
		SOLUTION GLU NO H (I)	SOLUTION GLU H (K)
		SOLUTION ASP NO H (J)	SOLUTION ASP H (L)