Re: [AMBER] ATOMIC NUMBER section in prmtop file - sulfur atom in ligand is -1?

From: He, Amy via AMBER <>
Date: Mon, 5 Sep 2022 17:21:17 +0000

Dear Dr. Case,

Thank you so much for your thoughtful reply! That was really helpful. I thought the problem was relevant to hydrogen mass partitioning, but as you correctly pointed out, it is likely more about the parametrization of the sulfur atom.
But I’m still not so sure what to do about the problem… And I would really appreciate it if you could give some more advice on it. My apologies for the lengthy message in advance.

1. New atom type…?
I did notice that the ligand bisulfite sulfur could be a new atom type. The sulfur has atom type SO in my prepin file, as assigned from Antchamber, which is not in $AMBERHOME/dat/leap/parm/parm10.dat, but in $AMBERHOME/dat/antechamber/ATOMTYPE_AMBER.DEF.

Therefore, Antechamber didn’t reject the atom type. In addition, I was able to use parmchk2 so that the relevant parameters are complemented with the GAFF2 equivalents.

However, without declaring the new atom type in Leap, there would be a message in the standard output:
# +Currently only Sp3-Sp3/Sp3-Sp2/Sp2-Sp2 are supported
# +---Tried to superimpose torsions for: *-C55-S58-*
# +--- With Sp3 - Sp0
# +--- Sp0 probably means a new atom type is involved
# +--- which needs to be added via addAtomTypes

Still, Leap was able to write a prmtop file for my system, though it seems like it could not recognize the hybridization state of the sulfur… I disregarded this message the first time I ran the calculation, because I thought this torsion had been taken care of in the GAFF2.frcmod I made:
H2-CT-SO-O 9 1.300 0.000 3.000 same as X -c3-s6-X , penalty score= 0.0

2. My Leap input file
cat << EOF >
source leaprc.protein.ff14SB
source leaprc.gaff2
source leaprc.water.tip3p
addAtomTypes {
    { "SO" "S" "sp3" }
set default PBRadii mbondi3
loadamberprep RES.prepin
loadamberparams RES_gaff2.frcmod
loadamberparams RES_ff14SB.frcmod
complex = loadpdb complex.pdb4amber.pdb
saveamberparm complex prmtop.prmtop inpcrd.inpcrd

3. Difference in prmtop files
I used a compare tool to examine the prmtop files. I found that addAtomTypes changed ATOMIC_NUMBER, RADII, and SCREEN in the prmtop.
ATOMIC_NUMBER: from -1 to 16
RADII: from 1.5 to 1.8
SCREEN: from 0.8 to 0.96

4. Results of the minimization
I found a big difference in the EEL term, with or without adding the new atom type in Leap!!!
(1) No Leap addAtomTypes, No ParmED HMR
   1000 -1.0844E+04 9.4968E-02 1.5259E+00 CE2 3311

BOND = 164.5286 ANGLE = 594.5273 DIHED = 3359.7041
VDWAALS = -2761.2412 EEL = -21251.9696 EGB = -4408.7354
1-4 VDW = 934.7608 1-4 EEL = 11987.3438 RESTRAINT = 0.0000
ESURF = 536.9040

(2) Leap addAtomTypes, No ParmED HMR
   1000 -7.2047E+04 2.9975E+00 1.0349E+02 O60 4654

BOND = 316.5836 ANGLE = 1018.8918 DIHED = 3472.6794
VDWAALS = -2755.0830 EEL = -83097.9364 EGB = -4361.2372
1-4 VDW = 865.3401 1-4 EEL = 11956.3905 RESTRAINT = 0.0000
ESURF = 537.265

With adding the new atom type in Leap, the endpoint structure appeared to be really distorted, just like it was with the ParmED prmtop file. And I thought that was not just because of the changes in RADII and SCREEN.
So I was kind of thinking, maybe ATOMIC NUMBER is -1 is the correct option (and it is necessary) for the sulfur atom, because it doesn’t have a type in parm10, and has been overridden by its equivalent in GAFF2?? Or maybe I did something that was serious wrong in the parameterization :’(

I’m sharing my input files and commands:
Ligand pdb as input for Antechamber
Complex pdb as input for Leap
Bash commands I used to execute programs and generate other inputs

I hope the formatting of this message is readable.. The outcome should be reproducible with the shared files. Once again, massive thanks to you for your time and kind advice.

Many Thanks,

Amy He
Chemistry Graduate Teaching Assistant
Hadad Research Group
Ohio State University
From: David A Case <>
Date: Monday, September 5, 2022 at 9:14 AM
To: He, Amy <>, AMBER Mailing List <>
Subject: Re: [AMBER] ATOMIC NUMBER section in prmtop file - sulfur atom in ligand is -1?
On Mon, Sep 05, 2022, He, Amy via AMBER wrote:
>One is the original prmtop I generated by Leap. In the ATOMIC NUMBER
>section, the ligand sulfur atom is -1.
I'm guessing(?) that you defined a new sulfur atom type for the ligand,
but neglected to have a addAtomTypes command in your tleap script:
addAtomTypes {
    { "xx"   "S" "sp3" }
Where "xx" is your new atom type.  Still, this should not make any
difference, since (I think!) Amber only really uses the atomic number for
QM/MM calculations, not for conventional force fields.
Also, the masses are irrelevant for minimization.  I think we would need
to have files to reproduce the problem you are seeing (for a small test
case).  It sounds like my assumption that the atomic number section is not
relevant is probably wrong.
You might be able to make some progress by looking the breakdown of energy
terms for the starting structure, using both prmtop files.  Are all fields
the same except for the VDW term?
If you did leave out the addAtomTypes entry, you could try to fix that, then
regenerate the prmtop file with tleap.  Maybe something else is different as
well.  Use the parmed commands printLJTypes and printLJMatrix to examine
your nonbonded interactions.  Look especially for things related to the
sulfur atom.
...good luck...dac
AMBER mailing list
Received on Mon Sep 05 2022 - 10:30:02 PDT
Custom Search