From: Peter Freddolino (petefred_at_ks.uiuc.edu)
Date: Fri Jun 05 2009 - 10:16:29 CDT
Hi Stephen,
what Langevin damping constant are you using, and have you tried a lower
value to see if it is affecting your results? I would be surprised if
the thermostat alone were responsible for 10x too small diffusion, but
it certainly could account for a lot of what you're seeing in the
pulling experiment:
If you have a Langevin damping constant of 5 ps^-1, for example, then if
you pulled your protein through vacuum at 0.8 m/s you would have (if
I've done my math right and your protein is about 8000 amu, which seems
right for 1AUM) a force opposing the velocity of 0.764 kcal / mol A.
That means you're not, effectively, pulling nearly as hard as you think
you are when you take the friction into account (and ignores all manner
of more complicated effects when you consider the first couple of
solvation shells as well). This is precisely why people generally do NVE
simulations with SMD.
So, you might want to try both of your experiments again with a much
weaker damping constant and see how that affects things. I'd also note
that obtaining a longer trajectory, maybe for the smaller system,
wouldn't hurt...
Best,
Peter
Stephen Hicks wrote:
> Hi,
>
> I'm attempting to measure the diffusion constant of a more or less
> spherical protein, the C-terminal domain of HIV capsid protein 1AUM.
> It has 70 amino acids and I estimate the radius to be about a=1.3nm.
> At 310K I estimate the viscosity of water to be anywhere from
> eta=.25cP to .7cP (dependent on the model - I understand TIP3P is much
> lower than experiment). So I ran NVT simulations in a large TIP3P
> water box (10nm to a side) with a 0.8fs timestep (langevin thermostat,
> rigid bonds) and measured D as follows:
>
> D = <(x(t+t')-x(t'))^2>/6t
>
> where x is the (3-vector) position of the protein and the <..> are
> averaging over times t'. For the t part, I calculate the average as a
> function of t and fit the part near zero to a line to compute the
> slope, which is my diffusion constant. But when I do this, my result
> is always on the order of about 3e-7 cm^2/s, while I expect, using
> Einstein/Stokes, that
>
> D = kT/(6*pi*eta*a) ~ 2.5e-6 cm^2/s (using .7cP) or 7e-6 cm^2/s (using .25cP).
>
> So at best (using the experimental viscosity rather than the TIP3P
> viscosity) I'm off by nearly an order of magnitude. I've done this at
> different box sizes (4.3nm sides) and gotten the same D, so I'm
> confident that it's not a finite-size effect. I also tried another
> approach, applying a constant force of F=1.093kcal/mol*A =
> (6*pi*.25cP*1.3nm)*12m/s in the +x direction, distributed evenly (by
> mass) among all the atoms in the protein (i.e. F_i = F * m_i/m_total),
> and found that the protein was drifting with a velocity of about
> 0.8m/s - again about an order of magnitude too small!
>
> (On the other hand, I've integrated a box with nothing but water and
> got more or less the correct experimental self-diffusion constant,
> 3.5e-5 cm^2/s -- which may be only a factor of 2-3 too small if TIP3P
> self-diffusion is supposed to be larger than experiment, but then my
> protein results are off even more, by the same factor. I've also
> tried to use TIP4P but so far I haven't been able to get NAMD to run
> properly with it)
>
> I've seen a variety of papers that claim to measure the translational
> diffusion of proteins with MD, so my question is, is there any known
> effect that would be throwing off my results? Why can't I reproduce
> reasonable numbers here? Am I missing something in my integration, or
> my water, or the way I'm measuring?
>
> Any help is greatly appreciated!
> Steve Hicks
> Ph.D. Candidate, Henley Group
> Laboratory of Atomic and Solid State Physics
> Cornell University
This archive was generated by hypermail 2.1.6 : Wed Feb 29 2012 - 15:52:54 CST