Calculate dihedral angles from xtc trajectory

[xiki-tempula]

Expected behaviour

I'm interested in calculating the time dependent change of a particular dihedral angle from a xtc trajectory and I expect to get a list of angles. I copied the example from https://github.com/MDAnalysis/MDAnalysisCookbook/blob/master/examples/backbone_dihedral.py

Actual behaviour

File "./dihedral_analysis.py", line 66, in
collection.compute(u.trajectory)

File "./mdanalysis/package/MDAnalysis/core/Timeseries.py", line 139, in compute
self.data = trj.correl(self, start, stop, step)

AttributeError: 'XTCReader' object has no attribute 'correl'

Code to reproduce the behaviour

from MDAnalysis import Universe, collection, Timeseries
collection.clear()
u = Universe(gro, xtc)
selection = 'resname {} and resid {} and name {}'
ALA2N = selection.format('ALA', 2, 'N')
ALA2CA = selection.format('ALA', 2, 'CA')
ALA2C = selection.format('ALA', 2, 'C')
TYR3N = selection.format('TYR', 3, 'N')
psi = [ALA2N, ALA2CA, ALA2C, TYR3N]
psi_angle = [u.select_atoms(atom) for atom in psi]
collection.addTimeseries(Timeseries.Dihedral(psi_angle))
collection.compute(u.trajectory)

I'm not quite sure why you need to call this correl function? I guess to calculate the dihedral angle, you only need four cartesian coordinates?
Something like

psi_list = []
for ts in u.trajectory:
    atoms = [u.select_atoms(atom) for atom in psi]
    cords = [atom.positions[0] for atom in atoms]
    psi_list.append(dihedral(cords))

Currently version of MDAnalysis:

py35, 0.16.0

[richardjgowers]

I've never really used timeseries myself, I think the following should work..

from MDAnalysis import Universe

u = Universe(gro, xtc)
selection = 'resname {} and resid {} and name {}'
ALA2N = selection.format('ALA', 2, 'N')
ALA2CA = selection.format('ALA', 2, 'CA')
ALA2C = selection.format('ALA', 2, 'C')
TYR3N = selection.format('TYR', 3, 'N')
psi = [ALA2N, ALA2CA, ALA2C, TYR3N]
psi_angle = sum([u.select_atoms(atom) for atom in psi])  # sum of Atoms creates an AtomGroup
psi_angle = psi_angle.dihedral  # convert AtomGroup to Dihedral object

psi_list = []
for ts in u.trajectory:
    psi_list.append(psi_angle.value())

• EDIT (@orbeckst): removed any reference to Timeseries or correl because it is not needed for the solution. Just to make sure that anyone else who might find this solution later will not become confused.

[xiki-tempula]

@richardjgowers Thank you. It solved my problem.

[kain88-de]

The timeseries feature only works for dcd files or a trajectory copied into memory.

[orbeckst]

As @kain88-de says, correl is not generally implemented. For DCD it contains code at the C-level to do simple trajectory analysis. You basically describe all the analysis that you want to do in your correl collection (the whole thing was inspired by CHARMM's CORREL functionality). Then you compute all your observables in one sweep with correl.compute(). In principle this is a good way to stack analysis and amortize the high I/O cost per frame but it relies on the reader implementing the functionality, which is not very flexible.

Typically we recommend using AnalysisBase to write general analysis functions. At some point we might figure out how to "stack" multiple analysis classes and basically do the same thing that correl does, just in a general way, and at the Python level.