[gmx-users] Large differences in calculated and experimental relative hydration free energy.
Abhishek Acharya
abhi117acharya at gmail.com
Sat Aug 11 09:10:38 CEST 2018
Hi,
No, I have not tried it for a simpler molecule. I thought initially to do
something like ethane to ethanol which would be appropriate I guess. But
decide to mail to the forum to know if I was missing something obvious. The
protocol is based on my understanding of a few papers from Mobley group and
uses the setup scripts provided by them.
But as you also suggest, I should try this for a simpler system and see if
it works.
Thank you Mark.
On Sat 11 Aug, 2018, 05:02 Mark Abraham, <mark.j.abraham at gmail.com> wrote:
> Hi,
>
> Sorry for not having time to dig in detail, but does your protocol work for
> a simpler molecule?
>
> Mark
>
> On Fri, Aug 10, 2018, 13:33 Abhishek Acharya <abhi117acharya at gmail.com>
> wrote:
>
> > Dear Gromacs Users,
> >
> > I am trying to calculate the relative hydration free energy of Catechol
> wrt
> > Benzene as a way of validating the parameters of Catechol, which I need
> for
> > further free energy calculations. I am using OPLS-AA for these
> calculations
> > and the parameters for both Catechol and Benzene were generated using the
> > LigParGen webserver. Subsequently, I used the alchemical_setup.py script
> (
> > https://github.com/MobleyLab/alchemical-setup) to generate the topology
> > for
> > transformation of catechol to benzene. This was followed by the necessary
> > steps to generate a solvated system.
> >
> > Relative hydration free energy was calculated in 3 steps:
> >
> > a) Discharging catechol using the following topology setup done (using 11
> > lambda values):
> >
> > ; nr type resnr res atom cgnr charge mass
> > typeB chargeB massB comments MCSS = Minimum Common Sub-Structure
> > 1 opls_800 1 LIG C1 1 -0.1361 12.0110
> > opls_800 0.0 12.0110 ; MCSS
> > 2 opls_801 1 LIG C2 2 -0.1361 12.0110
> > opls_801 0.0 12.0110 ; MCSS
> > 3 opls_802 1 LIG C3 3 -0.1693 12.0110
> > opls_802 0.0 12.0110 ; MCSS
> > 4 opls_803 1 LIG C4 4 0.2517 12.0110
> > opls_803 0.0 12.0110 ; MCSS
> > 5 opls_804 1 LIG C5 5 0.3382 12.0110
> > opls_804 0.0 12.0110 ; MCSS
> > 6 opls_805 1 LIG C6 6 -0.1693 12.0110
> > opls_805 0.0 12.0110 ; MCSS
> > 7 opls_810 1 LIG HA1 7 0.1508 1.0080
> > opls_808 0.0 1.0080 ; MCSS
> > 8 opls_810_dummy 1 LIG HA2 8 0.00000 1.0080
> > opls_810 0.0 1.0080 ; to be appeared
> > 9 opls_811 1 LIG HB1 9 0.1508 1.0080
> > opls_811 0.0 1.0080 ; MCSS
> > 10 opls_812 1 LIG HC1 10 0.4461 1.0080
> > opls_812_dummy 0.00000 1.0080 ; to be annihilated
> > 11 opls_813 1 LIG HD1 11 0.4422 1.0080
> > opls_813_dummy 0.00000 1.0080 ; to be annihilated
> > 12 opls_806 1 LIG O1 12 -0.7156 15.9990
> > opls_806_dummy 0.00000 15.9990 ; to be annihilated
> > 13 opls_807 1 LIG O2 13 -0.742 15.9990
> > opls_807_dummy 0.00000 15.9990 ; to be annihilated
> > 14 opls_808 1 LIG H1 14 0.1443 1.0080
> > opls_808 0.0 1.0080 ; MCSS
> > 15 opls_809 1 LIG H2 15 0.1443 1.0080
> > opls_809 0.0 1.0080 ; MCSS
> > 16 opls_809_dummy 1 LIG H3 16 0.00000 1.0080
> > opls_809 0.0 1.0080 ; to be appeared
> >
> > Here only coulomb-lambdas were varied as: 0.00 0.10 0.20 0.30 0.40 0.50
> > 0.60 0.70 0.80 0.90 1.00
> >
> > b) LJ-transformation from catechol to benzene (26 lambda values):
> >
> > ; nr type resnr res atom cgnr charge mass
> > typeB chargeB massB comments MCSS = Minimum Common Sub-Structure
> > 1 opls_800 1 LIG C1 1 0.0 12.0110
> > opls_800 0.0 12.0110 ; MCSS
> > 2 opls_801 1 LIG C2 2 0.0 12.0110
> > opls_801 0.0 12.0110 ; MCSS
> > 3 opls_802 1 LIG C3 3 0.0 12.0110
> > opls_802 0.0 12.0110 ; MCSS
> > 4 opls_803 1 LIG C4 4 0.0 12.0110
> > opls_803 0.0 12.0110 ; MCSS
> > 5 opls_804 1 LIG C5 5 0.0 12.0110
> > opls_804 0.0 12.0110 ; MCSS
> > 6 opls_805 1 LIG C6 6 0.0 12.0110
> > opls_805 0.0 12.0110 ; MCSS
> > 7 opls_808 1 LIG HA1 7 0.0 1.0080
> > opls_808 0.0 1.0080 ; MCSS
> > 8 opls_810_dummy 1 LIG HA2 8 0.0 1.0080
> > opls_810 0.0 1.0080 ; to be appeared
> > 9 opls_811 1 LIG HB1 9 0.0 1.0080
> > opls_811 0.0 1.0080 ; MCSS
> > 10 opls_812 1 LIG HC1 10 0.0 1.0080
> > opls_812_dummy 0.0 1.0080 ; to be annihilated
> > 11 opls_813 1 LIG HD1 11 0.0 1.0080
> > opls_813_dummy 0.0 1.0080 ; to be annihilated
> > 12 opls_806 1 LIG O1 12 0.0 15.9990
> > opls_806_dummy 0.0 15.9990 ; to be annihilated
> > 13 opls_807 1 LIG O2 13 0.0 15.9990
> > opls_807_dummy 0.0 15.9990 ; to be annihilated
> > 14 opls_808 1 LIG H1 14 0.0 1.0080
> > opls_808 0.0 1.0080 ; MCSS
> > 15 opls_809 1 LIG H2 15 0.0 1.0080
> > opls_809 0.0 1.0080 ; MCSS
> > 16 opls_809_dummy 1 LIG H3 16 0.0 1.0080
> > opls_809 0.0 1.0080 ; to be appeared
> >
> > Here only vdw-lambdas were varied as:
> > vdw_lambdas = 0.00 0.04 0.08 0.12 0.16 0.20 0.24 0.28 0.32
> > 0.36 0.40 0.44 0.48 0.52 0.56 0.60 0.64 0.68 0.72 0.76 0.80 0.84 0.88
> 0.92
> > 0.96 1.00
> >
> > Meanwhile, coulomb_lambdas were kept 0 throughout, as the charges of
> State
> > A and B have been set to zero in this topology.
> >
> > c) Recharging to full benzene charges (11 lambda values):
> >
> > ; nr type resnr res atom cgnr charge mass
> > typeB chargeB massB comments MCSS = Minimum Common Sub-Structure
> > 1 opls_800 1 LIG C1 1 0.0 12.0110
> > opls_800 -0.1474 12.0110 ; MCSS
> > 2 opls_801 1 LIG C2 2 0.0 12.0110
> > opls_801 -0.1476 12.0110 ; MCSS
> > 3 opls_802 1 LIG C3 3 0.0 12.0110
> > opls_802 -0.1474 12.0110 ; MCSS
> > 4 opls_803 1 LIG C4 4 0.0 12.0110
> > opls_803 -0.1472 12.0110 ; MCSS
> > 5 opls_804 1 LIG C5 5 0.0 12.0110
> > opls_804 -0.1470 12.0110 ; MCSS
> > 6 opls_805 1 LIG C6 6 0.0 12.0110
> > opls_805 -0.1470 12.0110 ; MCSS
> > 7 opls_808 1 LIG HA1 7 0.0 1.0080
> > opls_808 0.1474 1.0080 ; MCSS
> > 8 opls_810 1 LIG HA2 8 0.0 1.0080
> > opls_810 0.1470 1.0080 ; to be appeared
> > 9 opls_811 1 LIG HB1 9 0.0 1.0080
> > opls_811 0.1470 1.0080 ; MCSS
> > 10 opls_812_dummy 1 LIG HC1 10 0.0 1.0080
> > opls_812_dummy 0.0 1.0080 ; to be annihilated
> > 11 opls_813_dummy 1 LIG HD1 11 0.0 1.0080
> > opls_813_dummy 0.0 1.0080 ; to be annihilated
> > 12 opls_806_dummy 1 LIG O1 12 0.0 15.9990
> > opls_806_dummy 0.0 15.9990 ; to be annihilated
> > 13 opls_807_dummy 1 LIG O2 13 0.0 15.9990
> > opls_807_dummy 0.0 15.9990 ; to be annihilated
> > 14 opls_806 1 LIG H1 14 0.0 1.0080
> > opls_806 0.1474 1.0080 ; MCSS
> > 15 opls_807 1 LIG H2 15 0.0 1.0080
> > opls_807 0.1474 1.0080 ; MCSS
> > 16 opls_809 1 LIG H3 16 0.0 1.0080
> > opls_809 0.1474 1.0080 ; to be appeared
> >
> > In this case the topology was changed a bit to have the same LJ
> parameters
> > for both State A and B = Benzene.
> > Lambda values used for this transformation were:
> > coul_lambdas = 0.00 0.10 0.20 0.30 0.40 0.50 0.60 0.70 0.80
> > 0.90 1.00
> > vdw_lambdas = 1.00 1.00 1.00 1.00 1.00 1.00 1.00 1.00 1.00
> > 1.00 1.00 (Although, even a value of 0 throughout should work here)
> >
> > Finally, the free energy change for the these transformations were
> > calculated using the alchemical_analysis script (
> > https://github.com/MobleyLab/alchemical-analysis).
> > Result for the three steps [MBAR output only (in kcal/mol); statistical
> > error in brackets]: Discharging: 27.449 (0.028), LJ-transform: - 0.192
> > (0.001) and Recharging: - 0.348 (0.039)
> >
> > The problem is the the reported experimental value for relative hydration
> > free energy of catechol and benzene is around 8.5 kcal/mol which is very
> > different from the calculated value (27.6 kcal/mol). The simulations have
> > been performed using a similar (cutoff, electrostatics etc.) setup as
> > mentioned in the LigParGen paper. The default soft-core parameters were
> > used for all three steps.
> >
> > I also tried a one-step transformation (both vdw and coulomb lambda
> varied
> > simultaneously) and got a value of around 29 kcal/mol.
> >
> > Since, I am attempting a transmutation of this sort for the first time, I
> > wanted to be sure the setup is correct before I contact LigParGen
> authors.
> > Since this seems like a simple transformation, I expected differences
> > (between experimental and calculated values) of the order of 1.5-2
> kcal/mol
> > as seen in case of hydration free energies of test cases used in the
> > LigParGen paper. It would be nice to get some pointers as to what might
> be
> > going wrong with my setup.
> >
> > Thanks in advance.
> >
> > Abhishek Acharya
> > --
> > Gromacs Users mailing list
> >
> > * Please search the archive at
> > http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before
> > posting!
> >
> > * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
> >
> > * For (un)subscribe requests visit
> > https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or
> > send a mail to gmx-users-request at gromacs.org.
> >
> --
> Gromacs Users mailing list
>
> * Please search the archive at
> http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before
> posting!
>
> * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
>
> * For (un)subscribe requests visit
> https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or
> send a mail to gmx-users-request at gromacs.org.
>
More information about the gromacs.org_gmx-users
mailing list