[gmx-users] deltaG from PMF

[Justin Lemkul]

On 8/21/12 10:42 AM, Steven Neumann wrote:
> Please see the example of the plot from US simulations and WHAM:
>
> http://speedy.sh/Ecr3A/PMF.JPG
>
> First grompp of frame 0 corresponds to 0.8 nm - this is what was shown
> by grompp at the end.
>
> The mdp file:
>
> ; Run parameters
> define      = -DPOSRES_T
> integrator  = md        ; leap-frog integrator
> nsteps      = 25000000     ; 100ns
> dt          = 0.002     ; 2 fs
> tinit       = 0
> nstcomm     = 10
> ; Output control
> nstxout     = 0       ; save coordinates every 100 ps
> nstvout     = 0       ; save velocities every
> nstxtcout   = 50000        ; every 10 ps
> nstenergy   = 1000
> ; Bond parameters
> continuation    = no           ; first dynamics run
> constraint_algorithm = lincs    ; holonomic constraints
> constraints     = all-bonds     ; all bonds (even heavy atom-H bonds)
> constrained
> ; Neighborsearching
> ns_type     = grid      ; search neighboring grid cells
> nstlist     = 5         ; 10 fs
> vdwtype     = Switch
> rvdw-switch = 1.0
> rlist       = 1.4       ; short-range neighborlist cutoff (in nm)
> rcoulomb    = 1.4       ; short-range electrostatic cutoff (in nm)
> rvdw        = 1.2       ; short-range van der Waals cutoff (in nm)
> ewald_rtol  = 1e-5      ; relative strenght of the Ewald-shifted
> potential rcoulomb
> ; Electrostatics
> coulombtype     = PME       ; Particle Mesh Ewald for long-range electrostatics
> pme_order       = 4         ; cubic interpolation
> fourierspacing  = 0.12      ; grid spacing for FFT
> fourier_nx      = 0
> fourier_ny      = 0
> fourier_nz      = 0
> optimize_fft    = yes
> ; Temperature coupling is on
> tcoupl      = V-rescale                     ; modified Berendsen thermostat
> tc_grps     = Protein LIG_Water_and_ions   ; two coupling groups - more accurate
> tau_t       = 0.1   0.1                     ; time constant, in ps
> ref_t       = 318   318                     ; reference temperature,
> one for each group, in K
> ; Pressure coupling is on
> pcoupl      = Parrinello-Rahman             ; pressure coupling is on for NPT
> pcoupltype  = isotropic                     ; uniform scaling of box vectors
> tau_p       = 2.0                           ; time constant, in ps
> ref_p       = 1.0                           ; reference pressure, in bar
> compressibility = 4.5e-5                    ; isothermal
> compressibility of water, bar^-1
> ; Periodic boundary conditions
> pbc         = xyz       ; 3-D PBC
> ; Dispersion correction
> DispCorr    = EnerPres  ; account for cut-off vdW scheme
> ; Velocity generation
> gen_vel     = yes       ; assign velocities from Maxwell distribution
> gen_temp    = 318       ; temperature for Maxwell distribution
> gen_seed    = -1        ; generate a random seed
> ; These options remove COM motion of the system
> ; 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         = 500      ; kJ mol^-1 nm^-2
> pull_nstxout    = 1000    ; every 2 ps
> pull_nstfout    = 1000      ; every 2 ps
>
>

Based on these settings you're allowing grompp to set the reference distance to 
whatever it finds in the coordinate file.  It seems clear to me that the 
sampling indicates what I said before - you have an energy minimum somewhere 
other than where you "started" with.  What that state corresponds to relative to 
what you think is going on is for you to decide based on the nature of your 
system.  Whatever is occurring at 0.6 nm of COM separation is of particular 
interest, since the energy minimum is so distinct.

I hope you're doing a thorough analysis of convergence if you're generating 
velocities at the outset of each run, and removing unequilibrated frames from 
your analysis.

-Justin

-- 
========================================

Justin A. Lemkul, Ph.D.
Research Scientist
Department of Biochemistry
Virginia Tech
Blacksburg, VA
jalemkul[at]vt.edu | (540) 231-9080
http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin

========================================