GROMACS version: 2021.2
GROMACS modification: Yes: from Homebrew on mac
I have a problem when dealing with n terminal capping and/or c terminal capping for protein. I am using the charmm force field downloaded from MacKerell Lab
Here I have my pdb file with the capping added manually:
ATOM 1 CAY ACE U 0 -16.453 -12.961 -24.554 1.00 0.00 U4 C
ATOM 2 HY2 ACE U 0 -15.896 -12.727 -23.654 1.00 0.00 U4 H
ATOM 3 HY3 ACE U 0 -16.253 -13.969 -24.896 1.00 0.00 U4 H
ATOM 4 HY1 ACE U 0 -17.514 -12.782 -24.425 1.00 0.00 U4 H
ATOM 5 CY ACE U 0 -15.963 -12.008 -25.648 1.00 0.00 U4 C
ATOM 6 OY ACE U 0 -16.436 -12.058 -26.782 1.00 0.00 U4 O
ATOM 7 N NDC U 1 -15.004 -11.123 -25.308 1.00 0.00 U4 N
ATOM 8 CA NDC U 1 -14.418 -10.138 -26.223 1.00 0.00 U4 C
ATOM 9 HA2 NDC U 1 -13.998 -10.620 -26.992 1.00 0.00 U4 H
(More lines not shown, I have my user-defined NDC residue which works perfectly.)
When I run pdb2gmx, the ACE capping has the error:
"Identified residue ACE0 as a starting terminus.
(some line omitted here)
Start terminus ACE-0: NH3+ "
So I guess the ACE in the merged.rtp is matched with NH3+ in the merged.n.tdb.
and more lines of error(s):
" Fatal error:
atom N not found in buiding block 1ACE while combining tdb and rtp"
I have found that merged.rtp has the ACE record as a residue, but merged.*.tdb has no such record, so I guess pdb2gmx has some matching rules that ACE is matched with the closest capping entry, i.e. NH3+, to be combined with. Since the default NH3+ in the merged.n.tdb file has the line of N atom, which is to replace the N in the first residue, pdb2gmx returned the error in search of N.
[ NH3+ ]
[ replace ]
N N NH3 14.0070 -0.30
CA CA CT1 12.011 0.21
HA HA HB1 1.008 0.10
[ Add ]
3 4 H N CA C
HC 1.008 0.33 -1
[ delete ]
To resolve this, I renamed the ACE in pdb into ACT and create the ACT residue in the rtp file with the same structure as of ACE in merged.rtp. Then, I created an empty entry of [ACTOD] in tdb file, and pdb2gmx returned:
Start terminus ACT-0: ACTOD
which indicates that pdb2gmx has matched the ACT in rtp file with the ACTOD in n.tdb file and combined the two files to generate the topology.
Moreover, I tried pdb2gmx with the option -ter in order to assign the termini manually to ‘None’, since I had the capping included in the pdb file. It gave me:
"Start terminus ACT-0: None "
Therefore, I believe that when manually choose “None”, it will not search for the matching entry in *.tdb files, and thus avoid the confusing problem.
The problem seems resolved, but I hope someone can give the guidelines in developing the force field files. How pdb2gmx works with the capping? How the searching and matching works? Please advise about the priorities. Thank you in advance.
P.S. The default ACE entry in merged.rtp has the very unique atom types as
[CT3]-[HA3]*3 — C=O — ,
where the HA3-CT3-C=O dihedral parameters are missing. I guess the ACE is actually incomplete and hope you can fix this in the future. (I have used different atom types to avoid this problem.)