AMBER Archive (2009)

Subject: Re: [AMBER] Using idecomp=3 with igb=10 (PB)

From: Chris Whittleston (csw34_at_cam.ac.uk)
Date: Wed Oct 07 2009 - 11:57:25 CDT

Dear all,

Thanks for all your help on this already - I have just returned to this
problem now that I have an AMBER10 license, and am still having trouble
understanding the group specification. I have decomposed the energy two
ways, using idecomp=1 and idecomp=3 just to make sure they give the same
answer when I look at a single residue, in this case 159. The protein runs
from residue 1->329 and the ligand is residue 330.

Here is my idecomp=3 SANDER input file:

Interaction energy per residue input for SANDER
 &cntrl
  imin = 1,
  idecomp = 3,
  ncyc = 1,
  maxcyc = 0,
  igb = 1, saltcon=0.1,
  ntb = 0,
  cut = 999.0,
  rgbmax = 8.22,
 /
Protein
RRES 1 329
END
Ligand
LRES 330 330
END
Printing
RES 1 330
END
END

As expected, this gives a pairwise per-residue decomposition - so to look
for interactions to residue 330, I just use:

grep TDC sander.out | grep '\-> 330'

Here is a sample of the output from grep:
...
TDC 157-> 330 0.000 -0.061 0.892 -0.839 0.000
TDC 158-> 330 0.000 -0.024 -0.100 0.091 0.000
TDC 159-> 330 0.000 -1.146 -25.181 19.673 0.000
TDC 160-> 330 0.000 -1.158 -0.396 0.434 0.000
TDC 161-> 330 0.000 -0.087 -0.020 0.037 0.000
TDC 162-> 330 0.000 -0.008 0.218 -0.216 0.000
...

So - adding up the interaction energy from 159->330, you get -6.654kcal/mol

Now - as Hannes pointed out, I really want to be using idecomp=1, especially
when using PB as I waste a lot of time calculating interactions I just throw
away. Unfortunately, so far I have been unable to reproduce this value for
residue 159. Here is the idecomp=1 input I'm using to test this:

Interaction energy per residue input for SANDER
 &cntrl
  imin = 1,
  idecomp = 1,
  ncyc = 1,
  maxcyc = 0,
  igb = 1, saltcon=0.1,
  ntb = 0,
  cut = 999.0,
  rgbmax = 8.22,
 /
Protein
RRES 159 159
END
Ligand
LRES 330 330
END
Printing
RES 159
END
END

Here, I'm only printing for residue 159 and I get:

                    PRINT DECOMP - TOTAL ENERGIES

    resid |internal |vdw |eel |pol |sas
============================================================
TDC 159 -238.487 -2.889 130.226 -35.665 0.000

This is clearly a different total interaction! What am I missing here? Is it
something to do with how you specify groups as I've tried all sorts of
things (switching LRES<->RRES, LRES and RRES ->RES etc) and nothing can
produce that same total interaction between residue 159 and 330.

What is causing this difference?

Any help figuring this out is much appriciated!

Chris

2009/7/28 Ashish Runthala <ashish.runthala_at_gmail.com>

> Hello Chris,
>
> See this output, you will realize what is the problem out there.
>
>
> NSTEP = 0 TIME(PS) = 0.000 TEMP(K) = 0.00 PRESS = 0.0
> Etot = ************ EKtot = 0.0000 EPtot = ************
> BOND = 27342.0949 ANGLE = 978.1773 DIHED = 1354.4338
> 1-4 NB = 5277.1875 1-4 EEL = 7145.0848 VDWAALS = ************
> EELEC = -9360.5582 EHBOND = 0.0000 RESTRAINT = 0.0000
>
> ------------------------------------------------------------------------------
>
> vlimit exceeded for step 0; vmax = 24186095.7
>
> NSTEP = 1 TIME(PS) = 0.001 TEMP(K) = 0.00 PRESS = 0.0
> Etot = ************ EKtot = 0.0000 EPtot = ************
> BOND = 27342.0949 ANGLE = 978.1773 DIHED = 1354.4338
> 1-4 NB = 5277.1875 1-4 EEL = 7145.0848 VDWAALS = ************
> EELEC = -9360.5582 EHBOND = 0.0000 RESTRAINT = 0.0000
>
> ------------------------------------------------------------------------------
>
>
> A V E R A G E S O V E R 1 S T E P S
>
>
> NSTEP = 1 TIME(PS) = 0.001 TEMP(K) = 0.00 PRESS = 0.0
> Etot = ************ EKtot = 0.0000 EPtot = ************
> BOND = 27342.0949 ANGLE = 978.1773 DIHED = 1354.4338
> 1-4 NB = 5277.1875 1-4 EEL = 7145.0848 VDWAALS = ************
> EELEC = -9360.5582 EHBOND = 0.0000 RESTRAINT = 0.0000
>
> ------------------------------------------------------------------------------
>
>
> R M S F L U C T U A T I O N S
>
>
> NSTEP = 1 TIME(PS) = 0.001 TEMP(K) = 0.00 PRESS = 0.0
> Etot = 0.0000 EKtot = 0.0000 EPtot = 0.0000
> BOND = 0.0002 ANGLE = 0.0000 DIHED = 0.0000
> 1-4 NB = 0.0000 1-4 EEL = 0.0000 VDWAALS = 0.0000
> EELEC = 0.0001 EHBOND = 0.0000 RESTRAINT = 0.0000
>
> ------------------------------------------------------------------------------
>
> Total Energy is coming very high. So this system is totally unstable.
> You do one thing. This i checked to check the confirmation.
>
> You do one thing.Try minimizing the energy of the system. And then
> work your script out. Make the system with stable with lowest feasible
> energy.
> Ashish
>
>
> On Tue, Jul 28, 2009 at 6:17 PM, Hannes Kopitz<Hannes.Kopitz_at_gmx.de>
> wrote:
> > Hi Chris,
> >
> > Sorry, my mistake. It should be something like ...
> >
> > ...
> > RRES 1 318
> > END
> > LRES 319 319
> > END
> > RES 1 319
> > END
> > END
> >
> > The RRES and LRES cards define the residues to be considered for
> decomposition, whereas the RES card defines the residues for the output.
> > Moreover you should chose idecomp=1 for a per-residue decomposition. With
> idecomp=3 you would get a pairwise per-residue decomposition.
> > Also set ntmin to 2 and maxcyc to 0, cause you don't really want to
> perform any minimization step.
> > But anyway, for doing a PB decomposition you need AMBR10.
> > I hope that helps.
> >
> > Cheers!
> > Hannes
> >
> > -------- Original-Nachricht --------
> >> Datum: Tue, 28 Jul 2009 11:54:09 +0100
> >> Von: Chris Whittleston <csw34_at_cam.ac.uk>
> >> An: AMBER Mailing List <amber_at_ambermd.org>
> >> Betreff: Re: [AMBER] Using idecomp=3 with igb=10 (PB)
> >
> >> Dear Hannes,
> >>
> >> Thanks for the quick reply! I've tried changing to LRES and RRES and it
> >> seems to just suppress the output:
> >>
> >> FINAL RESULTS
> >>
> >>
> >>
> >> NSTEP ENERGY RMS GMAX NAME
> NUMBER
> >> 1 -1.1108E+04 1.3245E+00 1.4038E+01 HG 2021
> >>
> >> BOND = 191.2619 ANGLE = 696.4471 DIHED =
> >> 2808.3488
> >> VDWAALS = -2871.4868 EEL = -23739.7079 EPB =
> >> -3329.3383
> >> 1-4 VDW = 980.7815 1-4 EEL = 14155.9277 RESTRAINT =
> >> 0.0000
> >> ECAVITY = 0.0000 EDISPER = 0.0000
> >>
> >>
> >> CHECK DECOMP - TOTAL ENERGIES (w/ REST)
> >>
> >> INTERNAL= 0.0000
> >> VDWAALS = 0.0000 EEL = 0.0000
> >> EGB = 0.0000 ESURF = 0.0000
> >>
> >>
> >> CHECK DECOMP - SELF ENERGIES (w/o REST)
> >>
> >> INTERNAL= 0.0000
> >> VDWAALS = 0.0000 EEL = 0.0000
> >> EGB = 0.0000 ESURF = 0.0000
> >>
> >>
> >> CHECK DECOMP - INDIRECT ENERGIES (w/o REST)
> >>
> >> INTERNAL= 0.0000
> >> VDWAALS = 0.0000 EEL = 0.0000
> >> EGB = 0.0000 ESURF = 0.0000
> >>
> >>
> >> CHECK DECOMP - DIRECT ENERGIES (w/o REST)
> >>
> >> INTERNAL= 0.0000
> >> VDWAALS = 0.0000 EEL = 0.0000
> >> EGB = 0.0000 ESURF = 0.0000
> >>
> >>
> >> CHECK DECOMP - REST ENERGIES
> >>
> >> INTERNAL= 0.0000
> >> VDWAALS = 0.0000 EEL = 0.0000
> >> EGB = 0.0000 ESURF = 0.0000
> >>
> >>
> >>
> >>
> >> PRINT PAIR DECOMP - TOTAL ENERGIES
> >>
> >> resid1 ->resid2 |internal |vdw |eel |pol |sas
> >> ======================================================================
> >>
> >>
> >> PRINT PAIR DECOMP - SIDECHAIN ENERGIES
> >>
> >> resid1 ->resid2 |internal |vdw |eel |pol |sas
> >> ======================================================================
> >>
> >>
> >> PRINT PAIR DECOMP - BACKBONE ENERGIES
> >>
> >> resid1 ->resid2 |internal |vdw |eel |pol |sas
> >> ======================================================================
> >>
> >> It seems that the values for the internal, electrostatic and vdw
> energies
> >> are actually non-zero in the total output, just for the decomposition so
> I
> >> guess it could be that because the PB breakdown is not implemented in
> >> AMBER9
> >> - I'm just getting zeroes. Thanks for the help, I might try switching to
> >> AMBER10 at some point soon :)
> >>
> >> Best,
> >>
> >> Chris
> >>
> >> 2009/7/28 Hannes Kopitz <Hannes.Kopitz_at_gmx.de>
> >>
> >> > Dear Chris,
> >> >
> >> > A decomposition of the PB reaction field energy is implemented from
> >> AMBER10
> >> > but unfortunately not in AMBER9.
> >> >
> >> > Nevertheless the values for internal, vdw and eel shouldn't be zero.
> >> > Just a guess, you should try it with RRES and LRES instead of RES.
> >> >
> >> > Cheers!
> >> > Hannes
> >> >
> >> > -------- Original-Nachricht --------
> >> > > Datum: Mon, 27 Jul 2009 21:04:14 +0100
> >> > > Von: Chris Whittleston <csw34_at_cam.ac.uk>
> >> > > An: amber_at_ambermd.org
> >> > > Betreff: [AMBER] Using idecomp=3 with igb=10 (PB)
> >> >
> >> > > Dear AMBER users,
> >> > >
> >> > > I'm currently looking at decomposed per-residue interaction energies
> >> in a
> >> > > protein-ligand system with a GB solvent model (e.g. igb=1/5) but
> would
> >> > > like
> >> > > to be able to use Poisson Boltzmann (igb=10) to compare with some
> >> results
> >> > > from an MMPBSA calculation of the same system. Unfortunately, I am
> >> just
> >> > > getting a huge table of 0.000 for every interaction energy.
> >> > >
> >> > > Here is my SANDER input:
> >> > >
> >> > > Interaction energy per residue input for SANDER
> >> > > &cntrl
> >> > > imin = 1,
> >> > > idecomp = 3,
> >> > > ncyc = 1,
> >> > > maxcyc = 1,
> >> > > igb = 10,
> >> > > ntb = 0,
> >> > > cut = 0, <- have also tried 999.9
> >> > > /
> >> > > &pb
> >> > > npopt = 0,
> >> > > radiopt = 0
> >> > > /
> >> > > First set
> >> > > RES -1 318
> >> > > END
> >> > > Second set
> >> > > RES 319
> >> > > END
> >> > > END
> >> > >
> >> > > I am using modified GLYCAM parameters for a sugar ligand, hence have
> >> set
> >> > > npopt=0 and radiopt=0. SANDER seems to run ok (it doesn't produce
> any
> >> > > warnings or crash out), but after taking 20 seconds or so
> calculating
> >> > > things
> >> > > - it produces a table as follows:
> >> > >
> >> > >
> >> > > PRINT PAIR DECOMP - TOTAL ENERGIES
> >> > >
> >> > > resid1 ->resid2 |internal |vdw |eel |pol |sas
> >> > >
> ======================================================================
> >> > > TDC 1-> 1 0.000 0.000 0.000 0.000
> 0.000
> >> > > TDC 1-> 2 0.000 0.000 0.000 0.000
> 0.000
> >> > > TDC 1-> 3 0.000 0.000 0.000 0.000
> 0.000
> >> > > TDC 1-> 4 0.000 0.000 0.000 0.000
> 0.000
> >> > > TDC 1-> 5 0.000 0.000 0.000 0.000
> 0.000
> >> > > TDC 1-> 6 0.000 0.000 0.000 0.000
> 0.000
> >> > > TDC 1-> 7 0.000 0.000 0.000 0.000
> 0.000
> >> > > TDC 1-> 8 0.000 0.000 0.000 0.000
> 0.000
> >> > > TDC 1-> 9 0.000 0.000 0.000 0.000
> 0.000
> >> > > TDC 1-> 10 0.000 0.000 0.000 0.000
> 0.000
> >> > > TDC 1-> 11 0.000 0.000 0.000 0.000
> 0.000
> >> > > TDC 1-> 12 0.000 0.000 0.000 0.000
> 0.000
> >> > > TDC 1-> 13 0.000 0.000 0.000 0.000
> 0.000
> >> > > TDC 1-> 14 0.000 0.000 0.000 0.000
> 0.000
> >> > > TDC 1-> 15 0.000 0.000 0.000 0.000
> 0.000
> >> > > TDC 1-> 16 0.000 0.000 0.000 0.000
> 0.000
> >> > > TDC 1-> 17 0.000 0.000 0.000 0.000
> 0.000
> >> > > TDC 1-> 18 0.000 0.000 0.000 0.000
> 0.000
> >> > > TDC 1-> 19 0.000 0.000 0.000 0.000
> 0.000
> >> > > TDC 1-> 20 0.000 0.000 0.000 0.000
> 0.000
> >> > > TDC 1-> 21 0.000 0.000 0.000 0.000
> 0.000
> >> > > TDC 1-> 22 0.000 0.000 0.000 0.000
> 0.000
> >> > > TDC 1-> 23 0.000 0.000 0.000 0.000
> 0.000
> >> > > TDC 1-> 24 0.000 0.000 0.000 0.000
> 0.000
> >> > > TDC 1-> 25 0.000 0.000 0.000 0.000
> 0.000
> >> > > TDC 1-> 26 0.000 0.000 0.000 0.000
> 0.000
> >> > > TDC 1-> 27 0.000 0.000 0.000 0.000
> 0.000
> >> > > ...
> >> > > ...
> >> > > etc (they're all zero)
> >> > >
> >> > > The total energy however is not zero as can be seen in the output
> >> above
> >> > > the
> >> > > interaction energy table:
> >> > >
> >> > > NSTEP ENERGY RMS GMAX NAME
> >> NUMBER
> >> > > 1 -1.1108E+04 1.3245E+00 1.4038E+01 HG
> >> 2021
> >> > >
> >> > > BOND = 191.2619 ANGLE = 696.4471 DIHED =
> >> > > 2808.3488
> >> > > VDWAALS = -2871.4868 EEL = -23739.7079 EPB =
> >> > > -3329.3383
> >> > > 1-4 VDW = 980.7815 1-4 EEL = 14155.9277 RESTRAINT =
> >> > > 0.0000
> >> > > ECAVITY = 0.0000 EDISPER = 0.0000
> >> > >
> >> > >
> >> > > Maximum number of minimization cycles reached.
> >> > >
> >> > >
> >> > > FINAL RESULTS
> >> > >
> >> > >
> >> > >
> >> > > NSTEP ENERGY RMS GMAX NAME
> >> NUMBER
> >> > > 1 -1.1108E+04 1.3245E+00 1.4038E+01 HG
> >> 2021
> >> > >
> >> > > BOND = 191.2619 ANGLE = 696.4471 DIHED =
> >> > > 2808.3488
> >> > > VDWAALS = -2871.4868 EEL = -23739.7079 EPB =
> >> > > -3329.3383
> >> > > 1-4 VDW = 980.7815 1-4 EEL = 14155.9277 RESTRAINT =
> >> > > 0.0000
> >> > > ECAVITY = 0.0000 EDISPER = 0.0000
> >> > >
> >> > >
> >> > > CHECK DECOMP - TOTAL ENERGIES (w/ REST)
> >> > >
> >> > > INTERNAL= 0.0000
> >> > > VDWAALS = 0.0000 EEL = 0.0000
> >> > > EGB = 0.0000 ESURF = 0.0000
> >> > >
> >> > >
> >> > > CHECK DECOMP - SELF ENERGIES (w/o REST)
> >> > >
> >> > > INTERNAL= 0.0000
> >> > > VDWAALS = 0.0000 EEL = 0.0000
> >> > > EGB = 0.0000 ESURF = 0.0000
> >> > >
> >> > >
> >> > > CHECK DECOMP - INDIRECT ENERGIES (w/o REST)
> >> > >
> >> > > INTERNAL= 0.0000
> >> > > VDWAALS = 0.0000 EEL = 0.0000
> >> > > EGB = 0.0000 ESURF = 0.0000
> >> > >
> >> > >
> >> > > CHECK DECOMP - DIRECT ENERGIES (w/o REST)
> >> > >
> >> > > INTERNAL= 0.0000
> >> > > VDWAALS = 0.0000 EEL = 0.0000
> >> > > EGB = 0.0000 ESURF = 0.0000
> >> > >
> >> > >
> >> > > CHECK DECOMP - REST ENERGIES
> >> > >
> >> > > INTERNAL= 0.0000
> >> > > VDWAALS = 0.0000 EEL = 0.0000
> >> > > EGB = 0.0000 ESURF = 0.0000
> >> > >
> >> > >
> >> > > So - I was wondering if anyone could help me understand why I'm just
> >> > > getting
> >> > > zeroes! My one idea so far is that it is not possible to break down
> >> the
> >> > PB
> >> > > surface energy on a per-atom or per-residue basis and so SANDER is
> >> > > printing
> >> > > zeros. I'm using AMBER9 with all current bug fixes.
> >> > >
> >> > > Any and all help is greatly appriciated!
> >> > >
> >> > > Chris
> >> > >
> >> > > --
> >> > > Chris Whittleston
> >> > > Department of Chemistry
> >> > > University of Cambridge
> >> > > Lensfield Road, Cambridge, CB2 1EW
> >> > > Email: csw34_at_cam.ac.uk
> >> > > Tel: +44 (0)1223 336423
> >> > > _______________________________________________
> >> > > AMBER mailing list
> >> > > AMBER_at_ambermd.org
> >> > > http://lists.ambermd.org/mailman/listinfo/amber
> >> >
> >> > --
> >> > Neu: GMX Doppel-FLAT mit Internet-Flatrate + Telefon-Flatrate
> >> > für nur 19,99 Euro/mtl.!* http://portal.gmx.net/de/go/dsl02
> >> >
> >> > _______________________________________________
> >> > AMBER mailing list
> >> > AMBER_at_ambermd.org
> >> > http://lists.ambermd.org/mailman/listinfo/amber
> >> >
> >>
> >>
> >>
> >> --
> >> Chris Whittleston
> >> Department of Chemistry
> >> University of Cambridge
> >> Lensfield Road, Cambridge, CB2 1EW
> >> Email: csw34_at_cam.ac.uk
> >> Tel: +44 (0)1223 336423
> >> _______________________________________________
> >> AMBER mailing list
> >> AMBER_at_ambermd.org
> >> http://lists.ambermd.org/mailman/listinfo/amber
> >
> > --
> > Jetzt kostenlos herunterladen: Internet Explorer 8 und Mozilla Firefox 3
> -
> > sicherer, schneller und einfacher! http://portal.gmx.net/de/go/chbrowser
> >
> > _______________________________________________
> > AMBER mailing list
> > AMBER_at_ambermd.org
> > http://lists.ambermd.org/mailman/listinfo/amber
> >
>
>
>
> --
> Ashish Runthala,
> Faculty Division III,
> Lecturer, Biological Sciences,
> Birla Institute of Technology and Science,
> Pilani, Rajasthan- 333031
> INDIA
>
> _______________________________________________
> AMBER mailing list
> AMBER_at_ambermd.org
> http://lists.ambermd.org/mailman/listinfo/amber
> 

-- 
Chris Whittleston
Department of Chemistry
University of Cambridge
Lensfield Road, Cambridge, CB2 1EW
Email: csw34_at_cam.ac.uk
Tel: +44 (0)1223 336423
_______________________________________________
AMBER mailing list
AMBER_at_ambermd.org
http://lists.ambermd.org/mailman/listinfo/amber