[gmx-users] Converting Amber frcmod and mol2 files to Gromacs
Simon Dürr
simon.duerr at uni-konstanz.de
Wed Sep 30 19:37:56 CEST 2015
Hi,
we had the pleasure to convert the *.mol2 and *.frcmod for the HEME
and CYP prosthetic groups from the following paper for use with
GROMACS:
Shahrokh, K., Orendt, A., Yost, G. S. and Cheatham, T. E. (2012),
Quantum mechanically derived AMBER-compatible heme parameters for
various states of the cytochrome P450 catalytic cycle. J. Comput.
Chem., 33: 119–133. doi: 10.1002/jcc.21922
As the protocol for doing so is rather lengthy we want to share our
protocol with you guys.
If you can spare some time we would be grateful If some of you could
check our procedure.
If you find any errors let us know. However so far the trajectories
look good (angles and bond distances are good).
After further testing of the trajectories we want to release the
rtp/hdb/itp - files.
Where would be the best way to upload them?
The Protocol is as follows:
#################################
## Table of Contents ##
#################################
Part 1: RTP-file
Part 2: ffbondend-file
Part 3: ffnonbonded
Part 4: residuetypes.dat
Part 5: specbond.dat
Part 6: hdb-file
Part 7: atomtypes.atp
##########################
## PART 1 ##
##########################
##### BEGIN RTP-File #######
[ HEME ]
[ atoms]
NC ncr 0.00660 1
C1C ccr -0.01920 2
C4C ccr -0.01930 3
C2C ccr 0.01140 4
C3C ccr -0.03570 5
…
[ bonds ]
NC C1C
NC C4C
NC FE
[…]
###### END RTP-FILE #####
The information to fill this file is found in the HEM.mol2 file.
###### HEM.mol2 ######
@<TRIPOS>ATOM
1 NC 2.946000 -0.972000 -0.732000 nc 1 HEM 0.0066 ****
2 C1C 4.282000 -0.756000 -0.518000 cc 1 HEM -0.0192 ****
4 C2C 5.025000 -1.989000 -0.700000 cc 1 HEM 0.0114 ****
5 C3C 4.097000 -2.957000 -1.046000 cc 1 HEM -0.0357 ****
[…]
@<TRIPOS>BOND
1 1 2 1
2 1 3 1
3 1 28 1
###### END HEM.mol2 #####
To convert to Gromacs
1. take the uppercase name in column 2 of the mol2 which is the name
found in the PDB input file and put it in column 1 of the [atoms]
section of the .rtp
2. Add the lowercase identifier for GAFF (column 6 in the mol2 ) to
the second column of the atoms section.
3. As we don't want to overwrite any of the normal definitions from
AMBER99SB-ILDN with GAFF definitions we use a unique name for the GAFF
atoms. Here we just added "r" at the end of each atom name.
4. Add the column 6 of the mol2 with the charge as third column to the
atoms section
5 Group the atoms to charge groups in the fourth column of the
[atoms]-section. We put every atom in a new group.
Look at the bonds section in the mol2:
The first column is the index of the bond, the second column is the id
of the first atom of the bond, the third column the id of the second
atom of the bond. The last is the bond type (single or double bond)
7.Take the id of the first atom and replace it with the uppercase
identifier from the ATOM section of the mol2. Do the same for the id
of the second atom.
The columns 1,3,4,5,7 and 8 in the ATOM section are not relevant for Gromacs.
##########################
## PART 2 ##
##########################
##### BEGIN FFBONDED.ITP #####
[ bondtypes ]
; i j func b0 kb
fer ncr 1 0.2 83680 ; IC6.frcmod
[..]
[ angletypes ]
; i j k func th0 cth
ncr fer ndr 1 85.521 1104.576 ; IC6.frcmod
[..]
[ dihedraltypes ]
;i j k l func phase kd pn
X SH fer X 9 180.0 37.65600 2 ; IC6.frcmod
[..]
[ dihedraltypes ]
X orr crr orr 4 180.0 4.60240 2 ; GAFF.dat improper
#### END FFBONDED.ITP #####
The information for this file can be found in the IC6.frcmod file OR
in the gaff.dat
#### IC6.frcmod #####
BOND
fe-nc 100.000 2.000
[..]
ANGLE
nc-fe-nd 132.000 85.521 # average angle
[..]
DIHEDRAL
X -SH-fe-X 1 9.000 180.000 2.000
[..]
IMPROPER
NONBON
fe 1.8 0.01 1.0
#### END IC6.frcmod ####
#### BEGIN GAFF.dat ###
[…]
X -c -c -X 4 1.200 180.000 2.000
X -o -c -o 1.1 180. 2. JCC,7,(1986),230
#### END GAFF.dat #####
To convert to Gromacs
1. Add the bonds section:
a) Copy the two names of the atoms of the bond. Rename to new
identifier and separate by a single space.
b) Set func to 1 (single bond)
c) Convert the b0 value from the Amber b0 value(in angstroms) as follows
b0(Gro) = b0(amber) / 10 = [nm]
d) Convert the kb value from the Amber kb value as follows:
kb(Gro) = kb(amber) /10 *2 *4.184 = [joule]
2. Add the Angle section
a) Copy the three names of the atoms of the angle. Rename to new
identifier and separate by a single space.
b) Set func to 1 (single bond)
c) Copy the angle from amber
d) Calculate cth as follows
E(amber) = Emin/2
cth = E(amber) *2 * 4.184 = [joule]
3. Add the dihedrals.
First add all the dihedrals specified in the frcmod. Then use
Ambertools to generate a coordinate and topology file with the loaded
mol2 and frcmod files but don't load the GAFF parameters. Ambertools
will tell you which impropers and dihedrals you have to add from the
GAFF.dat because it could not find them in the frcmod.
For each dihedral/improper do this:
a) Copy the 4 names of the atoms of the dihedral. Rename the GAFF
residues to your new identifier. Don't rename X and residues in
uppercase
b) Set func to 9 for a dihedral, set func to 4 for improper
c) Take the phase from the fourth column of the Amber file
d) Calculate the kd value:
For dihedral (IDIVF is the second column in the Amber file)
kd(gro) = PK(amber) * 4.184 / IDIVF(amber) = [joule]
For impropers (IDIVF is not set)
kd(gro) = PK(amber) *4.184 = [joule]
d) Copy the pn value from the fifth column of the Amber file
##########################
## PART 3 ##
##########################
#### BEGINf ffnonbonded.ITP ####
[ atomtypes ]
; name at.num mass charge ptype sigma epsilon
ncr 7 14.01 0.0000 A 3.25000E-01 7.11280E-01 ; GAFF.dat
ccr 6 12.01 0.0000 A 3.39967E-01 3.59824E-01 ; GAFF.dat
[...]
#### END ffnonbonded.ITP ####
#### BEGIN Mol.2 ####
@<TRIPOS>ATOM
1 NC 2.946000 -0.972000 -0.732000 nc 1 HEM 0.0066 ****
2 C1C 4.282000 -0.756000 -0.518000 cc 1 HEM -0.0192 ****
#### END Mol2 ####
#### BEGIN GAFF.dat ####
MOD4 RE
cc 1.9080 0.0860 OPLS
#### END GAFF.dat
To convert to Gromacs:
a) Copy the lowercase name of the atom in column six from the mol2 and
rename to your new unique identifier.
b) Set atomic number from the periodic table. Easy ;)
c) Set mass from the periodic table or just copy it from other
definitions in this file from the same element
d) Set Charge to 0.000 which is the default value in Amber in Gromacs
in this file.
e) Set ptype to A
f) Calculate sigma from the Rmin/2 value (second column in Gaff.dat
MOD4 section)
sigma = Rmin(amber) * 1/(2^(1/6)) *2 * 1/10
g) Calculate Epsilon from the Epsilon value of Amber (third column in
Gaff.dat MOD4 section)
epsilon = epsilon(amber) *4.184 = [joule]
##########################
## PART 4 ##
##########################
#### BEGIN Residuetypes.dat ####
GLY Protein
HEME Protein
[..]
#### END residuetypes.dat ####
Add your new resname with the PDB Name here and assign it to the right group.
##########################
## PART 5 ##
##########################
#### BEGIN specbond.dat ####
[..]
resA atomA nbondsA resB atomB nbondsB length newresA newresB
CYS SG 1 HEM FE 2 0.2 CYP
HEME
[..]
#### END specbond.dat
Rename special residues like the CYS connected to the heme iron.
Gromacs detects what it has to rename by the bond length which is the
parameter in column 7 (in nanometers).
For each of your residues you need a line here.
##########################
## PART 6 ##
##########################
#### BEGIN HDB-File ####
HEME 16
;nH type Name C-Atom Controlatoms
1 1 HHD CHD C1D C4C
1 1 HAC CAC C3C CBC
2 3 HBC CBC CAC C3C
3 4 HMC CMC C2C C1C
[...]
#### END HDB FIle #####
There is no equivalent for this file in Amber. Just open the PDB file
for your residue in a viewer(like UCSF Chimera) and let it display
the atom names. What you now have to do is to add each hydrogen bond
to its C atom and specify two control atoms.
The index after the residue name(HEME 16) is the number of times a new
hydrogen is defined.
For each C atom you now have to add the hydrogens.
Type can either be:
1: Planar hydrogens
2: single hydrogen
3: two planar hydrogens -CR=CH2
4: two or three tetrahedral hydrogens
Others see: Gromacs Manual , chapter 5.6.4
##########################
## PART 7 ##
##########################
#### BEGIN ATOMTYPES.ATP ####
ncr 14.01000 ; GAFF.dat Sp2 N in non-pure aromatic systems
ccr 12.01000 ; GAFF.dat Sp2 carbons in non-pure aromatic systems
cgr 12.01000 ; GAFF.dat Inner Sp carbons in conjugated systems
Add the same mass for the new atomtypes from the ffbonded.itp
#### END atomtypes.atp ####
##########################
## Finalization ##
##########################
After having done all this make a new folder in the Gromacs share
folder to add your new FF as a duplicate of the old one with the new
additions added.
Closely examine the Gromacs log for errors or warnings.
More information about the gromacs.org_gmx-users
mailing list