[gmx-users] pdb2gmx picks up the wrong .tdb files?

[Justin Lemkul]

On 1/16/20 8:04 AM, András Ferenc WACHA wrote:
> Dear fellow Gromacs users,
>
> I have developed an extended version of the CHARMM36m force field for
> beta-amino acids (https://charmm-betaff.readthedocs.io,
> https://gitlab.com/awacha/charmm-beta.ff). For beta-amino acids I found
> pdb2gmx works well. However, if I try to use it for natural peptides and
> proteins (e.g. 6LYZ from the PDB), the following command:
>
> gmx pdb2gmx -f 6LYZ.pdb, -ff charmm-beta -ter -ignh
>
> fails to recognize the correct termini database files. The output is the
> following:
>
> (-------------------------- START of gmx output
> -----------------------------------)
>
> GROMACS:      gmx pdb2gmx, version 2020
> Executable:
> /opt/spack/opt/spack/linux-archrolling-skylake/gcc-9.2.0/gromacs-2020-ls53yrft6droj2ovhuzdwdufi4an4vir/bin/gmx
> Data prefix:
> /opt/spack/opt/spack/linux-archrolling-skylake/gcc-9.2.0/gromacs-2020-ls53yrft6droj2ovhuzdwdufi4an4vir
> Working dir:  /home/wachaandras/gromacs/forcefields
> Command line:
>    gmx pdb2gmx -f 6lyz.pdb -ff charmm-beta -ter -ignh
>
> Using the Charmm-beta force field in directory ./charmm-beta.ff
>
> Opening force field file ./charmm-beta.ff/watermodels.dat
>
> Select the Water Model:
>   1: TIP3P    TIP 3-point, recommended, by default uses CHARMM TIP3 with
> LJ on H
>   2: TIP4P    TIP 4-point
>   3: TIP5P    TIP 5-point
>   4: SPC        simple point charge
>   5: SPC/E    extended simple point charge
>   6: Methanol   Equilibrated methanol
>   7: Octanol    Equilibrated octanol
>   8: None
> 1
> going to rename ./charmm-beta.ff/merged.r2b
> Opening force field file ./charmm-beta.ff/merged.r2b
> Reading 6lyz.pdb...
> WARNING: all CONECT records are ignored
> Read '', 1001 atoms
> Analyzing pdb file
> Splitting chemical chains based on TER records or chain id changing.
> There are 1 chains and 0 blocks of water and 129 residues with 1001 atoms
>
>    chain  #res #atoms
>    1 'A'   129   1001
>
> All occupancies are one
> Opening force field file ./charmm-beta.ff/atomtypes.atp
> Opening force field file ./charmm-beta.ff/betaaminoacids.atp
> Opening force field file ./charmm-beta.ff/cyclicbeta.atp
> Reading residue database... (Charmm-beta)
> Opening force field file ./charmm-beta.ff/betaaminoacids.rtp
> Opening force field file ./charmm-beta.ff/cyclicbeta.rtp
> Opening force field file ./charmm-beta.ff/extra.rtp
> Warning: file does not end with a newline, last line:
> ;
> Opening force field file ./charmm-beta.ff/merged.rtp
> Opening force field file ./charmm-beta.ff/merged.hdb
> Opening force field file ./charmm-beta.ff/betaaminoacids.n.tdb
> Opening force field file ./charmm-beta.ff/cyclicbeta.n.tdb
> Opening force field file ./charmm-beta.ff/merged.n.tdb
> Opening force field file ./charmm-beta.ff/betaaminoacids.c.tdb
> Opening force field file ./charmm-beta.ff/cyclicbeta.c.tdb
> Opening force field file ./charmm-beta.ff/merged.c.tdb
> Processing chain 1 'A' (1001 atoms, 129 residues)
> Analysing hydrogen-bonding network for automated assignment of histidine
>   protonation. 213 donors and 184 acceptors were found.
> There are 262 hydrogen bonds
> Will use HISD for residue 15
> Identified residue LYS1 as a starting terminus.
> Identified residue LEU129 as a ending terminus.
> 8 out of 8 lines of specbond.dat converted successfully
> Special Atom Distance matrix:
>                      CYS6   MET12   HIS15   CYS30   CYS64   CYS76   CYS80
>                      SG48    SD87  NE2118   SG238   SG513   SG601   SG630
>     MET12    SD87   1.205
>     HIS15  NE2118   1.804   0.998
>     CYS30   SG238   1.450   1.070   2.063
>     CYS64   SG513   2.873   1.779   1.755   2.236
>     CYS76   SG601   2.740   1.544   1.418   2.132   0.778
>     CYS80   SG630   3.006   1.943   1.892   2.387   0.200   0.958
>     CYS94   SG724   2.576   1.388   1.329   1.976   0.678   0.211   0.871
>    MET105   SD799   1.866   0.933   1.679   0.896   1.858   1.483   2.047
>    CYS115   SG889   1.625   1.099   2.068   0.204   2.108   2.005   2.260
>    CYS127   SG981   0.203   1.089   1.733   1.332   2.808   2.631   2.950
>                     CYS94  MET105  CYS115
>                     SG724   SD799   SG889
>    MET105   SD799   1.389
>    CYS115   SG889   1.855   0.790
>    CYS127   SG981   2.474   1.701   1.505
> Linking CYS-6 SG-48 and CYS-127 SG-981...
> Linking CYS-30 SG-238 and CYS-115 SG-889...
> Linking CYS-64 SG-513 and CYS-80 SG-630...
> Linking CYS-76 SG-601 and CYS-94 SG-724...
> Select start terminus type for LYS-1
>   0: Beta3NH3+
>   1: Beta2NH3+
>   2: Beta23NH3+
>   3: B3_NH2
>   4: B2_NH2
>   5: B23_NH2
>   6: 5TER
>   7: None
>
> (-------------------------- END of gmx output
> -----------------------------------)
>
> Information in the force field relevant to the beta-peptides is stored
> in betaaminoacids.{arn,atp,c.tdb,n.tdb,rtp} files. I have kept the
> original ones (merged.*) intact. From the output above it seems that
> pdb2gmx recognizes and loads all the relevant files (both
> betaaminoacids.* and merged.*), but only shows a selection of N-terminal
> patches from the betaaminoacids.n.tdb file. More precisely, the 5TER
> terminus is from merged.n.tdb, but the other ones before it get
> overridden by entries from betaaminoacids.n.tdb. Is there some rule in
> pdb2gmx which considers some endgroups equivalent and keeps only one? I
> have encountered it neither in the main documentation (section 5.6.5:
> pdb2gmx input files), nor in the help of gmx pdb2gmx.
>
> Could anyone point me to how I should do this correctly?

What is the source of these force field files? Based on the .tdb 
entries, it seems they do not include the "standard" patches that should 
be included with the official CHARMM force field.

-Justin

-- 
==================================================

Justin A. Lemkul, Ph.D.
Assistant Professor
Office: 301 Fralin Hall
Lab: 303 Engel Hall

Virginia Tech Department of Biochemistry
340 West Campus Dr.
Blacksburg, VA 24061

jalemkul at vt.edu | (540) 231-3129
http://www.thelemkullab.com

==================================================