On Wed, Jul 27, 2011 at 8:26 AM, Peter Varnai <P.Varnai.sussex.ac.uk> wrote:
> Thanks, Jason for the response.
>
> >> I prepare a protein prmtop file using tleap with current bug-fixed
> AmberTools15.
> >> I obtain different results using the commands below:
> >
> > What are "results"?
>
> I obtain prmtop files that differ, but I had not tested the corresponding
> energy values by sander. I realise that the CT to CX change must have effect
> on the format of the prmtop sections. An energy calculation does show indeed
> identical results thus no change in the underlying force fields. What was
> the reason for the CT to CX type change for CA?
>
Yes, the prmtop files differ. As you've found out, this is not at all
indicative of a different, underlying 'force field'. Here's actually an
interesting exercise to try with tleap:
source leaprc.ff99SB
peptide = sequence {NALA ALA CALA}
peptide2 = copy peptide
saveamberparm peptide parm1 crd1
saveamberparm peptide2 parm2 crd2
Then diff parm1 and parm2 -- you will get what appear to be vastly different
topologies (although they give identical energies).
As for why they changed CT to CX, this is just as I understand it -- ff99SB
took the Amber99 force field and fixed some of its most glaring issues with
backbone dihedrals by implementing a couple torsion corrections (see
frcmod.ff99SB in $AMBERHOME/dat/leap/parm). An unexpected consequence of
this, however, was that several of the backbone torsion corrections actually
ended up affecting some sidechain torsion angles as well, since CT was a
common atom type for carbon atoms. Thus, by retyping the alpha carbon to CX
instead of CT, you can *specifically* target backbone dihedrals without fear
of affecting any sidechain dihedrals at the same time. That is, it adds
greater potential flexibility into force field definitions.
>
> >> source leaprc.ff99SBildn
> >
> > This is the "proper" way to load this force field.
>
> section 2.2.1 of the manual describes this differently, leaprc.ff99SBildn
> is not mentioned.
>
Perhaps it hasn't been documented yet. However, if a leaprc for a specific
force field exists, I would argue that it was authored in such a way that it
loads all of the files that were *meant* to be loaded exactly. Thus, the
leaprc.ff99SBildn should mimic the *proper* way of doing it as defined in
the manual. Since I didn't actually have any hand in
contributing/generating this force field, I don't consider myself the
definitive source, but if this *isn't* the case, it certainly doesn't follow
Amber's typical conventions.
>
>
> >> give almost the same prmtop except the modified Bondi radii - I wonder
> where these radii are stored for leap and why the new ildn atom types
> changed those.
> >
> > These are stored in the leap source code. I forget where off the top of
> my head, but it's easily enough greppable. Assignments are done based on
> atom type and bonded neighbors (see the source code for details).
> >
> > I've written a program that will (among other prmtop modifications) allow
> you to adjust the radii to whichever set you want. The logic leap uses to
> assign radii can be found in that source as well, which I can share with you
> if you want (it may be easier to read the logic in python with separate
> cases for each radii set rather than in C all jumbled together).
>
> probably the different mbondi assigned by leap are due to some perceived
> difference in the bonded neighbours depending on which way the prmtop was
> created.
>
It compares typing. See line 5268 of unitio.c in
$AMBERHOME/AmberTools/src/leap/src/leap/ for the section that defines these
radii. The hydrogen atom has the most complex rules for defining radii
(only H radii depend on bonded atoms), but it looks like you're looking for
Carbon rules. This is what I see:
5315 case 6:
5316 if ( strncmp(sType,"C1",2) &&
strncmp(sType,"C2",2) && strncmp (sType,"C3",2) )
5317 dGBrad = 1.7;
5318 else
5319 dGBrad = 2.2;
5320 break;
That is, C1*, C2*, and C3* (where * is a wildcard that can be any (or no)
character(s)) are given radii of 1.7, and all other carbons are given radii
of 2.2. Thus, if one FF assigns a type C1 and another assigns a type C4/C5,
then this control block will apply different radii. I'm not sure how
serious this is, since GB is rarely used (successfully) on nucleic acids
anyway, and it's mainly nucleic acids that use atom types C1, C2, C3, etc.
(I think).
HTH,
Jason
--
Jason M. Swails
Quantum Theory Project,
University of Florida
Ph.D. Candidate
352-392-4032
_______________________________________________
AMBER mailing list
AMBER.ambermd.org
http://lists.ambermd.org/mailman/listinfo/amber
Received on Wed Jul 27 2011 - 13:00:17 PDT
