Parameterising your ligands
Overview
Teaching: 15 min
Exercises: 0 minQuestions
How do I generate parameters for non-conventional residues?
Objectives
Generate parameters for GWS ligand using
antechamberandparmchk2.
In this tutorial, we will make use of the antechamber and parmchk2 tools from AmberTools package to generate AMBER force field parameters. We will use the “general AMBER force field 2 (GAFF2)” parameters. This force field was designed to provide atom types and atom parameters needed to parameterise most pharmaceutical molecules, and is compatible with the traditional AMBER force fields for proteins. To assign partial charges, we are going to use AM1-BCC2 charge method. It is an inexpensive and fast method for calculating partial charges and be we don’t need to worry too much about its limitations as partial charges won’t be used during production runs due to the QM treatment of the ligand in the QM/MM simulation.
Generating point charges and assigning atomtypes with antechamber
antechamber allows the rapid generation of topology files for use with the AMBER simulation programs. It is a higher-level wrapper around other AmberTools programs (sqm, divcon, atomtype, am1bcc, bondtype, espgen, respgen and prepgen) and is able to:
- Automatically identify bond and atom types
- Judge atomic equivalence
- Generate residue topology files
- Find missing force field parameters and supply reasonable suggestions
We will use antechamber to assign atom types to the GWS ligand and calculate a set of point charges. We will have to specify GAFF2 force field (-at gaff2, more information on the parameters here: $AMBERHOME/dat/leap/parm/gaff2.dat), the AM1-BCC2 charge method (-c bcc), the net charge of the molecule (-nc 0) and the name of the new residue generated (We are going to keep GWS: -rn GWS). We are going to output a mol2-type file (-fo mol2) containing the atomtypes and the point charges.
antechamber -i GWS.H.pdb -fi pdb -o GWS.mol2 -fo mol2 -c bcc -nc 0 -rn GWS -at gaff2
Welcome to antechamber 20.0: molecular input file processor.
acdoctor mode is on: check and diagnosis problems in the input file.
-- Check Format for pdb File --
Status: pass
-- Check Unusual Elements --
Status: pass
-- Check Open Valences --
Status: pass
-- Check Geometry --
for those bonded
for those not bonded
Status: pass
-- Check Weird Bonds --
Status: pass
-- Check Number of Units --
Status: pass
acdoctor mode has completed checking the input file.
Info: Total number of electrons: 118; net charge: 0
Running: /home/d118/d118/shared/amber20/bin/sqm -O -i sqm.in -o sqm.out
This step generates a lot of intermediate files (in capital letters). You can safely delete them as these are used in antechamber and we no longer require them. The sqm.xxx files are input and output from the sqm quantum mechanics code used by antechamber to calculate the atomic point charges. We are not interested in the data here except to check that the sqm calculation completed successfully.
tail sqm.out
--------- Calculation Completed ----------
Once we are sure that the point charges were generated The MOL2 file contains tridimensional information of the molecule as well as atom type and point charges. You can find more information on the MOL2 format here. The resulting MOL2 file is shown below:
<TRIPOS>MOLECULE
GWS
34 35 1 0 0
SMALL
bcc
@<TRIPOS>ATOM
1 N1 7.5340 0.3440 17.5810 nb 1001 GWS -0.643000
2 C4 13.4880 -0.5140 24.0620 c3 1001 GWS -0.075900
3 C5 13.9760 -1.6820 23.2210 c3 1001 GWS -0.078400
...
TIP
The GAFF and GAFF2 atom types are in lower case. This is the mechanism by which the GAFF force field is kept independent from the macromolecular AMBER force fields. All of the traditional AMBER force fields use uppercase atom types. In this way the GAFF and traditional force fields can be mixed in the same calculation.
Looking for missing force field parameters with parmchk2
While the most likely combinations of bond, angle and dihedral parameters are defined in the parameter file it is possible that our molecule might contain combinations of atom types for bonds, angles or dihedrals that have not been parameterised. If this is the case, we will have to specify any missing parameters before we can create our prmtop and inpcrd files in LEap.
We will use parmchk2 to test if all the parameters we require are available.
parmchk2 -i GWS.mol2 -f mol2 -o GWS.frcmod -s gaff2
parmchk2 generates a parameter file that can be loaded into Leap in order to add missing parameters. It contains all of the missing parameters. If possible, it will fill in these missing parameters by analogy to a similar parameter. Let’s look at the generated GWS.frcmod file:
Remark line goes here
MASS
BOND
ca-ns 315.20 1.412 same as ca- n, penalty score= 0.0
c -ns 356.20 1.379 same as c- n, penalty score= 0.0
hn-ns 527.30 1.013 same as hn- n, penalty score= 0.0
ANGLE
c -ns-ca 65.700 123.710 same as c -n -ca, penalty score= 0.0
ca-ns-hn 48.000 116.000 same as ca-n -hn, penalty score= 0.0
ns-c -o 113.800 123.050 same as n -c -o , penalty score= 0.0
c -ns-hn 48.700 117.550 same as c -n -hn, penalty score= 0.0
ca-ca-ns 85.600 120.190 same as ca-ca-n , penalty score= 0.0
c3-c -ns 84.300 115.180 same as c3-c -n , penalty score= 0.0
DIHE
o -c -ns-ca 4 10.000 180.000 2.000 same as X -c -n -X , penalty score= 0.0
c3-c -ns-ca 1 0.750 180.000 -2.000 same as c3-c -n -ca
c3-c -ns-ca 1 0.500 0.000 3.000 same as c3-c -n -ca, penalty score= 0.0
o -c -ns-hn 1 2.500 180.000 -2.000 same as hn-n -c -o
o -c -ns-hn 1 2.000 0.000 1.000 same as hn-n -c -o , penalty score= 0.0
ca-ca-ns-c 1 0.950 180.000 2.000 same as c -n -ca-ca, penalty score= 0.0
ns-c -c3-c3 1 0.000 180.000 -4.000 same as n -c -c3-c3
ns-c -c3-c3 1 0.710 180.000 2.000 same as n -c -c3-c3, penalty score= 0.0
ca-ca-ns-hn 4 1.800 180.000 2.000 same as X -ca-n -X , penalty score= 0.0
c3-c -ns-hn 4 10.000 180.000 2.000 same as X -c -n -X , penalty score= 0.0
IMPROPER
ca-ca-ca-ns 1.1 180.0 2.0 Using the default value
ca-ca-ca-ha 1.1 180.0 2.0 Using general improper torsional angle X- X-ca-ha, penalty score= 6.0)
c3-ns-c -o 10.5 180.0 2.0 Using general improper torsional angle X- X- c- o, penalty score= 6.0)
c -ca-ns-hn 1.1 180.0 2.0 Same as X -X -n -hn, penalty score= 6.0 (use general term))
ca-h4-ca-nb 1.1 180.0 2.0 Same as X -X -ca-ha, penalty score= 44.3 (use general term))
NONBON
TIP
You should check these parameters carefully before running a simulation. If antechamber can’t empirically calculate a value or has no analogy it will either add a default value that it thinks is reasonable or alternatively insert a place holder (with zeros everywhere) and the comment “ATTN: needs revision”. In this case you will have to manually parameterise this yourself. See the links at the beginning of the section.
FURTHER INFORMATION:
This step is not trivial and there are many ways to do it. In this example, as we are aiming to use a QM description of the ligand, we can use a simple approach. If you have trouble parameterising your ligands here are some useful links:
Key Points
Generating parameters for non-conventinal residues is not easy and there are different ways to parameterise them.
The quality of the parameters is important to obtain valid results.
antechambercreates parameters for a non-conventional residues.
parmchk2checks for missing parameters in the forcefield.ALWAYS check the parameters for your system and test that they behave as expected.