[gmx-users] Free energy of conformation

David Mobley dmobley at gmail.com
Wed Jun 7 00:37:56 CEST 2006

I just e-mailed suggesting umbrella sampling. Since you put some more
detail in your question, I'll give you a bit more detail in terms of
how you might go about doing that.

> I didn't get a response to this question before, so here I go again: I have
> a protein that has a residue modeled as a 2-state rotamer in the crystal
> structure. The transition occurs in ms, so I cannot see it in the MD
> simulations. I'd like to calculate the free energy difference between both
> states, possibly using the perturbation method. But I am not sure how to do
> this, as examples are usually done for mutations and hydration. Has anyone
> an idea of how to do this? So far I have thought of changing the dihedral of
> the side chain in the top file (which goes from avge -90 to 90), but I am
> thinking that maybe there are other things that I need to account for. Any
> ideas …?

In my opinion, the easiest way to do this is to run a series of
separate calculations where the dihedral is restrained to a particular
value using simple harmonic restraints (make sure you are using 3.3.1
or 3.3 with relevant bugfixes, as dihedral restraints were broken
previously). You can do this by just putting a [ dihedral_restraints ]
section in your topology (see the documentation). You restrain to a
series of different angles to cover the whole 360^o of interest, and
then use the weighted histogram analysis method (WHAM) to take the
potential energies and calculate a potential of mean force (PMF) --
basically the free energy profile for rotating the dihedral. The free
energy difference is then the appropriately weighted difference of the
log probability of being in each basin. (As a first approximation you
can just take the difference in well depths, but it isn't quite right,
especially if they have different shapes).

There is in principle a WHAM code in gromacs (g_wham), but I haven't
used it and can't vouch for it, and it probably expects some input
other than just the results of a bunch of separate simulations with
different dihedral restraints, so if you use my route, it probably
means you need your own WHAM code (or ask us for ours). Alternatively,
there is probably some way you can do the whole umbrella sampling and
PMF routine in GROMACS alone using g_wham, but I can't help you with

Things to watch out for:
- Make sure if you look at the binned dihedral angles, you have good
enough overlap in the distributions, otherwise your free energy
differences will be wrong
- Make sure your spring constant is big enough to ensure you sample
across the barriers (and you equilibrate long enough), but not so big
that you get delta function-type probability distributions.

The overlap thing is key. Otherwise you will get a PMF that looks OK
(maybe), but is completely wrong.


> Thanks!!
> Esther Caballero-Manrique
> Graduate Teaching Fellow
> Chemistry Department
> 1253 University of Oregon
> Eugene, OR 97403-1253
> (541) 346-2485
> _______________________________________________
> gmx-users mailing list    gmx-users at gromacs.org
> http://www.gromacs.org/mailman/listinfo/gmx-users
> Please don't post (un)subscribe requests to the list. Use the
> www interface or send it to gmx-users-request at gromacs.org.
> Can't post? Read
> http://www.gromacs.org/mailing_lists/users.php

More information about the gromacs.org_gmx-users mailing list