MechStiff tutorial¶

This is an example how to use Mechanical Stiffness Calculations implemented in ProDy package.

The theory of MechStiff has been describe in:

Eyal E, Bahar I (2008) Toward a molecular understanding of the anisotropic response of proteins to external forces: insights from elastic network models Biophys J 94: 3424-35. PMID: 18223005; PMC2292382

http://www.ncbi.nlm.nih.gov/pubmed/18223005?ordinalpos=5&itool=EntrezSystem2.PEntrez.Pubmed.Pubmed_ResultsPanel.Pubmed_RVDocSum

In [1]:
%matplotlib inline
from prody import *
import matplotlib.pylab as plt
In [2]:
gfp, header = parsePDB('1gfl', header=True)
@> PDB file is found in the local folder (/data/prody_tutorials/.../1gfl.pdb.gz).
@> 3950 atoms and 1 coordinate set(s) were parsed in 0.06s.
In [3]:
gfp
Out[3]:
<AtomGroup: 1gfl (3950 atoms)>
In [4]:
calphas = gfp.select('protein and chain A and name CA')
In [5]:
calphas
Out[5]:
<Selection: 'protein and chain A and name CA' from 1gfl (230 atoms)>

Building Hessian using ANM model¶

We instantiate an ANM instance and we are building Hessian matrix by passing selected atoms.

In [6]:
anm = ANM('GFP ANM analysis')
In [7]:
anm.buildHessian(calphas, cutoff=13.0)
@> Hessian was built in 0.08s.
In [8]:
anm.getHessian().round(3)
Out[8]:
array([[ 9.017,  2.741, -0.979, ...,  0.   ,  0.   ,  0.   ],
       [ 2.741,  5.32 ,  1.181, ...,  0.   ,  0.   ,  0.   ],
       [-0.979,  1.181,  5.663, ...,  0.   ,  0.   ,  0.   ],
       ..., 
       [ 0.   ,  0.   ,  0.   , ...,  4.794, -0.052, -3.557],
       [ 0.   ,  0.   ,  0.   , ..., -0.052,  3.885,  1.827],
       [ 0.   ,  0.   ,  0.   , ..., -3.557,  1.827,  7.321]])

Mechanical Stiffness Matrix Calculations¶

We can perfom MechStiff calculations for selected group of 230 Ca's atoms.

In [9]:
anm.buildMechStiff(calphas)
@> 690 modes were calculated in 0.45s.
@> Calculating stiffness matrix.
@> Stiffness matrix calculated in 0.57s.
@> The range of effective force constant is: 5.27141824578 to 17.3849094465.
In [10]:
anm.getStiffness()
Out[10]:
array([[ 0.        ,  8.54744466,  8.55133411, ...,  8.73270017,
         8.36529536,  7.42692504],
       [ 8.54744466,  0.        ,  8.89990284, ...,  8.98097054,
         8.61849698,  7.77257692],
       [ 8.55133411,  8.89990284,  0.        , ...,  8.87863489,
         8.55675667,  7.84521673],
       ..., 
       [ 8.73270017,  8.98097054,  8.87863489, ...,  0.        ,
         8.87731778,  8.36590599],
       [ 8.36529536,  8.61849698,  8.55675667, ...,  8.87731778,
         0.        ,  9.52997339],
       [ 7.42692504,  7.77257692,  7.84521673, ...,  8.36590599,
         9.52997339,  0.        ]])

To obtain matrix with effective spring constant values use:

In [11]:
showMechStiff(anm, calphas, 'jet_r')
@> 690 modes were calculated in 0.62s.
@> Calculating stiffness matrix.
@> Stiffness matrix calculated in 0.24s.
@> The range of effective force constant is: 5.27141824578 to 17.3849094465.
Out[11]:
(<matplotlib.image.AxesImage at 0x7f7be0830590>,
 <matplotlib.colorbar.Colorbar at 0x7f7bc988ba90>)