[gmx-users] Pulling ligand - Different Profiles (Force vs time)
Justin A. Lemkul
jalemkul at vt.edu
Wed Jun 27 15:39:54 CEST 2012
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
========================================
More information about the gromacs.org_gmx-users
mailing list