Skip to content

[BUG?] Incompatible AMBER atom names for protein cap residues #365

@jmichel80

Description

@jmichel80

I have come across something annoying with AMBER templates.
 
So if you export a PDB from Flare version 8.0.2 in 'AMBER PDB' mode it writes ACE and NME residues
with atom names consistent with those found in
 

$AMBERHOME/dat/leap/lib/aminont10.lib and 
$AMBERHOME/dat/leap/lib/aminoct10.lib

 
In these lib files we have for ACE
 

!entry.ACE.unit.atoms table  str name  str type  int typex  int resx  int flags  int seq  int elmnt  dbl chg
 "HH31" "HC" 0 1 131072 1 1 0.112300
 "CH3" "CT" 0 1 131072 2 6 -0.366200
 "HH32" "HC" 0 1 131072 3 1 0.112300
 "HH33" "HC" 0 1 131072 4 1 0.112300
 "C" "C" 0 1 131072 5 6 0.597200

 for NME

!entry.NME.unit.atoms table  str name  str type  int typex  int resx  int flags  int seq  int elmnt  dbl chg
 "N" "N" 0 1 131072 1 7 -0.415700
 "H" "H" 0 1 131072 2 1 0.271900
 "CH3" "CT" 0 1 131072 3 6 -0.149000
 "HH31" "H1" 0 1 131072 4 1 0.097600
 "HH32" "H1" 0 1 131072 5 1 0.097600
 "HH33" "H1" 0 1 131072 6 1 0.097600

 
 
However ambertools (here 23.6 from conda-forge) now seems to load lib files for all their forcefields from
 

$AMBERHOME/dat/leap/lib/aminont12.lib and 
$AMBERHOME/dat/leap/lib/aminoct12.lib

 
In these entries ACE is specified as
 

!entry.ACE.unit.atoms table  str name  str type  int typex  int resx  int flags  int seq  int elmnt  dbl chg
 "H1" "HC" 0 1 131072 1 1 0.112300
 "CH3" "CT" 0 1 131072 2 6 -0.366200
 "H2" "HC" 0 1 131072 3 1 0.112300
 "H3" "HC" 0 1 131072 4 1 0.112300
 "C" "C" 0 1 131072 5 6 0.597200
 "O" "O" 0 1 131072 6 8 -0.567900
(…)

 
And NME as
 

!entry.NME.unit.atoms table  str name  str type  int typex  int resx  int flags  int seq  int elmnt  dbl chg
 "N" "N" 0 1 131072 1 7 -0.415700
 "H" "H" 0 1 131072 2 1 0.271900
 "C" "CT" 0 1 131072 3 6 -0.149000
 "H1" "H1" 0 1 131072 4 1 0.097600
 "H2" "H1" 0 1 131072 5 1 0.097600
 "H3" "H1" 0 1 131072 6 1 0.097600
(…)

 
As a result  trying to parameterise a protein exported from Flare
 
protein = BSS.Parameters.ff14SB(protein, work_dir='tmp').getMolecule()
 
fails with this leap error message
 

FATAL:  Atom .R<ACE 1>.A<HH32 7> does not have a type.
FATAL:  Atom .R<ACE 1>.A<HH33 8> does not have a type.
FATAL:  Atom .R<ACE 1>.A<HH31 9> does not have a type.
FATAL:  Atom .R<NGLU 3>.A<H 18> does not have a type.
FATAL:  Atom .R<NME 164>.A<CH3 7> does not have a type.
FATAL:  Atom .R<NME 164>.A<HH31 8> does not have a type.
FATAL:  Atom .R<NME 164>.A<HH32 9> does not have a type.
FATAL:  Atom .R<NME 164>.A<HH33 10> does not have a type.

 

I wonder if there is a sensible way to automatically deal with inconsistent amber naming schemes behind the scenes when calling
 
protein_xtal = BSS.IO.readPDB('inputs/3STD_prep.pdb', pdb4amber=True)[0]
 
Note how pdb4amber doesn’t seem to catch this issue.

To reproduce the issue

3STD.zip

import BioSimSpace as BSS

protein_xtal = BSS.IO.readPDB('inputs/3STD_prep.pdb', pdb4amber=True)[0]

protein = protein_xtal.extract(
   [atom.index() for atom in protein_xtal.search("not resname HOH").atoms()]
)

protein = BSS.Parameters.ff14SB(protein, work_dir='tmp').getMolecule()

cat tmp/leap.log

tagging @mark-mackey-cresset as this is perhaps more a Flare bug than a BioSimSpace bug

Metadata

Metadata

Assignees

No one assigned

    Labels

    bugSomething isn't working

    Type

    No type
    No fields configured for issues without a type.

    Projects

    No projects

    Milestone

    No milestone

    Relationships

    None yet

    Development

    No branches or pull requests

    Issue actions