Author: | Oliver Beckstein |
---|---|
Year: | 2012 |
Copyright: | GNU Public License v2 |
New in version 0.7.7.
The module contains code to analyze root mean square quantities such as the RMSD or RMSF (not implemented yet).
This module uses the fast QCP algorithm [Theobald2005] to calculate the root mean square distance (RMSD) between two coordinate sets (as implemented in MDAnalysis.core.qcprot.CalcRMSDRotationalMatrix()).
When using this module in published work please cite [Theobald2005].
See also
In this example we will globally fit a protein to a reference structure and investigate the relative movements of domains by computing the RMSD of the domains to the reference. The example is a DIMS trajectory of adenylate kinase, which samples a large closed-to-open transition. The protein consists of the CORE, LID, and NMP domain.
The trajectory is included with the test data files. The data in RMSD.rmsd is plotted with matplotlib.pyplot.plot():
import MDAnalysis
from MDAnalysis.tests.datafiles import PSF,DCD,CRD
u = MDAnalysis.Universe(PSF,DCD)
ref = MDAnalysis.Universe(PSF,DCD) # reference closed AdK (1AKE) (with the default ref_frame=0)
#ref = MDAnalysis.Universe(PSF,CRD) # reference open AdK (4AKE)
import MDAnalysis.analysis.rms
R = MDAnalysis.analysis.rms.RMSD(u, ref,
select="backbone", # superimpose on whole backbone of the whole protein
groupselections=["backbone and (resid 1-29 or resid 60-121 or resid 160-214)", # CORE
"backbone and resid 122-159", # LID
"backbone and resid 30-59"], # NMP
filename="rmsd_all_CORE_LID_NMP.dat")
R.run()
R.save()
import matplotlib.pyplot as plt
rmsd = R.rmsd.T # transpose makes it easier for plotting
time = rmsd[1]
fig = plt.figure(figsize=(4,4))
ax = fig.add_subplot(111)
ax.plot(time, rmsd[2], 'k-', label="all")
ax.plot(time, rmsd[3], 'k--', label="CORE")
ax.plot(time, rmsd[4], 'r--', label="LID")
ax.plot(time, rmsd[5], 'b--', label="NMP")
ax.legend(loc="best")
ax.set_xlabel("time (ps)")
ax.set_ylabel(r"RMSD ($\AA$)")
fig.savefig("rmsd_all_CORE_LID_NMP_ref1AKE.pdf")
Returns RMSD between two coordinate sets a and b.
a and b are arrays of the coordinates of N atoms of shape N*3 as generated by, e.g., MDAnalysis.core.AtomGroup.AtomGroup.coordinates().
An implicit optimal superposition is performed, which minimizes the RMSD between a and b although both a and b must be centered on the origin before performing the RMSD calculation so that translations are removed.
One can use the center = True keyword, which subtracts the center of geometry (for weights = None) before the superposition. With weights, a weighted average is computed as the center.
The weights can be an array of length N, containing e.g. masses for a weighted RMSD calculation.
The function uses Douglas Theobald’s fast QCP algorithm [Theobald2005] to calculate the RMSD.
>>> u = Universe(PSF,DCD)
>>> bb = u.selectAtoms('backbone')
>>> A = bb.coordinates() # coordinates of first frame
>>> u.trajectory[-1] # forward to last frame
>>> B = bb.coordinates() # coordinates of last frame
>>> rmsd(A,B)
6.8342494129169804
Class to perform RMSD analysis on a trajectory.
Run the analysis with RMSD.run(), which stores the results in the array RMSD.rmsd:
frame time (ps) RMSD (A)
This class uses Douglas Theobald’s fast QCP algorithm [Theobald2005] to calculate the RMSD.
New in version 0.7.7.
Setting up the RMSD analysis.
The RMSD will be computed between select and reference for all frames in the trajectory in universe.
Arguments: |
|
---|
New in version 0.7.7.
Changed in version 0.8: groupselections added
Results are stored in this N×3 numpy.ndarray array, (frame, time (ps), RMSD (Å)).
Perform RMSD analysis on the trajectory.
A number of parameters can be changed from the defaults. The result is stored as the array RMSD.rmsd.
Keywords: |
|
---|
Save RMSD from RMSD.rmsd to text file filename.
If filename is not supplied then the default provided to the constructor is used.
The data are saved with numpy.savetxt().