4.7.4. Dihedral and Ramachandran analysis — MDAnalysis.analysis.dihedrals

Author:Henry Mull
Year:2018
Copyright:GNU Public License v2

New in version 0.18.1.

This module calculates the dihedral angles phi and psi in degrees for a given list of residues (or atoms corresponding to the residues). This can be done for selected frames or whole trajectories.

A list of time steps that contain phi and psi angles for each residue is generated, and a basic Ramachandran plot can be generated using the method Ramachandran.plot(). This plot is best used as a reference, but it also allows for user customization.

See also

MDAnalysis.lib.distances.calc_dihedrals()
function to calculate dihedral angles from atom positions

4.7.4.1. Example application

This example will show how to calculate the phi and psi angles of adenylate kinase (AdK) and generate a basic Ramachandran plot. The trajectory is included within the test data files:

import MDAnalysis as mda
from MDAnalysisTests.datafiles import GRO, XTC
u = mda.Universe(GRO, XTC)
r = u.select_atoms("protein")    # selection of residues

from MDAnalysis.analysis.dihedrals import Ramachandran
R = Ramachandran(r).run()

import matplotlib.pyplot as plt
fig, ax = plt.subplots(figsize=plt.figaspect(1))
ax.set_title("Ramachandran Plot (AdK)")
R.plot(ax=ax, color='k', marker='s')

Alternatively, if you wanted to plot the data yourself, the angles themselves can be accessed using Ramachandran.angles:

fig, ax = plt.subplots(figsize=plt.figaspect(1))
ax.axis([-180,180,-180,180])
ax.axhline(0, color='k', lw=1)
ax.axvline(0, color='k', lw=1)
ax.set(xticks=range(-180,181,60), yticks=range(-180,181,60),
       xlabel=r"$\phi$ (deg)", ylabel=r"$\psi$ (deg)")
for ts in R.angles:
    ax.scatter(ts[:,0], ts[:,1], color='k', marker='s')

4.7.4.2. Analysis Class

class MDAnalysis.analysis.dihedrals.Ramachandran(atomgroup, **kwargs)[source]

Calculate phi and psi dihedral angles of selected residues.

Phi and psi angles will be calculated for each residue corresponding to atomgroup for each time step in the trajectory. A ResidueGroup is generated from atomgroup which is compared to the protein to determine if it is a legitimate selection.

Note

If the residue selection is beyond the scope of the protein, then an error will be raised. If the residue selection includes the first or last residue, then a warning will be raised and they will be removed from the list of residues, but the analysis will still run.

Parameters:atomgroup (AtomGroup or ResidueGroup) – atoms for residues for which phi and psi are calculated
Raises:ValueError – If the selection of residues is not contained within the protein
angles

Contains the time steps of the phi and psi angles for each residue as an n_frames×n_residues×2 numpy.ndarray with content [[[phi, psi], [residue 2], ...], [time step 2], ...].

plot(ax=None, **kwargs)[source]

Plots data into standard ramachandran plot. Each time step in Ramachandran.angles is plotted onto the same graph.

Parameters:ax (matplotlib.axes.Axes) – If no ax is supplied or set to None then the plot will be added to the current active axes.
Returns:ax – Axes with the plot, either ax or the current axes.
Return type:matplotlib.axes.Axes
run()

Perform the calculation