From: Jan Saam (saam_at_charite.de)
Date: Wed May 31 2006 - 17:27:03 CDT
Thanks a lot for your comments, Jerome!
I see now that it is more complicated than I initially thought.
Jerome Henin wrote:
>Dear Jan and NAMDers,
>>>Another problem with this, though, is to define an absolute direction. Is
>>>your protein able to tumble or rotate in some way?
>>Unfortunately yes, it will tumble. We should urge Jim to provide a
>>function to remove the COM motion constantly throughout the simulation,
>>not just once at the beginning. That would help a lot defining RCs.
>The problem is, conserving the total momentum of the simulation box will not
>prevent the angular momentum of a solvated protein from fluctuating. I can
>think of no way to prevent that.
That is correct, one has to clamp only the protein to the coordinate
axes. I'm not entirely sure about it's thermodynamic effect, but I
assume as soon as you have a thermostat running, it wouldn't make a
I've taken a look into all possibilities that NAMD's TCL interface
provides. The momenta can be retrieved easily through the callback
command. The problem is that its a pain to reassign all velocities every
N steps so that the momenta are compensated (especially for the angular
momentum because each atom get its own additional velocity term).
As a simple minded countermeasure I just added weak harmonic constraints
on a bunch of helix backbone atoms far from the region of interest. This
will keep the system in place.
>>>>c) If I'm partitioning the channel into several shorter, almost straight
>>>>portions, I can stitch the PMF together and probably improve on these
>>>>problems. But I'll have to know the RC very well before and it is not
>>>>easy to define the directions in a system where everything moves.
>>>I would say the only useful coordinates would then be defined relative to
>>>the protein itself. The available RCs are probably not adapted for that,
>>>but we can talk about what would suit your needs. It would probably be
>>>helpful to many other people, too.
>>I think it would be cool to spiff up ABF (which I really have to say is
>I usually say that the advantage of ABF is that the idea is extremely simple.
>But I am sure some who are into more elaborate methods would rather call this
In 400 years of physics the most elegant and ssimplest methods that
described a system with the requested accuracy alway won ;-)
>>so that it automatically find the reaction coordinate and hence transition
>Sorry, but I think you are asking a bit too much from ABF. Smart people have
>been working hard on sophisticated ways to do that (e.g. transition path
>sampling), and I do not believe ABF is really fit for that purpose.
Maybe (see below).
>>This might not be possible for all types of reactions but in my case I
>>could envision the following:
>>You probably know "conformational flooding" due to H. Grubmüller or the
>>similar method "metadynamics" by Parinello. I'll quickly outline
>>metadynamics since it's suited better to our problem: The basic idea is
>>that you obtain information about your free energy profile along a
>>certail RC by adding a a little gaussian function to the total force
>>field potential for the current state. Each sample adds one more little
>>gaussian contribution to a potential that builds up slowly filling the
>>energy minimum you are sitting in. Figuratively speaking you are
>>dropping litter into your well until you are carried out of the pit by
>>all that trash. The potential that you built up by then is something
>>like the inverse of the PMF along the RC (very elegant method, too!).
>>I think this gaussian accumulation technique could also be used to
>>actually *find* the most likely RC. At least in such a simple case as a
>>diffusiing particle in a channel this should be possible.
>I would say things similar in spirit have been done (see for instance Ensing
>and Klein, PNAS 2005).
Thanks for the reference, a nice example of using metadynamics.
In a way metadynamics or "the hills method" as Ensing and Klein term it
is even more elegant since the math is fantastically simple as compared
to the stuff worked out in Darve and Pohorille work.
My question is do you have an opinion what are the advantages and
disadvantages of the two methods? When would you rather suggest using ABF?
>>But may be there's even a simpler way without running a simulation
>>before to get the RC, since ABF has the same effect as the gaussian
>>Wouldn't it be possible to replace Xi by a 3-dimensional analog Xi_j ;
>>j=1,2,3 so we would get:
>>dA(xi_j)/dXi_j = <partial(V(x))/partial(Xi_j) -
>>1/beta*partial(ln|J_j|)/partialXi_j>_(Xi_j) = -<F_Xi_j>_Xi_j
>>(Underscores denote subscripts)
>Things are a bit more complicated than this, because we would need
>partial(ln|J|)/partialXi_j, where J is a Jacobian for a coordinate transform
>involving all the RCs. In some cases, it may be trivial, and very painful in
Hmm. As a starting point 3D cartesian coords would be sufficient and I
guess that's rather simple?
>>We would have j^3 bins but we would only use a small fraction of them,
>>namely the ones that have more than a certain number of samples. So we
>>would get a 3D potential map of the accessible conformations and the
>>longer the simulation is run the more the 3D subspace is explored.
>>Since the bins with no samples in them are even in the 1D case treated
>>as having the highest energy while all other bins have a lower energy
>>this is ideal to be visualize using a volume denstity isosurface in VMD
>>(you'd have to invert the potential first). Then, depending on the
>>chosen free energy isovalue you'd see blobs or tubes with an energy
>>lower than that.
>>The problem here is again, that one would have to make the RC grid
>>relative to a protein coordinate system or constrain the COM motions and
>>rotations of the system.
>>I also thought about case a) again. When the system has spent some time
>>ina dead-end even if that's in a potential well this minimum will be
>>leveled out after enough samples have been taken and the particle can
>>travel back onto the main road it will finally first break the lowest
>>potential barrier wether it's on the branch 1 or 2 of the way.
>Or not. This is a timescale issue: how many samples is "enough samples" ?
Of couse. Enough samples here means enough for the particle to get
untrapped. I would just require a mimimum number of samples in each bin
which depends on the pmf. Has to be decided by the user...
>>What do you think about this. Do you think it could work? Would it be
>>terribly hard to implement this?
>Implementing a multidimensional version of ABF would be feasible, but probably
>at the expense of generality and flexibility with respect to the choice of
>reaction coordinates. The problem mostly arises because you'd have to handle
>a Jacobian term arising from all the coordinates. Again, it can be very
>simple in predefined cases, but being general could be a nightmare. A more
>specialized version would be a distant branch from the current state of the
>code, so I don't know if I want to do that. Maybe I can come up with an
Special cases would be sufficient and a separate branch would be ok,
too. You probably have a lot of other things to do, but I would find it
useful and would use it. Otherwise I'd try to program the hills method
-- --------------------------- Jan Saam Institute of Biochemistry Charite Berlin Monbijoustr. 2 10117 Berlin Germany +49 30 450-528-446 saam_at_charite.de
This archive was generated by hypermail 2.1.6 : Wed Feb 29 2012 - 15:42:07 CST