Did you flip the sign of the reference energy? That would make your
equivalence point off by dG/(kTpKa) pH units, which would probably be
quite a few.
Brian
On 10/04/2015 09:01 PM, Jason Swails wrote:
> On Sun, Oct 4, 2015 at 5:42 PM, Rahul Ramesh <raramesh.umich.edu> wrote:
>
>> Hello Jason
>> I'm trying to simulate a cellulosic molecule with a succinic group. ( This
>> is my model compound as I am going to simulate a larger system in future )
>> I assigned the charges using the AMBER 99SB FF. I performed Thermodynamic
>> Integration of this structure and obtained a value of -25 kcal/mol. Since I
>> wanted to confirm the value of this del G ref against the pKa , I performed
>> a CpHMD at pH = pKa. I expected the system to be protonated/ deprotonated
>> 50 % of the time. But I end up getting 100% protonation. I used the same
>> set of parameters for both TI and CpHMD.
>> I went through the mdout file and the EEL , EGB and 1-4 EEL are not
>> convincing.
>>
> If you are getting 100% protonation, then your reference energies must be
> wrong.
>
>
>> NSTEP = 125000 TIME(PS) = 2250.000 TEMP(K) = 237.22 PRESS =
>> 0.0
>>
>> Etot = -262.6830 EKtot = 24.2773 *EPtot =
>> -286.9603*
>>
>> BOND = 21.5394 ANGLE = 26.4520 DIHED =
>> 24.9527
>>
>> 1-4 NB = 6.7956 *1-4 EEL = -89.1130* VDWAALS =
>> -3.9488
>>
>> *EELEC = -25.1198 EGB = -248.5186* RESTRAINT =
>> 0.0000
>>
>> The electrostatic energy is positive if you subtract 1-4 EEL from EELEC.
>>
> Why are you subtracting them? The total electrostatic energy is the *sum*
> of the two.
>
>
>
>> I'm thinking of a reason as to why the electrostatic energy is really low
>> and 1-4 EEL and EGB is really high?
>>
>> I compared my system (cellulosic molecule attached to a succinic group) to
>> a succinic acid (COOHCH2CH2COOH) system. Basically I simulated just the
>> side group of the cellulosic molecule using the AMBER99SB FF. I followed
>> the same procedure as above and I was able to get 50 % protonation at pH =
>> pKa and a del G ref as -22 kcal/mol ( pH = 4.2).
>>
>> I tried comparing the mdout file for succinic acid to the above mdout file
>>
>> *Protonated State for Succinic Acid*
>>
>> NSTEP = 999600 TIME(PS) = 3999.200 TEMP(K) = 545.62 PRESS =
>> 0.0
>>
>> Etot = -63.8361 EKtot = 19.5165 *EPtot =
>> -83.3526*
>>
>> BOND = 2.1652 ANGLE = 8.1213 DIHED =
>> 7.6514
>>
>> 1-4 NB = 0.8964 *1-4 EEL = 19.4734* VDWAALS =
>> -0.8920
>>
>> *EELEC = -93.5325 EGB = -27.2358* RESTRAINT =
>> 0.0000
>> *Deprotonated State for Succinic Acid*
>>
>> NSTEP = 999550 TIME(PS) = 3999.100 TEMP(K) = 438.95 PRESS =
>> 0.0
>>
>> Etot = -88.1985 EKtot = 15.7011 * EPtot =
>> -103.8997*
>>
>> BOND = 2.2672 ANGLE = 10.6449 DIHED =
>> 8.7225
>>
>> 1-4 NB = 0.6721 *1-4 EEL = 21.5060 * VDWAALS =
>> -1.0602
>>
>> *EELEC = -55.1840 EGB = -91.4680 *RESTRAINT =
>> 0.0000
>> For the Succinic Acid case , EELEC is negative and different between the
>> EPtot of the protonated and deprotonated is approximately the Del G ref. I
>> was expecting the same behavior for my cellulosic molecule attached to
>> succinic group.
>> Could you comment on the reason for EELEC being positive in the first case
>> and high value of 1-4 EEL and EGB
>>
> EELEC can be positive, and if you have a system that carries a net charge,
> you should expect that the EELEC *would* be positive (like charges repel).
> EGB screens ELEC, so it is often around the same size, but opposite in
> sign, of EELEC. As for the 1-4 EEL interactions -- those can be very large
> if dihedrals connect atoms with large charges.
>
> The difference in electrostatic energy between two different charge
> vectors can be surprising and counterintuitive -- fixed-charge force fields
> are not *intended* to have changing charges, so the free energy of changing
> charges is unpredictable. That's why we have the reference energy -- to
> correct for everything the force field gets "wrong".
>
> Another, potentially easier, way of calculating reference energies is to
> calculate a single-point energy difference for a typical conformation
> between the two charge states and use that as a starting guess. Titrate
> the model compound and then use the pKa you fit to the titration curve to
> correct the reference energy (this can be done in a single step).
>
> HTH,
> Jason
>
--
Brian Radak
Postdoctoral Scholar
Gordon Center for Integrative Science, W323A
Department of Biochemistry & Molecular Biology
University of Chicago
929 E. 57th St.
Chicago, IL 60637-1454
Tel: 773/834-2812
email: radak.uchicago.edu
_______________________________________________
AMBER mailing list
AMBER.ambermd.org
http://lists.ambermd.org/mailman/listinfo/amber
Received on Mon Oct 05 2015 - 07:00:04 PDT
