(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 2008

3.1) Building a Leap library file for FAM5

The third section of this tutorial describes how to build a leap library file for our non-standard residue (FAM5) using the charges derived in section 1. The FAM5 residue used in this tutorial was optimized in two fragments, however, because we "combined" these fragments using the Sirius visualization program we will only need to build one library file for the entire FAM5 molecule. This library file will contain all of the parameters and charges needed to run MD on our system. A brief overview of the following steps include loading FAM5 into the xleap program, defining the structure's topology and atom types, as well as inputting previously derived charges into the library file.

3.1.1) Defining the FAM5 Topology

Because we deleted the connectivity information in our pdb, we will need to manually tell xleap what atoms are bonded together. Fortunately, this is made easy with xleap's edit GUI, which allows us to draw in the bonding structure. We will need to navigate to the directory that houses our pdb file and execute the xleap program.

$AMBERHOME/exe/xleap &

The Following GUI should appear:

Load the fam5_edit.pdb file into xleap:

FAM5 = loadpdb "fam5_edit.pdb"

The following output should appear on the screen:

Loading PDB file: ./fam5_edit.pdb
Unknown residue: FAM5   number: 0   type: Terminal/last
..relaxing end constraints to try for a dbase match
  -no luck
Creating new UNIT for residue: FAM5 sequence: 1
Created a new atom named: C within residue: .R<FAM5 1>
Created a new atom named: O within residue: .R<FAM5 1>
Created a new atom named: C1 within residue: .R<FAM5 1>
Created a new atom named: C2 within residue: .R<FAM5 1>
Created a new atom named: H2 within residue: .R<FAM5 1>
Created a new atom named: C3 within residue: .R<FAM5 1>
Created a new atom named: C4 within residue: .R<FAM5 1>
Created a new atom named: O41 within residue: .R<FAM5 1>
Created a new atom named: O42 within residue: .R<FAM5 1>
Created a new atom named: C5 within residue: .R<FAM5 1>
Created a new atom named: C6 within residue: .R<FAM5 1>
Created a new atom named: H6 within residue: .R<FAM5 1>
Created a new atom named: C7 within residue: .R<FAM5 1>
Created a new atom named: H7 within residue: .R<FAM5 1>
Created a new atom named: C8 within residue: .R<FAM5 1>
Created a new atom named: C9 within residue: .R<FAM5 1>
Created a new atom named: C10 within residue: .R<FAM5 1>
Created a new atom named: H10 within residue: .R<FAM5 1>
Created a new atom named: C11 within residue: .R<FAM5 1>
Created a new atom named: H11 within residue: .R<FAM5 1>
Created a new atom named: C12 within residue: .R<FAM5 1>
Created a new atom named: O12 within residue: .R<FAM5 1>
Created a new atom named: C13 within residue: .R<FAM5 1>
Created a new atom named: H13 within residue: .R<FAM5 1>
Created a new atom named: C14 within residue: .R<FAM5 1>
Created a new atom named: O15 within residue: .R<FAM5 1>
Created a new atom named: C16 within residue: .R<FAM5 1>
Created a new atom named: C17 within residue: .R<FAM5 1>
Created a new atom named: C18 within residue: .R<FAM5 1>
Created a new atom named: H18 within residue: .R<FAM5 1>
Created a new atom named: C19 within residue: .R<FAM5 1>
Created a new atom named: H19 within residue: .R<FAM5 1>
Created a new atom named: C20 within residue: .R<FAM5 1>
Created a new atom named: O20 within residue: .R<FAM5 1>
Created a new atom named: C21 within residue: .R<FAM5 1>
Created a new atom named: H21 within residue: .R<FAM5 1>
Created a new atom named: N within residue: .R<FAM5 1>
Created a new atom named: H within residue: .R<FAM5 1>
Created a new atom named: C22 within residue: .R<FAM5 1>
Created a new atom named: H221 within residue: .R<FAM5 1>
Created a new atom named: H222 within residue: .R<FAM5 1>
Created a new atom named: C23 within residue: .R<FAM5 1>
Created a new atom named: H231 within residue: .R<FAM5 1>
Created a new atom named: H232 within residue: .R<FAM5 1>
Created a new atom named: C24 within residue: .R<FAM5 1>
Created a new atom named: H241 within residue: .R<FAM5 1>
Created a new atom named: H242 within residue: .R<FAM5 1>
Created a new atom named: C25 within residue: .R<FAM5 1>
Created a new atom named: H251 within residue: .R<FAM5 1>
Created a new atom named: H252 within residue: .R<FAM5 1>
Created a new atom named: C26 within residue: .R<FAM5 1>
Created a new atom named: H261 within residue: .R<FAM5 1>
Created a new atom named: H262 within residue: .R<FAM5 1>
Created a new atom named: C27 within residue: .R<FAM5 1>
Created a new atom named: H271 within residue: .R<FAM5 1>
Created a new atom named: H272 within residue: .R<FAM5 1>
Created a new atom named: O3* within residue: .R<FAM5 1>
  total atoms in file: 57

What happened?
Xleap did not recognize the FAM5 residue in its database, therefore it began creating a "UNIT" or new residue to define this unknown molecule. As part of this process xleap added all of our atoms to the UNIT. Because we are dealing with a non-standard residue this is exactly what we want. The UNIT created by xleap will allow the program to recognize FAM5 once we save our library file. It is important, however to double check that xleap "created" the correct atoms. If there is any discrepancy, you should check the pdb file you loaded into leap. The easiest discrepancy to spot is the total number of atoms in the file, the value for FAM5 should be 57.

We have now loaded the fluorescein molecule into xleap. However, as stated earlier, the fam5_edit.pdb file only contains information about atom coordinates. Therefore, xleap does not have information concerning fluorescein's parameters and connectivity, and we need to manually enter this information using the edit GUI provided with xleap.

Execute the edit GUI using xleap

edit FAM5

Before we begin, you should know how to navigate around the edit interface. The LEFT mouse button controls atom selection, the MIDDLE mouse button rotates the molecule and the RIGHT mouse button translates it. If the MIDDLE button does not work you can simulate it by holding down the Ctrl key and the LEFT mouse button. To make the image larger or smaller simultaneously hold down the MIDDLE and RIGHT buttons (or Ctrl and RIGHT) and move the mouse up to zoom in or down to zoom out. To help you out: carbon atoms are colored green, oxygen-red, hydrogen-white, and nitrogen-blue. Atoms that are currently selected appear in purple (you can deselect atoms by single clicking in open space, or holding down the shift key while rubber banding / selecting). You can use the erase tool to undo bonds/atoms accidentally created during the drawing process. Moreover, you can view atom names by selecting Display -> names. Once you have entered "draw" mode, you can bond the atoms together. To connect atoms, LEFT click the first atom and HOLD. While holding, drag the cursor to the second atom. The bond should "snap to" when you get close. Release the mouse button to complete the bond. The following images display the edit window before and after the topology information has been defined. **CAUTION... It is also important to note that when using xleap do not pess the X button to close a window. This action will terminate the session and any unsaved data will be lost. This will be especially unfortunate with the amount of data we will need to input. Instead use the File -> Close method. For the main edit window this can be accomplished using Unit -> Close.

3.2.1) Defining atom type and previously derived atomic charges

Now that we have defined the FAM topology, it is time to move on to the library file. As stated earlier we will need to define the atom types of each atom in our molecule. The list of atom types supported in the ff99SB force field is described in the force field papers and the files in the AMBER dat directory. We will define our atom types by analogy, comparing the atoms present in the manuscripts with those in our actual system. This analogy will take into account each atom with respect to its bonding and dihedral interactions. The parm99.dat file in $AMBERHOME/dat/leap/parm does a good job of describing the environment of each atom type and will allow us to make accurate comparisons. As we work through this process we will encounter an oxygen that is not properly defined by any of the parm99 atom types. So, we will create a new atom type describing this oxygen's bonding/dihedral interactions and define it in the library.

To define these parameters we will use the atom description table given to us in xleap. Select the entire FAM molecule by double clicking in the black region of the edit window. Once all the atoms have been selected select, Edit > Edit selected atoms. The following GUI should appear. Note: Mac users should select Window > Zoom from the X11 toolbar if you are having problems viewing the entire "Edit selected atoms" window.

Figure 3.4

We will input all of our atom types and charges in this window. Due to the large number of atoms in our system, we will not go through every single atom comparison in the tutorial. We will, however, focus on a few atoms; giving you the option to analyze the rest on your own. For the purpose of this tutorial we will provide you with a file that contains all the atom types and charges ready for you to input into the "edit selected atoms" window.

Atom Comparison #1:

In our first atom comparison we will look at C10, which is atom #17 in the fam5_edit.pdb file. This atom is part of the xanthene ring, and resembles the average carbon atom present in a benzene ring, therefore we assign an atom type appropriate for an aromatic carbon, CA.

Atom Comparison#2:

In our second comparison we will look at atom H10 which is atom #18 in the fam5_edit.pdb file. This atom is a hydrogen atom bonded to a aromatic carbon. As such we will assign it an atom type appropriate for this type of hydrogen, in this case, HA.

Atom Comparison#3:

For the third example consider atoms O41 and O42. These are equivalent oxygen atoms that share a resonance structure between COO- and CO-O. These atoms are carbonyl oxygens and are similar to the oxygens found in the protein backbone but actually are a better match to the oxygens in a deprotonated aspartic acid sidechain. In this case we assign these atoms to be of type O2.

Section 3.2.2) Defining the atom type of a non-analogous Oxygen atom

The oxygen atom in the xanthene O15 is aromatic in character and therefore does not correspond well with any of the oxygen atom types that already exist in FF99SB as they are all either SP2 type carbonyl or SP3 type ether oxygens. Deriving parameters for a new atom type is beyond the scope of this tutorial but fortunately parameters for this particular atom type have been published (VanBeek et al, Biophys. J. 2007, 92, 4168-4178 [pdf] [supp mat]). Therefore we will just make use of these parameters. At this stage all we need worry about is assigning the correct atom type which in this case simply needs to be a unique type that does not exist in the current force field files. Following the terminology in the VanBeek paper we will assign it type OA.

The remaining atoms types were assigned as shown in the figure below.

Figure 3.5

Section 3.3.1) Saving our new library file

Now that we have inputted all of the parameters necessary for xleap to adequately recognize our FAM5 residue, we can save a library file which will enable xleap to recognize this residue in the future. This is very important so that we do not have to repeat all of the above steps each time we import FAM5 into the xleap program. To do this type the following command. It is also important to save a new pdb at this step

saveoff FAM5 fam5.lib
savepdb FAM5 fam5_leap.pdb

Section 3.4.1) Creating a frcmod file.

A frcmod file is required to define the mass and VDW parameters for the new OA atom type we added and also to provide all of the bond, angle and dihedral parameters that are not present in the standard FF99SB force field.  We will use leap to help us identify the missing parameters.  This will require a couple of passes as the leap check function only identifies missing bond and angle parameters and to look for missing dihedrals we will need to try to save a prmtop and inpcrd file after we have provided the missing bond and angle parameters.  If you have leap open you should close it now (being sure that you saved your lib and pdb files). Since we don't need the graphical front end for this we can just use tleap which is faster and more convenient.

$AMBERHOME/exe/tleap -s -f $AMBERHOME/dat/leap/cmd/oldff/leaprc.ff99SB

We can then load the lib file we created above:

>loadoff fam5.lib

Now we can check the fam5 unit: 

>check FAM5

and Leap will return the missing bond and angle parameters:

Checking 'FAM5'....
ERROR: The unperturbed charge of the unit: -2.307860 is not integral.
WARNING: The unperturbed charge of the unit: -2.307860 is not zero.
ERROR: The perturbed charge: -2.307860 is not integral.
WARNING: The perturbed charge: -2.307860 is not zero.
Checking parameters for unit 'FAM5'.
Checking for bond parameters.
Could not find bond parameter for: CA - O
No bond parameter for: CA - O
Could not find bond parameter for: CA - OA
No bond parameter for: CA - OA
Could not find bond parameter for: OA - CA
No bond parameter for: OA - CA
Could not find bond parameter for: CA - O
No bond parameter for: CA - O
Checking for angle parameters.
Could not find angle parameter: O - C - CA
Can't find angle parameter: O - C - CA
Could not find angle parameter: CA - C - N
Can't find angle parameter: CA - C - N
Could not find angle parameter: CA - C - O2
Can't find angle parameter: CA - C - O2
Could not find angle parameter: CA - C - O2
Can't find angle parameter: CA - C - O2
Could not find angle parameter: CA - CA - OA
Can't find angle parameter: CA - CA - OA
Could not find angle parameter: CA - CA - O
Can't find angle parameter: CA - CA - O
Could not find angle parameter: O - CA - CA
Can't find angle parameter: O - CA - CA
Could not find angle parameter: CA - CA - OA
Can't find angle parameter: CA - CA - OA
Could not find angle parameter: CA - OA - CA
Can't find angle parameter: CA - OA - CA
Could not find angle parameter: OA - CA - CA
Can't find angle parameter: OA - CA - CA
Could not find angle parameter: OA - CA - CA
Can't find angle parameter: OA - CA - CA
Could not find angle parameter: CA - CA - O
Can't find angle parameter: CA - CA - O
Could not find angle parameter: O - CA - CA
Can't find angle parameter: O - CA - CA
There are missing parameters.
check: Errors: 2 Warnings: 2

The first 2 errors and warnings came about because our FAM5 unit does not have an integer charge.  As discussed earlier, this is what we intended since 5' residues in AMBER are not supposed to have an integer charge.  The list that follows are all of our missing parameters though there is quite a bit of redundancy.  For instance the last two parameters that are reported are CA-CA-O and O-CA-CA.  These are equivalent angles, differing only in the order of the atoms.  So, removing these redundancies leaves us with the following list of parameters that we need to define:

CA - O
CA - OA
CA - C - O  (re-ordered from O - C - CA)
CA - C - N
CA - C - O2
CA - CA - OA
CA - CA - O
CA - OA - CA

Coming up with values for these parameters is non-trivial.  If we're lucky, these missing parameters describe atoms in environments that are very similar to existing parameters and we can assign them 'by anology.'  For instance there exists a ca - ca - o parameter in the generalized amber force field that we could use for the CA - CA - O parameter.  This kind of parameter assignment must be done with care and a thorough description is beyond the scope of this tutorial.  Taking parameters from the VanBeek paper mentioned earlier (VanBeek et al, Biophys. J. 2007, 92, 4168-4178 [pdf] [supp mat]) we come up with the following first pass at a frcmod file, fam5_incomplete.frcmod.

From VanBeek et al. Biophys J. (2007) 92, 4168-4178
MASS
OA            16.00        0.465    same as gaff os and parm99

BOND
CA-OA     372.40    1.373      same as gaff ca-os
CA- O       648.00    1.214      same as gaff c -o

ANGLE
CA-CA-OA   69.800   119.200    same as gaff ca-ca-os
CA-CA-O      70.900   123.430   same as gaff ca-ca-o
CA-C -O        68.700   123.440   same as ca-c -o in gaff
CA-C -O2      68.700   123.440   same as ca-c -o in gaff
CA-OA-CA   63.600    118.960   same as GAFF ca-os-ca
CA-C -N        69.400    112.03    same as GAFF ca-c -n

DIHE

IMPROPER

NONBON
OA 1.6837 0.1700 OPLS ether e.g. same as gaff os

Note that these parameters provide a reasonable description of the molecule, but if highly accurate behavior of the fluorescein is needed, these parameters would need to be refined.  Note also that we have not included any dihedral parameters as we have not yet figured out which are missing.  Finally, even though the check command did not provide a warning message about the missing atom type OA, we know that this will be a problem eventually, so we have included the MASS and NONBON parameters for this atom type as well.

We now load the fam5_incomplete.frcmod file and then find out what dihedral parameters need to be defined:

>loadamberparams fam5_incomplete.frcmod

Unfortunately the check command does not check dihedrals so instead we will just attempt to save a prmtop and inpcrd file which will fail listing the missing dihedral parameters.

>saveamberparm FAM5 prmtop inpcrd

and Leap will return the missing bond and angle parameters:

 ** No torsion terms for  CA-CA-OA-CA
 ** No torsion terms for  CA-CA-OA-CA
 ** No torsion terms for  CA-OA-CA-CA
 ** No torsion terms for  CA-OA-CA-CA

In this case there is actually only one missing dihedral (CA-CA-OA-CA) since the above four all describe the same dihedral angle.  So we just need to add one entry to the DIHE section of our frcmod file:

CA-CA-OA-CA   4   14.500   180.000   2.000   Same as gaff X-ca-ca-X and because M. Zwier found ~18 when hand fitting

Use a text editor to update the fam5_incomplete.frcmod file and save it as fam5.frcmod.  Now load this into Leap and attempt the saveamberparm command again just as a final check; we won't actually use the resulting prmtop and inpcrd file.

We have successfully saved our FAM5 library file and created the associated frcmod file and are now ready to attach it to our polyAT decamer.

---

| INDEX | SECTION 1 | SECTION 2 | SECTION 3 | SECTION 4 |

---

*Funding and computational support for the creation of this tutorial was provided by NSF-CIEG (BDI0726924), NSF-REU, NSF-MRI, HHMI and ACS-PRF.

(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 2008