Thanks for your feedback Ray.
I was already using radiopt=0 so I tried radiopt=1 but pbsa couldnt
define radii for the extra points. This could easily be fixed by
defining radii for EP and recompiling but I thought I should look at
the radii in the prmtop files that I was using for both the simulation
and the MM-PBSA calculations first. It turns out that the radii for
TIP5P in both parmtop files are incorrect with one hydrogen set to
0.8, one to 1.2 and the two EP have radii of 1.5 which is the same as
the oxygen.
Searching through the code I found that radii are defined in unitio.c
and it is here where a problems lie. The issue with the EP are that
they are not recognised so they are given the default value of 1.5.
This therefore effects both the simulation and MM-PBSA prmtop files. I
dont think this will effect the simulation as the radii are not used
in sander for PMEMD simulations if I understand things correctly.
However the hydrogen radii issues are different and I cant quite
figure out why we get different radii. The default hydrogen radii is
1.2 but if bonded to an Oxygen it should be 0.8. Therefore it seems
that only one hydrogen is being recognised as being bonded to oxygen
with maybe the other one being recognised as being bound to a hydrogen
as multiply bonded hydrogen used the first bond to define which radii
to use. The only problem with this theory is that for the TIP3P
simulation parmtop, the hydrogens are recognized correctly and given
the correct radii. However for the TIP3P MM-PBSA parmtop, where only
the structural water is included, the hydrogens are not recognized
properly and are defined radii as seen for the TIP5P parmtop. It is at
this point I am confused as to why the TIP3P radii are not recognized
properly for the structural water parmtop.
Given that I thought I had identified why I was getting highly
variable PBCAL values, I hand edited the parmtop files to define
appropriate radii (EP=0.0, hydrogens=0.8) and ran the MM-PBSA
calculations again. However, I still got the same increased variance
but this time, the GB energy could not be calculated (value of
-infinity was given in *all.out files).
Complex SD Receptor SD
Ligand SD Delta SD
PBCAL -2634.95 341.17 -2032.19 66.55
-45.77 1.23 -556.98 347.61
GB 0.00 0.00 0.00 0.00
-39.91 1.26 39.91 1.26
I then tried Ray's suggestion of using the default hydrogen radius of
1.3A (even though in unitio.c the default is 1.2).
Complex SD Receptor SD
Ligand SD Delta SD
PBCAL -2206.90 193.60 -2126.48 97.35
-45.77 1.23 -34.64 216.70
So defining the radii "correctly" didnt help but in the end has posed
more question. Why are the radii defined incorrectly in some parmtop's
and is there any way of reducing the variability in the PBCAl values
in the MM-PBSA calculations?
Any help with these issues would be greatly appreciated.
Daniel Oehme
PhD Student
La Trobe University
On Tue, Dec 2, 2008 at 12:55 PM, Ray Luo <rayhuangluo.gmail.com> wrote:
> Daniel,
>
> Just my 2c ...
>
> If you use radiopt=1, I think the radius optimized for TIP3P may be too
> small for TIP5P, i.e. the two extra points are too close to the surface of
> the solvent probe sphere boundary, which may cause highly fluctuating PB
> energies.
>
> Maybe use of radiopt=0 could resolve the problem since the radius for the
> two hydrogen atoms (1.30A) may be large enough to cover the two EP's.
>
> All the best,
> Ray
>
> ==========================================
> Ray Luo, Ph.D.
> Associate Professor
> Dept Molecular Biology & Biochemistry
> University of California, Irvine, CA 92697
> USPS: PO Box 3900 Email: rluo.uci.edu
> Phones: (949) 824-9528, 9562
> Web: http://rayl0.bio.uci.edu/
> ==========================================
>
>
> -----Original Message-----
> From: owner-amber.scripps.edu [mailto:owner-amber.scripps.edu] On Behalf Of
> Daniel Oehme
> Sent: Monday, December 01, 2008 3:45 PM
> To: amber
> Subject: AMBER: Variable PBCAL values in MM-PBSA calculations
>
> Hi,
>
> I am running some MM_PBSA calculations using AMBER 10 and have been
> getting some highly variable PBCAL values depending on what water type
> I assign a structural water molecule in the active site to be. The
> issue comes when trying to use the TIP5P water type. For TIP3P
> simulations I get the following results:
> Complex SD Receptor SD
> Ligand SD Delta SD
> PBCAL -2066.21 39.76 -2103.88 39.57
> -45.30 1.39 82.96 7.24
>
> while for the TIP5P I get:
> Complex SD Receptor SD
> Ligand SD Delta SD
> PBCAL -2549.08 369.70 -2188.62 295.53
> -45.77 1.23 314.69 243.66
>
> Please note that two different simulations were run to get these
> results. As can be seen, the Complex and Receptor have extremely high
> SD's but the Ligands values are fairly similar. If however I run the
> MM-PBSA calculations without the structural water molecules in the
> TIP5P simulation I get the following results:
> Complex SD Receptor SD
> Ligand SD Delta SD
> PBCAL -1982.99 66.10 -2024.60 66.13
> -45.77 1.23 87.38 4.61
>
> I have also run a calculation on the TIP5P simulation by making the
> structural water molecule TIP3P instead of TIP5P and got these
> results:
> Complex SD Receptor SD
> Ligand SD Delta SD
> PBCAL -1982.75 66.07 -2027.20 65.99
> -45.77 1.23 90.23 4.65
>
> As can be seen, the last two sets of results give SD's are a lot lower
> in this case yet still a bit higher than in the inital TIP3P
> simulations. My theory is that the PBCAL calculations are getting
> confused by the extra points (EP) that are included in the TIP5P water
> model compared to TIP3P. My question is can anybody think of a reason
> as to why the PBCAL values would get confused with the EP's and how I
> can change the setup of the calculation to be able to run the MM-PBSA
> calculations with a TIP5P structural water molecule?
>
> MM-PBSA calculations were run without a solvent sphere and with
> standard input except for fillratio which was increased to 4.
>
> Daniel Oehme
> PhD Student
> La Trobe University
> -----------------------------------------------------------------------
> The AMBER Mail Reflector
> To post, send mail to amber.scripps.edu
> To unsubscribe, send "unsubscribe amber" (in the *body* of the email)
> to majordomo.scripps.edu
>
> -----------------------------------------------------------------------
> The AMBER Mail Reflector
> To post, send mail to amber.scripps.edu
> To unsubscribe, send "unsubscribe amber" (in the *body* of the email)
> to majordomo.scripps.edu
>
-----------------------------------------------------------------------
The AMBER Mail Reflector
To post, send mail to amber.scripps.edu
To unsubscribe, send "unsubscribe amber" (in the *body* of the email)
to majordomo.scripps.edu
Received on Fri Dec 05 2008 - 18:40:43 PST
