[gmx-users] Pulling ligand - Different Profiles (Force vs time)

[Justin A. Lemkul]

On 6/27/12 7:48 AM, Steven Neumann wrote:
> Dear Gmx Users,
>
> I obtained a protein-ligand complex from 100ns simulation. Now I am
> pulling my ligand away from the protein after the energy minimzation
> in water and equilibration of 100ps (two coupling baths: Protein,
> LIG_Water_and_ions).
> Then I proceed my pulling :
>
> grompp -f pull.mdp -c npt.gro -p topol.top -n index.ndx -t npt.cpt -o pull.tpr
>
> mdrun -s pull.tpr -deffnm pull
>
>
> title       = Umbrella pulling simulation
> define      = -DPOSRES
> ; Run parameters
> integrator  = md
> dt          = 0.002
> tinit       = 0
> nsteps      = 500000    ; 1 ns
> nstcomm     = 10
> ; Output parameters
> nstxout     = 0
> nstvout     = 0
> nstfout     = 500
> nstxtcout   = 1000       ; every 1 ps
> nstenergy   = 500
> ; Bond parameters
> constraint_algorithm    = lincs
> constraints             = all-bonds
> continuation            = yes       ; continuing from NPT
> ; Single-range cutoff scheme
> nstlist     = 5
> ns_type     = grid
> rlist       = 1.4
> rcoulomb    = 1.4
> rvdw        = 1.2
> vdwtype     = Switch
> rvdw-switch = 1.0
> ; PME electrostatics parameters
> coulombtype     = PME
> fourierspacing  = 0.12
> fourier_nx      = 0
> fourier_ny      = 0
> fourier_nz      = 0
> pme_order       = 4
> ewald_rtol      = 1e-5
> 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       = 298   298                                  ; reference
> temperature, one for each group, in K
> ; Pressure coupling is on
> Pcoupl          = Parrinello-Rahman
> pcoupltype      = isotropic
> tau_p           = 2.0
> compressibility = 4.5e-5
> ref_p           = 1.0
> ; Generate velocities is off
> gen_vel     = no
> ; Periodic boundary conditions are on in all directions
> pbc     = xyz
> ; Long-range dispersion correction
> DispCorr    = EnerPres
> ; Pull code
> pull            = umbrella
> pull_geometry   = distance  ; simple distance increase
> pull_dim        = N N Y
> pull_start      = yes       ; define initial COM distance > 0
> pull_ngroups    = 1
> pull_group0     = Protein
> pull_group1     = LIG
> pull_rate1      = 0.004      ; 0.004 nm per ps = 4 nm per ns
> pull_k1         = 500      ; kJ mol^-1 nm^-2
>
> I run 3 pulling simulations with the same mdp  and I obtain 3
> different profiles (Force vs time). Then I used 2xlonger pulling with
> the same pulling distance and I run 3 simulations again. Each time I
> obtain different profile. Can anyone explain me this? I am using
> velocities from npt simulation as above (gen_vel = no and continuation
> = yes) so I presume the output should be similar. Please, advice.
>

I assume you're passing a checkpoint file to grompp?  If you're relying on 
velocities from the .gro file, they are of insufficient precision to guarantee 
proper continuation.

Small variations are inherent in any simulation set, and in the case of pulling, 
small changes (though intentional) are the basis for Jarzynski's method.  In any 
case, all MD simulations are chaotic and so it depends on what your definition 
of "different" is in the context of whether or not there are meaningful changes 
imparted through the course of each simulation.  Also note that in the absence 
of the -reprod flag, the same .tpr file may result in a slightly different 
outcome.  The implications of these outcomes are limited by sampling; the 
ensemble should converge with sufficient time and/or replicates.  For 
non-equilibrium processes like pulling, convergence is probably harder, but 
again you have to ask whether the differences are meaningful.

http://www.gromacs.org/Documentation/Terminology/Reproducibility

-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

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