[gmx-users] Does gmx covar/gmx anaeig give <dS> or T<dS> for ligand binding?
Billy Williams-Noonan
billy.williams-noonan at monash.edu
Wed Jun 22 08:22:40 CEST 2016
Yes, still with two ligands... So I assume that <dS> should be halved to
about -0.25 kJ/mol K, which would give T<dS>= -75.
I found out recently that we have access to a node with 1TB of RAM.
So the solvent is still relevant? If the amount of memory I need for the
entropy calculation is given by *( 3*N )^2 * 4 *where N is the number of
atoms, I should be able to calculate the entropy of systems (including
solvent) on that node. If I have to multiply by the number of frames then
I can't.
Is the formula above correct? Because according to it, with about 50,000
atoms in my system, I would need 90 GB of RAM, which is odd, because GMX
covar crashed when I used a node with 128GB RAM and one core.
Billy
On 22 June 2016 at 16:05, David van der Spoel <spoel at xray.bmc.uu.se> wrote:
> On 22/06/16 06:44, Billy Williams-Noonan wrote:
>
>> I re-did the calculation. When considering the entire biomolecule of
>> each ensemble:
>>
>> <S(P.L)>' = 51760 J/mol K
>>
>> <S(P)> = 50640 J/mol K
>>
>> <S(L)> = 814 J/mol K
>>
>> Resulting in <dS>= -0.51 kJ/mol K
>>
> Still with two ligands? This still corresponds to a binding entropy change
> of -150 kJ/mol. Of course you are ignoring the entropy change of the water
> which is probably almost the same magnitude and with opposite sign. If you
> want more quantitative results you could consider doing a PMF but your
> ligand is very large so that will be difficult to converge as well.
>
>
>
>>
>> And when just considering Protein-H, I got:
>>
>> <S(P.L)>' = 31941 J/mol K
>>
>> <S(P)> = 31340 J/mol K
>>
>> <S(L)> = 640 J/mol K
>>
>> Resulting in <dS>= -0.69 kJ/mol K
>>
>>
>> These values make more sense given my enthalpy calculation with
>> g_mmpbsa
>> is likely not converged. Thank you for your time and patience. :)
>>
>> Billy
>>
>>
>>
>>
>>
>>
>>
>>
>>
>> On 22 June 2016 at 13:40, Billy Williams-Noonan <
>> billy.williams-noonan at monash.edu> wrote:
>>
>> Sorry that was the ATB, not the ATP
>>>
>>> On 22 June 2016 at 13:39, Billy Williams-Noonan <
>>> billy.williams-noonan at monash.edu> wrote:
>>>
>>> Hi David,
>>>>
>>>> Thanks again for responding... Sorry if I came across the wrong way.
>>>> I'm not trying to disprove the code, but simply understand why my values
>>>> don't make sense I trust your knowledge on this subject too, since I
>>>> suspect you're one of the geniuses who helped to develop
>>>> g_covar/g_anaeig.
>>>> :) I know the units for entropy too.
>>>>
>>>> I should explain that I have previously performed relative FEP
>>>> calculations of ligands binding to the site of interest, and reproduced
>>>> experimental binding affinities within 1.4 kcal/mol of experiment.
>>>> Ligand
>>>> topologies came from the ATP using a GROMOS united atom force-field.
>>>> So I
>>>> know that the protocol I use for system equilibration is working.
>>>>
>>>> Using the same equilibration protocol as with the FEP protocol, and
>>>> having tried an absolute FEP calculation with restraints that failed
>>>> dismally, I have a cyclic peptide that has mM affinity for the same
>>>> binding
>>>> site as the aforementioned ligands (see above paragraph). So, using the
>>>> same protein as a model and placing the cyclic peptide in the correct
>>>> orientation as determined by the crystal structure, I am trying to use
>>>> g_mmpbsa to get an absolute binding affinity. Of course the entropic
>>>> term
>>>> from the g_mmpbsa calculation is missing, so I am using g_covar and
>>>> g_anaeig to determine the entropy.
>>>>
>>>> You're right about the size of my ligand too of course. The cyclic
>>>> peptide is 54 atoms in size and moves quite a lot in solution. I am
>>>> used a
>>>> Parrinello-Rahman/V-rescale NPT ensemble, set to 300K and 1 bar, for the
>>>> entirety of the 100ns simulation. And my protein is a symmetrical dimer
>>>> (two of the same protein bound to each other) so there is one ligand for
>>>> each monomer, forming 108 atoms between the two ligands.
>>>>
>>>> When I initially made this thread, the variables I was talking about
>>>> were:
>>>>
>>>> <S(P.L)> = 128,886 J/mol/K
>>>> = Entropy of one ligand bound to one side of the protein dimer,
>>>> despite another ligand being bound on the other side. a_1-3071 was
>>>> selected (twice) in g covar to represent the P.L complex, despite there
>>>> being 3125 atoms in total
>>>>
>>>> <S(P)> = 153,548 J/mol/K
>>>> = Entropy of the protein
>>>>
>>>> <S(L)> = 4137 J/mol/K
>>>> = Entropy of the cyclic peptide
>>>>
>>>> So I redefined <S(P.L)>, as <S(P.L)>', and selected a_1-3125 this
>>>> morning, to get the entropy of the dimer complexed with two cyclic
>>>> peptides, and got a value of 51,759.8 J/mol/K. I substituted this into
>>>> the
>>>> equation (1)
>>>>
>>>> <dS> = <S(P.L)>' - <S(P)> - 2*<S(L)> --------------(1)
>>>>
>>>> I multiplied the entropy of the ligand by two to account for the fact
>>>> that the beginning state now has the two ligands and the protein in
>>>> solution, while the end state has the protein dimer complexed with those
>>>> two ligands. And, the answer was -110.072 kJ/mol/K.
>>>>
>>>> So I am clearly doing something wrong and I'd like some advice on
>>>> what
>>>> it is... I doubt it's an equilibration problem, since my FEP
>>>> calculations
>>>> previously worked with the same equilibration protocol. And I doubt
>>>> this
>>>> is a convergence issue too, since a <dS> this high should prohibit
>>>> binding
>>>> in most cases, and my ligand definitely binds as seen by viewing the
>>>> complex simulation on VMD.
>>>>
>>>> Advice? Thoughts? I am about to try it with just the Protein-H atoms
>>>> from the index files to see if that changes anything...
>>>>
>>>> Billy
>>>>
>>>> On 21 June 2016 at 21:51, David van der Spoel <spoel at xray.bmc.uu.se>
>>>> wrote:
>>>>
>>>> On 21/06/16 11:26, Billy Williams-Noonan wrote:
>>>>>
>>>>> Hi Gromacs Users,
>>>>>>
>>>>>> I have used gmx covar and gmx anaeig to generate three ensemble
>>>>>> average
>>>>>> entropies over 100ns: first for a ligand in solution (<S(L)>), second
>>>>>> for a
>>>>>> protein in solution (<S(P)>) and third for their respective complex in
>>>>>> solution (<S(P.L)>).
>>>>>>
>>>>>> My understanding is that the change in entropy upon binding is
>>>>>> given
>>>>>> by:
>>>>>>
>>>>>> <dS> = <S(P.L)> - <S(P)> - <S(L)> -----------(1)
>>>>>>
>>>>>> Using gmx covar/gmx anaeig I got Quasi-Harmonic entropy estimates
>>>>>> of:
>>>>>>
>>>>>> <S(P.L)> = 128,886 J/mol/K
>>>>>>
>>>>>> <S(P)> = 153,548 J/mol/K
>>>>>>
>>>>>> <S(L)> = 4137 J/mol/K
>>>>>>
>>>>>> As stated, these values were generated using gmx covar/anaeig by
>>>>>> selecting for the relevant biomolecule in each ensemble and ignoring
>>>>>> the
>>>>>> effect of solvent movement.
>>>>>>
>>>>>> The unit printed by the program is J/mol K, which is the normal unit
>>>>> for
>>>>> entropy in all handbooks. You can not prove or disprove the
>>>>> correctness of
>>>>> the code by an example, you will have to check the code yourself if you
>>>>> doubt it.
>>>>>
>>>>> Looking at your numbers, they are huge and the difference is huge too.
>>>>> You should probably make sure first that all you simulations are in
>>>>> equilibrium. A ligand entropy och 4137 I would expect for an organic
>>>>> molecules with close to 100 carbon atoms.
>>>>>
>>>>>
>>>>>
>>>>> By subbing the above-described values into (1), I got about -28
>>>>>> kJ/mol/K
>>>>>> for <dS>, which is the right answer if the units are actually kJ/mol,
>>>>>> and
>>>>>> not kJ/mol/K. Strangely, upon multiplying by T, I got a value of
>>>>>> -8640
>>>>>> kJ/mol, which is quite obviously wrong.
>>>>>>
>>>>>> So does (1) yield a value for <S> or T<dS> ? Is anyone able to
>>>>>> explain
>>>>>> this to me?
>>>>>>
>>>>>> Kind regards,
>>>>>>
>>>>>> Billy
>>>>>>
>>>>>>
>>>>>>
>>>>>>
>>>>> --
>>>>> David van der Spoel, Ph.D., Professor of Biology
>>>>> Dept. of Cell & Molec. Biol., Uppsala University.
>>>>> Box 596, 75124 Uppsala, Sweden. Phone: +46184714205.
>>>>> spoel at xray.bmc.uu.se http://folding.bmc.uu.se
>>>>> --
>>>>> 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.
>>>>>
>>>>>
>>>>
>>>>
>>>> --
>>>> Billy Noonan* | *PhD Student *|* Bsci ( *Adv* ), IA Hon
>>>>
>>>> *LinkedIn Profile
>>>> <
>>>> http://www.linkedin.com/profile/preview?locale=en_US&trk=prof-0-sb-preview-primary-button
>>>> >
>>>> **|* +61420 382 557
>>>>
>>>> Monash Institute for Pharmaceutical Sciences ( *MIPS* )
>>>> Royal Parade, Parkville, 3052
>>>>
>>>>
>>>>
>>>
>>> --
>>> Billy Noonan* | *PhD Student *|* Bsci ( *Adv* ), IA Hon
>>>
>>> *LinkedIn Profile
>>> <
>>> http://www.linkedin.com/profile/preview?locale=en_US&trk=prof-0-sb-preview-primary-button
>>> >
>>> **|* +61420 382 557
>>>
>>> Monash Institute for Pharmaceutical Sciences ( *MIPS* )
>>> Royal Parade, Parkville, 3052
>>>
>>>
>>>
>>
>>
>
> --
> David van der Spoel, Ph.D., Professor of Biology
> Dept. of Cell & Molec. Biol., Uppsala University.
> Box 596, 75124 Uppsala, Sweden. Phone: +46184714205.
> spoel at xray.bmc.uu.se http://folding.bmc.uu.se
> --
> 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.
>
--
Billy Noonan* | *PhD Student *|* Bsci ( *Adv* ), IA Hon
*LinkedIn Profile
<http://www.linkedin.com/profile/preview?locale=en_US&trk=prof-0-sb-preview-primary-button>
**|* +61420 382 557
Monash Institute for Pharmaceutical Sciences ( *MIPS* )
Royal Parade, Parkville, 3052
More information about the gromacs.org_gmx-users
mailing list