[gmx-users] Fwd: Question on gmx_wham with two reaction coordinates.

[Yuanchao Liu (MSU)]

Hi all

Sorry to post this again. I will appreciate if anybody can help with this

I have a question on dealing with two reaction coordinates by gmx_wham. I
did umbrella sampling on glucose 6-phosphate molecule on a poly-Lys peptide
surface. Because the peptide is neither cylindrical nor spherical
symmetric, I have to restraint each window around a specific point, using
two reaction coordinates (distance and direction). As results, the ligand's
trajectory looks pretty well. It was restraint around a specific points,
instead of going around a sphere surface when I used 'distance' alone.

However, I got some issue using WHAM for analysis. That is, GROMACS seems
not to allow me to use WHAM with two pull geometry (distance and
direction). So I have to pick only one reaction coordinate to do the
analysis. When I run regular wham, I got this issue:

Fatal error:

umbrella0.tpr: Your pull coordinates have different pull geometry
(coordinate

1: distance, coordinate 2: direction)

If you want to use only some pull coordinates in WHAM, please select them
with

option gmx wham -is


I refer this link to pick one of the two reaction coordiantes to work on
WHAM.
https://urldefense.proofpoint.com/v2/url?u=http-3A__manual.gromacs.org_documentation_2016_onlinehelp_gmx-2Dwham.html&d=DwIFaQ&c=nE__W8dFE-shTxStwXtp0A&r=gCy395jfSTIKIGBjdJ8SwSg8yN7wywRqF0ilx8ZoefpMNPZKIUY0jXzHKcFUwcQV&m=erg0OZLT3d1xG1FDq6J357V6VeRaFwNylkFDDHWoibM&s=yaLzHjzylI66bC_grW7Y2MuN-XHXdno8fH3jtTGa5AI&e= 

See all the figures in this link:
https://urldefense.proofpoint.com/v2/url?u=http-3A__www.evernote.com_l_AAMsNYa6HXZFeqNDTbsiV4HO0ta2KY-2Dnz7s_&d=DwIFaQ&c=nE__W8dFE-shTxStwXtp0A&r=gCy395jfSTIKIGBjdJ8SwSg8yN7wywRqF0ilx8ZoefpMNPZKIUY0jXzHKcFUwcQV&m=erg0OZLT3d1xG1FDq6J357V6VeRaFwNylkFDDHWoibM&s=JmmnE853X5gvKEdb7Ru78xmCvsQ5foSjLDnAhhlitLY&e= 

When I only selected 'distance' for WHAM, the histogram and PMF looks like
this (Figure 1). I can somewhat understand the insufficient overlap at
0.5-0.75 nm could result in some abnormal segment in PMF, but what really
confused me is the global decreasing on PMF. Then histogram looks just
right.

Here is the result when only selecting 'direction' as the reaction
coordinate for wham (Figure 2). Gromacs's WHAM still output the distance as
x-axis, and the histogram is different form the above one.

Then I just simply add these two PMFs together and got following plots.
(Figure 3) The summation looks like what I expected. I can not understand
this. Because the restraints by 'distance' and 'direction' are
perpendicular to each other, the coordinate-1's WHAM result should be
something like what I got from single reaction coordinate. I expected the
coordinate-1's WHAM result should look like the summation result. I do not
have experience on dealing with two reaction coordinates. So, I will
appreciate it if anybody has suggestions on how to figure this out.

Here is the pulling parameter for my umbrella sampling: Pull code

pull                    = yes
pull_ngroups            = 2              ; there are two pulling groups
pull_ncoords            = 2             ; two reaction coordiante to
restraint the molecule around its initial position
pull_group1_name        = Other          ; ligand
pull_group2_name        = ref_group      . reference group, it is a part of
my peptide
; first reaction coordinate
pull_coord1_type        = umbrella      ; harmonic biasing force
pull_coord1_geometry    = distance      ; simple distance increase
pull_coord1_groups      = 1 2           ; two groups the pulling will be
applied
pull_coord1_dim         = Y Y Y
pull_coord1_rate        = 0.0          ; 0.01 nm per ps = 10 nm per ns
pull_coord1_k           = 1000          ; kJ mol^-1 nm^-2
pull_coord1_start       = yes           ; define initial COM distance > 0
; second reaction coordinate
pull_coord2_type        = umbrella      ; harmonic biasing force
pull_coord2_geometry    = direction      ; simple distance increase
pull_coord2_vec          = 11.018  -20.969  13.824
pull_coord2_groups      = 1 2
pull_coord2_dim         = Y Y Y
pull_coord2_rate        = 0.0          ; 0.01 nm per ps = 10 nm per ns
pull_coord2_k           = 1000          ; kJ mol^-1 nm^-2
pull_coord2_start       = yes           ; define initial COM distance > 0



Best

Yuanchao Liu (刘元超)
Ph.D. Candidate
Department of Chemical Engineering and Materials Science
428 South Shaw Lane, Room 3263
Michigan State University
East Lansing, Michigan 48824-4437

Email: ycLiu1987 at gmail.com; liuyuan6 at msu.edu; liuyuan6 at egr.msu.edu
Cell:+1 (517)-763-1806 <+1%20517-763-1806>