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

Steven Neumann s.neumann08 at gmail.com
Wed Jun 27 15:53:18 CEST 2012 

Thank you Justin.

On Wed, Jun 27, 2012 at 2:39 PM, Justin A. Lemkul <jalemkul at vt.edu> wrote:
>
>
> On 6/27/12 9:36 AM, Steven Neumann wrote:
>>
>> On Wed, Jun 27, 2012 at 1:51 PM, Justin A. Lemkul <jalemkul at vt.edu> wrote:
>>>
>>>
>>>
>>> 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.
>>
>>
>> Thank you Justin. I am using according to your tutorial:
>>
>> 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
>>
>> Would you suggest:
>>
>> grompp -f pull.mdp -c npt.gro -p topol.top -n index.ndx -t npt.cpt -o
>> pull.tpr
>> mdrun -s pull.tpr -cpi npt.cpt -deffnm pull ??
>>
>
> No, I would not, especially if the NPT run uses position restraints, in
> which case the two phases are different.  I missed the command line in the
> earlier message.  What you are doing makes sense.
>
>
>> Profiles do not vary slightly - the maximum pulling force (breaking
>> point) varies from 290 to 500 kJ/mol nm2 which is really a lot.
>>
>
> Consult the points below and watch your trajectories.  If you're getting
> different forces, your ligands are experiencing different interactions.  SMD
> is a path-dependent, non-equilibrium process.  Good sampling and a
> justifiable path are key.
>
> -Justin
>
>
>>>
>>> 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
>>>
>>> ========================================
>>>
>>>
>>> --
>>> gmx-users mailing list    gmx-users at gromacs.org
>>> http://lists.gromacs.org/mailman/listinfo/gmx-users
>>> * Only plain text messages are allowed!
>>> * Please search the archive at
>>> http://www.gromacs.org/Support/Mailing_Lists/Search before posting!
>>> * 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/Support/Mailing_Lists
>>
>>
>
> --
> ========================================
>
> 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
>
> ========================================
>
>
> --
> gmx-users mailing list    gmx-users at gromacs.org
> http://lists.gromacs.org/mailman/listinfo/gmx-users
> * Only plain text messages are allowed!
> * Please search the archive at
> http://www.gromacs.org/Support/Mailing_Lists/Search before posting!
> * 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/Support/Mailing_Lists

More information about the gromacs.org_gmx-users mailing list