[gmx-users] Problem in Free Energy Calculations
Justin Lemkul
jalemkul at vt.edu
Tue Jan 5 13:37:47 CET 2016
On 1/5/16 1:34 AM, shagun krishna wrote:
> Hi Justin..
> Thanks a ton to you. Finally I have completed the simulation of ligand in
> solvent. In your earlier mail you have written that to calculate the
> binding energy of protein-ligand I will need the nonbonded interaction
> energy between protein-ligand coming from protein-ligand simulation and the
> ligand-water interaction energy coming from the second simulation. For
> final ΔG binding I should substract these two energies, if I am not wrong.
> PLs suggest me.
>
No, you don't simply subtract them, but this approach requires LIE, which I
already said I don't use and requires some empirical parameters that you have to
figure out.
-Justin
> Thanx in advance..
>
> Best regards,
>
> SHAGUN
>
> On Mon, Jan 4, 2016 at 4:58 PM, Justin Lemkul <jalemkul at vt.edu> wrote:
>
>>
>>
>> On 1/4/16 6:25 AM, shagun krishna wrote:
>>
>>> Hii Justin,
>>>
>>> Thanks for the suggestion. When I am performing the ligand simulation in
>>> water I am getting the following error in nvt equilibration:
>>>
>>> Program mdrun, VERSION 4.6.5
>>> Source code file: /build/buildd/gromacs-4.6.5/src/mdlib/domdec_top.c,
>>> line:
>>> 393
>>>
>>> Fatal error:
>>> 56 of the 13118 bonded interactions could not be calculated because some
>>> atoms involved moved further apart than the multi-body cut-off distance
>>> (0.841804 nm) or the two-body cut-off distance (1.7 nm), see option -rdd,
>>> for pairs and tabulated bonds also see option -ddcheck
>>> For more information and tips for troubleshooting, please check the
>>> GROMACS
>>> website at http://www.gromacs.org/Documentation/Errors
>>>
>>>
>>
>> http://www.gromacs.org/Documentation/Terminology/Blowing_Up#Diagnosing_an_Unstable_System
>>
>> Your ligand topology is probably incorrect in some way. Turn off the free
>> energy code, use sensible cutoffs as required by the parent force field,
>> and follow the above link.
>>
>> -Justin
>>
>>
>> My nvt.mdp file is as follows:
>>>
>>> ; Run control
>>> integrator = sd ; Langevin dynamics
>>> tinit = 0
>>> dt = 0.002
>>> nsteps = 50000 ; 100 ps
>>> nstcomm = 100
>>> ; Output control
>>> nstxout = 500
>>> nstvout = 500
>>> nstfout = 0
>>> nstlog = 500
>>> nstenergy = 500
>>> nstxtcout = 0
>>> xtc-precision = 1000
>>> ; Neighborsearching and short-range nonbonded interactions
>>> nstlist = 10
>>> ns_type = grid
>>> pbc = xyz
>>> rlist = 1.7
>>> ; Electrostatics
>>> coulombtype = PME
>>> rcoulomb = 1.7
>>> ; van der Waals
>>> vdw-type = switch
>>> rvdw-switch = 0.8
>>> rvdw = 0.9
>>> ; Apply long range dispersion corrections for Energy and Pressure
>>> DispCorr = EnerPres
>>> ; Spacing for the PME/PPPM FFT grid
>>> fourierspacing = 0.12
>>> ; EWALD/PME/PPPM parameters
>>> pme_order = 6
>>> ewald_rtol = 1e-06
>>> epsilon_surface = 0
>>> optimize_fft = no
>>> ; Temperature coupling
>>> ; tcoupl is implicitly handled by the sd integrator
>>> tc_grps = system
>>> tau_t = 1.0
>>> ref_t = 300
>>> ; Pressure coupling is off for NVT
>>> Pcoupl = No
>>> tau_p = 0.5
>>> compressibility = 4.5e-05
>>> ref_p = 1.0
>>> ; Free energy control stuff
>>> free_energy = yes
>>> init_lambda = 0.0
>>> delta_lambda = 0
>>> foreign_lambda = 0.05
>>> sc-alpha = 0.5
>>> sc-power = 1.0
>>> sc-sigma = 0.3
>>> couple-moltype = 1M_R ; name of moleculetype to decouple
>>> couple-lambda0 = vdw ; only van der Waals interactions
>>> couple-lambda1 = none ; turn off everything, in this case
>>> only vdW
>>> couple-intramol = no
>>> nstdhdl = 10
>>> ; Generate velocities to start
>>> gen_vel = yes
>>> gen_temp = 300
>>> gen_seed = -1
>>> ; options for bonds
>>> constraints = h-bonds ; we only have C-H bonds here
>>> ; Type of constraint algorithm
>>> constraint-algorithm = lincs
>>> ; Do not constrain the starting configuration
>>> continuation = no
>>> ; Highest order in the expansion of the constraint coupling matrix
>>> lincs-order = 12
>>>
>>>
>>> On Mon, Jan 4, 2016 at 4:20 PM, Justin Lemkul <jalemkul at vt.edu> wrote:
>>>
>>>
>>>>
>>>> On 1/4/16 1:09 AM, shagun krishna wrote:
>>>>
>>>> Hi Justin,
>>>>> You mean to say that I don't need to follow your Free Energy
>>>>> Calculations
>>>>> :
>>>>> Methane in Water.
>>>>> Actually I want to calculate the binding energy of my ligands with my
>>>>> protein. All I should do is to run the one simulation of my ligands
>>>>> simply
>>>>> in water and another is the protein-ligand complex simulation, both for
>>>>> 10ns. I am really very confused because I have gone through literature
>>>>> and
>>>>> all of them have calculated the value of α, β and γ, which i guess can
>>>>> be
>>>>> calculated by using free energy codes.
>>>>>
>>>>>
>>>>> I don't use the LIE algorithm, so I don't know all the details. AFAIK
>>>> these are empirical parameters to weight various contributions to the
>>>> interaction energy. LIE is ultimately an estimate of the free energy,
>>>> and
>>>> its accuracy depends on the quality of the fitted parameters.
>>>>
>>>> A rigorous free energy calculation is the best method for calculating
>>>> free
>>>> energy of binding. For instance, LIE (again, AFAIK) relies solely on
>>>> enthalpic terms, so it does not account for entropy unless you do
>>>> additional work to figure this out separately.
>>>>
>>>>
>>>> -Justin
>>>>
>>>> --
>>>> ==================================================
>>>>
>>>> Justin A. Lemkul, Ph.D.
>>>> Ruth L. Kirschstein NRSA Postdoctoral Fellow
>>>>
>>>> Department of Pharmaceutical Sciences
>>>> School of Pharmacy
>>>> Health Sciences Facility II, Room 629
>>>> University of Maryland, Baltimore
>>>> 20 Penn St.
>>>> Baltimore, MD 21201
>>>>
>>>> jalemkul at outerbanks.umaryland.edu | (410) 706-7441
>>>> http://mackerell.umaryland.edu/~jalemkul
>>>>
>>>> ==================================================
>>>> --
>>>> 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.
>>>>
>>>>
>> --
>> ==================================================
>>
>> Justin A. Lemkul, Ph.D.
>> Ruth L. Kirschstein NRSA Postdoctoral Fellow
>>
>> Department of Pharmaceutical Sciences
>> School of Pharmacy
>> Health Sciences Facility II, Room 629
>> University of Maryland, Baltimore
>> 20 Penn St.
>> Baltimore, MD 21201
>>
>> jalemkul at outerbanks.umaryland.edu | (410) 706-7441
>> http://mackerell.umaryland.edu/~jalemkul
>>
>> ==================================================
>> --
>> 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.
>>
--
==================================================
Justin A. Lemkul, Ph.D.
Ruth L. Kirschstein NRSA Postdoctoral Fellow
Department of Pharmaceutical Sciences
School of Pharmacy
Health Sciences Facility II, Room 629
University of Maryland, Baltimore
20 Penn St.
Baltimore, MD 21201
jalemkul at outerbanks.umaryland.edu | (410) 706-7441
http://mackerell.umaryland.edu/~jalemkul
==================================================
More information about the gromacs.org_gmx-users
mailing list