Pdb2gmx: Generating an ATP topology fails with gromos54a7.ff

GROMACS version: 2022.3
GROMACS modification: No

Hi,

I want to simulate ATP using the GROMOS gromos54a7.ff forcefield that comes shipped with the GROMACS installation.

First of all, I realized that the ATP defined in the gromos54a7.ff/aminoacids.rtp is the 3- charged molecule, but I want to use the 4- charged one. Is there a reason that the 4- charged one is not defined and what can I do to build a topology for the 4- one?

Secondly, even if I try to go for the 3- version as defined in gromos54a7.ff, pdb2gmx fails with the following error

atom AH3PG is missing in residue ATP 653 in the pdb file

atom AH3PG is missing in residue ATP 653 in the pdb file

-------------------------------------------------------
Program:     gmx pdb2gmx, version 2022.3
Source file: src/gromacs/gmxpreprocess/pdb2top.cpp (line 1584)

Fatal error:
There were 1 missing atoms in molecule Other, if you want to use this
incomplete topology anyhow, use the option -missing

My call was: gmx pdb2gmx -f atp.gro -water spce -v
So somehow the Hydrogen atom at the gamma phosphate oxygen is not generated although defined in the rtp file.

Interestingly, when I use the older forcefields gromos53a5.ff or gromos53a6.ff, the topology can be generated successfully.

These are my .gro files that I use as input. Note that the atom names are different between gromos 53 and gromos 54:

atp.gro:
I do the call ‘gmx pdb2gmx -f atp.gro -water spce -v’ and select ‘14: GROMOS96 54a7 force field’.
It throws the Error I described above.

Gnomes, ROck Monsters And Chili Sauce
   31
  653ATP    C4'    1  13.945   1.888   5.103
  653ATP    O4'    2  14.043   1.890   5.210
  653ATP    C1'    3  13.989   1.831   5.327
  653ATP     C5    4  14.117   1.995   5.626
  653ATP     N7    5  14.062   2.102   5.560
  653ATP     C8    6  14.015   2.051   5.449
  653ATP     N9    7  14.035   1.917   5.437
  653ATP     N1    8  14.224   1.867   5.785
  653ATP     C2    9  14.202   1.764   5.708
  653ATP     N3   10  14.143   1.757   5.590
  653ATP     C4   11  14.102   1.879   5.555
  653ATP     C6   12  14.182   1.990   5.749
  653ATP     N6   13  14.197   2.104   5.824
  653ATP    C2'   14  13.838   1.845   5.314
  653ATP    O2'   15  13.771   1.735   5.374
  653ATP    C3'   16  13.815   1.844   5.166
  653ATP    O3'   17  13.789   1.710   5.127
  653ATP    C5'   18  13.937   2.026   5.037
  653ATP    O5'   19  13.900   2.130   5.129
  653ATP     PA   20  13.886   2.281   5.078
  653ATP    O1A   21  14.013   2.323   5.013
  653ATP    O2A   22  13.831   2.367   5.186
  653ATP    O3A   23  13.777   2.263   4.964
  653ATP     PB   24  13.619   2.211   4.964
  653ATP    O1B   25  13.613   2.060   4.988
  653ATP    O2B   26  13.537   2.304   5.055
  653ATP    O3B   27  13.580   2.240   4.805
  653ATP     PG   28  13.487   2.357   4.734
  653ATP    O1G   29  13.563   2.488   4.757
  653ATP    O2G   30  13.479   2.316   4.587
  653ATP    O3G   31  13.354   2.349   4.810
   0.00000   0.00000   0.00000

atp_53.gro:
I do the call ‘gmx pdb2gmx -f atp_53.gro -water spce -v’ and select ‘13: GROMOS96 53a6 force field (JCC 2004 vol 25 pag 1656)’
It generates the topology successfully

Gnomes, ROck Monsters And Chili Sauce
   31
  653ATP   AC4*    1  13.945   1.888   5.103
  653ATP   AO4*    2  14.043   1.890   5.210
  653ATP   AC1*    3  13.989   1.831   5.327
  653ATP    AC5    4  14.117   1.995   5.626
  653ATP    AN7    5  14.062   2.102   5.560
  653ATP    AC8    6  14.015   2.051   5.449
  653ATP    AN9    7  14.035   1.917   5.437
  653ATP    AN1    8  14.224   1.867   5.785
  653ATP    AC2    9  14.202   1.764   5.708
  653ATP    AN3   10  14.143   1.757   5.590
  653ATP    AC4   11  14.102   1.879   5.555
  653ATP    AC6   12  14.182   1.990   5.749
  653ATP    AN6   13  14.197   2.104   5.824
  653ATP   AC2*   14  13.838   1.845   5.314
  653ATP   AO2*   15  13.771   1.735   5.374
  653ATP   AC3*   16  13.815   1.844   5.166
  653ATP   AO3*   17  13.789   1.710   5.127
  653ATP   AC5*   18  13.937   2.026   5.037
  653ATP   AO5*   19  13.900   2.130   5.129
  653ATP    APA   20  13.886   2.281   5.078
  653ATP   AO1A   21  14.013   2.323   5.013
  653ATP   AO2A   22  13.831   2.367   5.186
  653ATP   AO3A   23  13.777   2.263   4.964
  653ATP    APB   24  13.619   2.211   4.964
  653ATP   AO1B   25  13.613   2.060   4.988
  653ATP   AO2B   26  13.537   2.304   5.055
  653ATP   AO3B   27  13.580   2.240   4.805
  653ATP    APG   28  13.487   2.357   4.734
  653ATP   AO1G   29  13.563   2.488   4.757
  653ATP   AO2G   30  13.479   2.316   4.587
  653ATP   AO3G   31  13.354   2.349   4.810
   0.00000   0.00000   0.00000

Why is pdb2gmx failing for the 54a7 forcefield?

Thanks a lot for your time and consideration in advance!

Best
Julian

If you don’t have the hydrogen in your structure you need to get pdb2gmx to build it by using the -missing option and having an appropriate .hdb file entry.

More importantly though, is there a good reason to use this force field? There are both AMBER and CHARMM ATP parameters available with the more common protonation state available that will behave better than this old GROMOS model. Plus if you’re simulating something like a protein with the ATP this will likely be better with ones of these all-atom force field families too.

Cheers

Tom

Hi Tom,

thanks for your answer and pointing to the more appropriate force fields.
The reason why I am interested in GROMOS is that I want to do a benchmark study comparing different force fields, and I was interested how GROMOS performs.

Regarding the hydrogens:
In principle, I was wondering why pdb2gmx worked without issues for 53a6 but not for 54a7.
Thanks for mentioning the .hdb files, I checked them and I think I found now the root cause why the gamma phosphate hydrogen is not added:
While (some) atom names are changed in the .rtp definition from 53a6 to 54a7 (eg, AH3G → AH3PG), the atom names stay the same in the .hdb files (eg, AH3G → AH3G). I am quite new to GROMACS, but I would assume that this atom name mismatch causes the AH3PG to not be generated. As only some atom names are different in the .rtp files, this would also explain why some hydrogens are added as expected but the gamma phosphate one not.

Here the excerpts from the respective files:
…/top/gromos53a6.ff/aminoacids.rtp

...
[ ATP ]
 [ atoms ]
  AN9    NR    -0.20000     0
  AC4     C     0.20000     0
  AN3    NR    -0.36000     1
  AC2   CR1     0.36000     1
  AN1    NR    -0.36000     2
  AC6     C     0.36000     2
  AN6    NT    -0.83000     3
 AH61     H     0.41500     3
 AH62     H     0.41500     3
  AC5     C     0.00000     4
  AN7    NR    -0.36000     4
  AC8   CR1     0.36000     4
 AC1*   CH1     0.20000     5
 AO4*    OA    -0.36000     5
 AC4*   CH1     0.16000     5
 AC2*   CH1     0.15000     6
 AO2*    OA    -0.54800     6
 AH2*     H     0.39800     6
 AC3*   CH1     0.15000     7
 AO3*    OA    -0.54800     7
 AH3*     H     0.39800     7
 AC5*   CH2     0.00000     8
 AO5*    OA    -0.36000     9
  APA     P     0.70500     9
AO1A    OM    -0.63500     9
AO2A    OM    -0.63500     9
AO3A    OA    -0.36000    10
  APB     P     0.70500    10
AO1B    OM    -0.63500    10
AO2B    OM    -0.63500    10
AO3B    OA    -0.36000    11
  APG     P     0.63000    11
AO1G    OM    -0.63500    11
AO2G    OM    -0.63500    11
AO3G    OA    -0.54800    11
AH3G     H     0.39800    11
...

…/top/gromos54a7.ff/aminoacids.rtp:

...
[ ATP ]
 [ atoms ]
  AN9    NR    -0.20000     0
  AC4     C     0.20000     0
  AN3    NR    -0.36000     1
  AC2   CR1     0.36000     1
  AN1    NR    -0.36000     2
  AC6     C     0.36000     2
  AN6    NT    -0.83000     3
 AH61     H     0.41500     3
 AH62     H     0.41500     3
  AC5     C     0.00000     4
  AN7    NR    -0.36000     4
  AC8   CR1     0.36000     4
 AC1*   CH1     0.20000     5
 AO4*    OA    -0.36000     5
 AC4*   CH1     0.16000     5
 AC2*   CH1     0.15000     6
 AO2*    OA    -0.54800     6
 AH2*     H     0.39800     6
 AC3*   CH1     0.15000     7
 AO3*    OA    -0.54800     7
 AH3*     H     0.39800     7
 AC5*   CH2     0.00000     8
 AO5*    OA    -0.36000     9
  APA     P     0.70500     9
AO1PA    OM    -0.63500     9
AO2PA    OM    -0.63500     9
AO3PA    OA    -0.36000    10
  APB     P     0.70500    10
AO1PB    OM    -0.63500    10
AO2PB    OM    -0.63500    10
AO3PB    OA    -0.36000    11
  APG     P     0.63000    11
AO1PG    OM    -0.63500    11
AO2PG    OM    -0.63500    11
AO3PG    OA    -0.54800    11
AH3PG     H     0.39800    11
...

…/top/gromos53a6.ff/aminoacids.hdb:

...
ATP    4
2      3       AH6     AN6     AC6     AC5
1      2       AH2*    AO2*    AC2*    AC3*
1      2       AH3*    AO3*    AC3*    AC2*
1      2       AH3G    AO3G    APG     AO2G
...

…/top/gromos54a7.ff/aminoacids.hdb:

...
ATP    4
2      3       AH6     AN6     AC6     AC5
1      2       AH2*    AO2*    AC2*    AC3*
1      2       AH3*    AO3*    AC3*    AC2*
1      2       AH3G    AO3G    APG     AO2G
...

And indeed, when I edit the …/top/gromos54a7.ff/aminoacids.hdb to

...
ATP    4
2      3       AH6     AN6     AC6     AC5
1      2       AH2*    AO2*    AC2*    AC3*
1      2       AH3*    AO3*    AC3*    AC2*
1      2       AH3PG    AO3PG    APG     AO2PG
...

, the topology can be generated successfully.