[gmx-users] Umbrella Sampling and flat-bottom position restraints

Oliver Stueker ostueker at gmail.com
Sun Nov 16 23:45:02 CET 2014 

Hi Justin,

Thanks for the feedback.

On Sun, Nov 16, 2014 at 2:51 PM, Justin Lemkul <jalemkul at vt.edu> wrote:

>
>
> On 11/16/14 1:04 PM, Oliver Stueker wrote:
>
>> Dear fellow Gromacs Users,
>>
>> I'm trying to calculate the binding free energy between a protein and a
>> small molecule using Umbrella Sampling.
>>
>> I made my first attempt last year using Gromacs 4.5 and ran into the same
>> convergence problems that have been pointed out on the mailing list
>> several
>> times.
>> To overcome these problems I now want to give it another try utilizing the
>> flat-bottomed position restraints that have been introduced in Gromacs
>> 5.0.
>>
>> My plan is to perform the Umbrella sampling along the Z-axis and to
>> restrict the protein as well as the ligand with flat-bottomed position
>> restraints along X and Y into a cylindrical space (as suggested e.g. by
>> Chris Neale in [1]).
>>
>> But I'm not quite sure what radius the cylinder should have ( or how wide
>> the flat bottom should be).
>> The protein is roughly spherical with a diameter of ~4nm. I was initially
>> thinking to use an r_{fb} of 2nm for the protein and  for the ligand 1/2
>> of
>> it's longest dimension to allow both of them to re-orient freely. But the
>> more I think about, I'm afraid that it's too large.
>>
>>
> Why is it necessary to pull only along the z-axis?  For simple ligand
> binding to a globular protein, you can define the pull vector in any or all
> of the three spatial dimensions.  Is there some reason you presuppose
> binding only along z?
>

​Hmm, indeed it's probably not necessary to ​pull only along Z. However the
box would need to be much larger if I can't control the direction.
In my first attempt I had already aligned the vector between the COMs of
Protein and ligand with the Z-axis so that I can pull along that direction.
That time I've used a box of 10 x 10 x 15 nm.

Last time I've used the following pull settings:
; Pull code
pull            = umbrella
pull_geometry   = distance
pull_dim        = N N Y
pull_start      = yes
pull_ngroups    = 1
pull_group0     = Protein
pull_group1     = LIG
pull_init1      = 0
pull_rate1      = 0.0
pull_k1         = 1000      ; kJ mol^-1 nm^-2
pull_nstxout    = 1000      ; every 2 ps
pull_nstfout    = 1000      ; every 2 ps

​I also had the problem that during the equillibration of the Umbrella
windows the ligand ​moved a bit away from it's initial position so that the
initial pull distance was a bit off. This time I want to prevent that by
setting pull_init for each window directly and leave pull_start = off.

​Say if I ​use pull_dim = Y Y Y and a cubic box of 15 x 15 x 15 nm. That
would restrict the ligand to the surface of a sphere of radius d (COM
reference distance of the umbrella window) around the protein. Wouldn't
that mean that I would need to sample the ligand conformations on the
complete surface of each sphere to achieve sufficient conformational
overlap between neighbouring windows? Those surfaces would become quite big
once I reach to high COM distances (6 nm).
Say if in two windows n-1 and n+1 the ligand moves towards +X / +Y and in
window n (in between) the ligand moves to -X / -Y. Then while the
histograms might overlap nicely (according the sampled COM distances) the
sampled conformations might be very different.

​It seems that simply using pull_dim = Y Y Y would also not really solve
the problem.​



>  Has anyone have experience with that? Maybe even a published paper?
>> So far I've unsuccessfully tried several times to find a published method.
>>
>>
>> Also to I need to consider the flat-bottomed potential when doing the
>> g_wham analysis?
>>
>>
> It's an additional bias to the Hamiltonian, so it should probably be
> accounted for in some way.  As to how, I'm not sure.  But it's an
> additional degree of freedom that is being restrained.
>

​Yes, you are probably right.

Would this work:

   - generate conformations by pulling along Z and select conformations for
   windows
   - equillibrate each window
   - run umbrella-sampling (to sample conformations/coordinates) along Z
   with pull_init = d(window) and use flat-bottomed PosRes to restrict Protein
   and Ligand in X/Y movement.
      - as reference coordinates, I use the COM of Protein and ligand
      respectively
      - as r_{fb} I can use the radius of the sphere that encloses the
      protein (or ligand) plus a small buffer.
      - store coordinates in XTC and discard the pullx- / pullf- files.
   - Do a re-run umbrella-sampling (to generate distances/forces) with pull_dim
   = Y Y Y and without the PosRes to generate the pullx / pullf files that
   will be used for g_wham.

​Instead of the last re-run I could maybe also run g_dist to generate my
pullx.xvg files and then pass g_wham tpr files with pull_dim = Y Y Y and
without the PosRes.


Best,
Oliver



> -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.
>

More information about the gromacs.org_gmx-users mailing list