AMBER Archive (2004)

Subject: RE: AMBER: what is averaged RMS fluctuation?

From: Ross Walker (ross_at_rosswalker.co.uk)
Date: Mon Jul 12 2004 - 15:40:23 CDT

 

> I still have one more question. The purpose of this RMS
> calculation is to
> find out which residue fluctuates most during the MD simulation. So my
> original idea was that: "x axis is the residue number; y axis is RMS
> fluctuation (bigger value means more mobile residue)". And I wanted to
> calculate the fluctuation like this: sqrt(<(r_i(k)-<r_i>)^2>).

In principle yes, but you won't be able to plot this on a single graph since
you have 3 variables and not 2. (residue number, rms fluctuation and time).

You you can either generate a graph for each point in time, a 3D graph that
may be difficult to understand (or may not, not sure on that) or a line
graph that has a different line for each residue and time on the x axis and
rmsd on the y axis (this is probably good for about 10 residues, then it
will start to look like a bit of a mess).

Now, what you suggest above is indeed a RMSd calculation that either ptraj
or carnal will happily do for you. Although be sure that you want the RMSd
per residue rather than atomic fluctuations. My understanding is that the
RMSd for a residue is the average RMSd for each atom in that residue, but I
may be mistaken. In my original suggestion r_i would have been the initial
structure rather than the average structure but I guess for fluctuations
rmsd from the average structure is more informative.
 
> I want to know if there is a better way to calculate this. Or
> is there any
> program can do this kind of calculation (like ptraj)?

Ptraj and carnal will. Probably best to do it in two steps though so you can
be sure you are doing what you think you are doing. First of all calculate
your average structure:

trajin job1.mdcrd
average average_struc.pdb pdb nowrap

Then use the average structure as the structure you want to perform the rmsd
against. I haven't tried doing a fit residue by residue in ptraj so I don't
have a script handy for this. Here is an example for carnal though:

FILES_IN
  PARM p1 1LDY_Charmm_Pavelites_NICH_full_params.prmtop;
  STATIC t1 1LDY_Charmm_Pavelites_NICH_full_params.prmcrd;
  STREAM s1 ../1LDY_Pavelites_NICH_Full_prod_300K_10ps.mdcrd.bz2;

FILES_OUT
  TABLE tab1 1LDY_pav_full_rms_per_res_1_to_756.tab;
  
DECLARE
  GROUP solute (RES 1-756);
  RMS fit1 FIT solute s1 t1;
  
OUTPUT
  TABLE tab1 fit1%residues;
 
END

All the best
Ross

/\
\/
|\oss Walker

| Department of Molecular Biology TPC15 |
| The Scripps Research Institute |
| Tel:- +1 858 784 8889 | EMail:- ross_at_rosswalker.co.uk |
| http://www.rosswalker.co.uk/ | PGP Key available on request |

-----------------------------------------------------------------------
The AMBER Mail Reflector
To post, send mail to amber_at_scripps.edu
To unsubscribe, send "unsubscribe amber" to majordomo_at_scripps.edu