Re: SMD and CPT? (fwd)

From: Gianluca Interlandi (gianluca_at_u.washington.edu)
Date: Fri Feb 09 2007 - 13:52:23 CST

I'm sending Marcos' answer to my last e-mail to the whole list as many
people might be interested in our discussion.

Marcos, Thank you very much for your help.

Best,

   Gianluca

---------- Forwarded message ----------
Date: Fri, 9 Feb 2007 13:43:32 -0600 (CST)
From: Marcos Sotomayor <sotomayo_at_ks.uiuc.edu>
To: Gianluca Interlandi <gianluca_at_u.washington.edu>
Subject: Re: namd-l: SMD and CPT?

Dear Gianluca,

Two comments on your plots:

- I see you are using a 2fs time step. Are you using rigid bonds and/or multiple
time stepping? (I guess you are using rigid bonds and not multiple time stepping
since the temperature seems to be quite under control, good!)

- Check the energy of your system (TOTAL3 in namd output) and compute the
difference between the starting value and the final value (due to intrinsic
fluctuations of the energy you will be off by ~ +/- 10 kcal/mol). Then plot
Force vs distance (I assume the plot you sent is F vs t). The area under Force
vs distance (you can integrate F vs d numerically using xmgrace, be sure to use
the appropriate units) should be similar (not identical since you are applying a
time dependent force) to the energy difference you computed before. That's a
good way to check that the simulation is OK in terms of the physics behind it.

The temperature rises at the breaking point because that's when the external
force starts to do work (F \times d). If you plot the energy of your system,
likely you will see that increases at the breaking point as well.

Although I feel more comfortable with NVE simulations, do not completely
disregard your CPT simulation, stiffer springs will give you a lot of large
fluctuations! Likely that's the reason why you see forces on the order of 1200pN
and probably you will see something similar in your NVE simulation. Plot F vs t
and a running average (just like you did for the temperature) to see at what
level is the force.

Finally, the choice of spring constant is a tough one. AFM people use really
soft cantilevers, while if you want to get a potential of mean force (PMF) you
need to use a stiff spring (have a look at Park and Schulten, Journal of
Chemical Physics, 120:5946-5961, 2004.)

Hope this helps and don't hesitate to contact me again if you have further
questions.

Regards,
Marcos
PS: perhaps you could resend this e-mail (without the attachment) to namd-l,
since I think many people would benefit from our discussion even without the
plot.

On Fri, 9 Feb 2007, Gianluca Interlandi wrote:

> Dear Marcos,
>
> Thanks a lot for your answer and for offering me your help.
>
> I have already run a bunch of pulling simulations in the NVE ensemble. I
> attach a plot of a constant velocity simulation where I'm pulling two
> proteins apart. At the point of breaking the temperature raises of a few
> degrees. The force peak is at around 670 pN. Recently, I have started a
> new simulation of the same system but with a stiffer spring, i.e., 140
> pN/A instead of 14 pN/A, and in the CPT (Langevin as thermostat and
> Langevin-Nose Hoover as barostat). After 5 ns the force has reached 1200
> pN but the two proteins are still bound and the total RMSD is still around
> 2 A. I started wondering whether this is due to the stiffer spring only or
> whether CPT has introduced some artifacts. To investigate this I have
> started a new simulation with exactly the same parameters but in the NVE.
> The latter simulation has just started so I'm waiting for the results.
>
> In any case, now I know that I have to stick with NVE.
>
> Many thanks,
>
> Gianluca
>
> On Fri, 9 Feb 2007, Marcos Sotomayor wrote:
>
> >
> > Actually, I have three reasons to disagree with Sterling.
> >
> > 1 - Performing simulations in the NVE ensemble permits you to monitor the
> > change in energy of your system, and check that it is at least similar to
> > the
> > work done by the external applied force (and exactly the same when using
> > constant forces). I've monitored temperature changes of 0.5 K in simulations
> > that do not use multiple time stepping, which matches well the work done on
> > the system by external forces.
> >
> > 2 - Algorithms used to control temperature usually are not designed to
> > handle
> > systems in which external forces induce motions in a preferred direction.
> > The
> > Langevin thermostat for instance, will apply a net (average) zero force to
> > atoms that do not move on a preferred direction, but likely will introduce
> > an
> > artificial viscous drag to atoms that are being pulled. Thus, your force
> > peaks
> > will be larger.
> >
> > 3 - Some temperature control methods induced center of mass motion, which is
> > not good when you have fixed reference points like in SMD.
> >
> > The best solution is to perform simulations in both ensembles and compare,
> > focusing on those results that are "ensemble independent".
> > However, with a water box big enough I suggest to perform SMD simulations in
> > the NVE ensemble.
> >
> > Ginaluca, that's the reason why we performed simulations in the NVE ensemble
> > in our Structure paper. Let me know if you have more questions, I'll be glad
> > to discuss them with you.
> >
> > Regards,
> > Marcos
> >
> >
> > On Fri, 9 Feb 2007, Sterling Paramore wrote:
> >
> > > I would definitely suggest you do SMD with a thermostat. Otherwise, the
> > > viscous heating that occurs when you pull on your system will increase the
> > > temperature. In any real system, the heat produced would be transferred
> > > to
> > > the surrounding thermal reservoir. The artificial thermostat is just a
> > > way
> > > to model this effect.
> > >
> > > -Sterling
> > >
> > >
> > > On Feb 8, 2007, at 11:49 PM, Gianluca Interlandi wrote:
> > >
> > > > I have a question concerning steered molecular dynamics simulations
> > > > (constant force and constant velocity). Is it appropriate to use a
> > > > thermostat and barostat (CPT) while performing a constant force or
> > > > constant velocity pulling simulation? I have seen that many people
> > > > prefer
> > > > NVE, i.e., no thermostat.
> > > >
> > > > I expect my protein to undergo large conformational changes during
> > > > pulling. Does a thermostat slow down the sequence of events, since part
> > > > of
> > > > the applied force is converted into heat?
> > > >
> > > > Many thanks,
> > > >
> > > > Gianluca
> > > >
> > > > -----------------------------------------------------
> > > > Dr. Gianluca Interlandi gianluca_at_u.washington.edu
> > > > +1 (206) 685 4435
> > > > +1 (206) 714 4303
> > > > http://biocroma.unizh.ch/gianluca/
> > > >
> > > > Postdoc at the Department of Bioengineering
> > > > at the University of Washington, Seattle WA U.S.A.
> > > > -----------------------------------------------------
> >
>
> -----------------------------------------------------
> Dr. Gianluca Interlandi gianluca_at_u.washington.edu
> +1 (206) 685 4435
> +1 (206) 714 4303
> http://biocroma.unizh.ch/gianluca/
>
> Postdoc at the Department of Bioengineering
> at the University of Washington, Seattle WA U.S.A.
> -----------------------------------------------------

This archive was generated by hypermail 2.1.6 : Wed Feb 29 2012 - 15:46:01 CST