Re: [AMBER] Calculating lipid diffusion

From: Daniel Roe <>
Date: Wed, 20 Aug 2014 08:57:14 -0600


On Tue, Aug 19, 2014 at 8:07 PM, Vijyan Achari <> wrote:
> Since I have simulated for around 300 ns long time, I have used options like nscm=1000 ,iwrap=1 during the simulation in MD-input file. Upon finishing the simulation, I unwrap the trajectories. (It has been indicated in the manual that, unwrapping is important for calculating thediffusion properties).

Note that if you really want a good estimate of diffusion in your
system you should simulate in the NVE ensemble to avoid effects from
the thermostat etc.

> Then I used script (as below) to calculate the diffusion of each lipid throughout 45000 frames for all the 128 lipids I used. Finally I averaged the MSD vs time data for all the 128 lipids, and I get a plot as I attached in this mail.

You don't have to do the averaging yourself since the 'diffusion'
command does it for you. Just specify a mask that includes atoms of
interest. For example, if you were calculating the diffusion of water
you could just specify ':WAT.O'. This will improve your statistics
(and results) quite a bit.

> My questions:
> 1) Is my MSD vs time plot correct? If is not than could I get some explanation what might have went wrong in my calculation?

That's really difficult to say without knowing exactly how you did
your averaging and what you are considering "correct". The line you
obtain does not have much of a slope so it does seem like it would be
difficult to obtain a good estimate of the diffusion constant from it.
If you haven't already I highly suggest you consult previous work in
which diffusion constants were calculated from simulations of lipid
systems, see e.g.

> 2) There was no option given for mass weighted diffusion calculation as in the GROMACS (g_msd). But option average is given in AMBER. Does this average option represents Mass weighted MSD?

No, the 'diffusion' command only calculates displacement of atoms.
Instead of mass-weighted MSD (which doesn't make sense to me) I think
you mean calculating the MSD for the center of mass of a group of
atoms, right? The 'stfcdiffusion' command can do this via the 'com'
keyword for a single mask. However if you do this you will have to do
the averaging over each molecule manually.

> Aditionally, visually I have seen the unwrapped trajectory file using VMD and I noticed the lipids remain in bilayer arrangement but the water molecules separate and go out of the box and filled vast space. (I dont see the lipid diffuse far apart from each other).
> I am little worry of my data here. Could I get some explanation on the MSD plot that I got (attached here ) and possible errors might occurred if the plot is wrong?

Assuming your bilayer fills the XY plane it "feels" itself via PBC in
those directions, so the lipids should not separate (if they did that
would indicate problems with your system). Water on the other hand is
comparatively free to diffuse where it likes, so I think the behavior
you observe is reasonable. I think you may need to repeat your
simulation under better conditions, recalculate the MSD and then
compare that to experimental values.

Hope this helps,


> Regards
> Vijay
> _______________________________________________
> AMBER mailing list

Daniel R. Roe, PhD
Department of Medicinal Chemistry
University of Utah
30 South 2000 East, Room 307
Salt Lake City, UT 84112-5820
(801) 587-9652
(801) 585-6208 (Fax)
AMBER mailing list
Received on Wed Aug 20 2014 - 08:00:03 PDT
Custom Search