Want to make distance vs. time and torsion vs. time graphs from MD trajectory

[Andrew Kleist]

Hello!

I am interesting in doing two analyses with MD trajectories with Bio3d:

1) Distance between two atoms versus simulation time: specifically, I want to measure the distance between two helices of a GPCR over time using the alpha carbon of residues at the base of each helix 2) Torsion angle fluctuation versus simulation time; specifically with chi2 angles of Trp, for instance: fluctuations in a core Trp are important for activation of my GPCR, so I want to use the chi2 angle fluctuation as a proxy for activation state.

It looks like dist.xyz can be used to measure distances, but I'm confused on the implementation. In the index example: dist.xyz( c(1,1,1, 3,3,3), c(3,3,3, 2,2,2, 1,1,1)) ...what do the two vectors represent? I'm hoping what I can do is do dist.xyz(atom1, atom2), but this doesn't seem to work. Also I want the measurement to be iterated for all MD frames, and I'm not sure how to implement this.

As for the torsion angles, there is a torsion.xyz fucntion, but I have the same question: where do I indicate the atom numbers for the torsion angle, and how do I iterate the measurement over all trajectory frames? In addition to single atom analysis, I also want to do a cross-correlation matrix analysis as in the trajectory analysis tutorial except use torsion angles instead of calpha rmsd.

I am new to R, so I apologize if these are basic questions! Thanks for the help!

-Andy

[Xinqiu Yao]

Hi,

The first two arguments of dist.xyz() are 'xyz' coordinates of atoms, represented as a vector in a form like 'atom1_x, atom1_y, atom1_z, atom2_x, atom2_y, atom2_z, ...'. When input are matrices, rows of each matrix represent coordinates of the same atom in different frames. In this case, the number of columns between 'a' and 'b' must be the same (usually 3). For your purpose, your input 'a' should be a matrix of N-by-3 containing the coordinates of atom1, where N is the number of frames, and input 'b' is the same format as 'a' but for atom2. You should also set all.pairs=FALSE to calculate distance in the same frame only.

For torsion angles, just provide 'xyz' coordinates of four successive atoms for which you want to calculate the angle. Note that you can use the function atom.select and your topology file of the MD (e.g. a pdb file) to locate specific atoms in the 'xyz' matrix. (You might also want to check the torsion.pdb() function to calculate pre-defined angles for a whole protein).

For example,

# Suppose your topology and trajectory are 'pdb' and 'xyz', respectively.
#
# To calculate the distance between C-alpha atoms of residues #1 and #100.
atom1 <- atom.select(pdb, resno=1, elety='CA')
atom2 <- atom.select(pdb, resno=100, elety='CA')
d <- dist.xyz(xyz[, atom1$xyz], xyz[, atom2$xyz], all.pairs=FALSE)

# To calculate the torsion angle defined by the first four successive atoms.
atoms <- atom.select(pdb, eleno=1:4)
tor <- torsion.xyz(xyz[, atoms$xyz])

Not sure about torsional cross-correlations. I guess you could directly use the R function cor(), taking a matrix of angle values as input where rows of the matrix representing frames and columns the angles. You have to take care of the periodic issue of angles probably...

[Andrew Kleist]

Great, thanks for the quick response! I'll give this a try and let you know if I run into any other problems.