# Anisotropic Network Model (ANM)¶

This example shows how to perform ANM calculations, and retrieve
normal mode data. An `ANM`

instance that stores Hessian matrix
(and also Kirchhoff matrix) and normal mode data describing the intrinsic
dynamics of the protein structure will be obtained. `ANM`

instances
and individual normal modes (`Mode`

) can be used as input to functions
in `dynamics`

module.

See [Doruker00] and [Atilgan01] for more information on the theory of ANM.

[Doruker00] | Doruker P, Atilgan AR, Bahar I. Dynamics of proteins predicted
by molecular dynamics simulations and analytical approaches: Application to
a-amylase inhibitor. Proteins 2000 40:512-524. |

[Atilgan01] | Atilgan AR, Durrell SR, Jernigan RL, Demirel MC, Keskin O,
Bahar I. Anisotropy of fluctuation dynamics of proteins with an
elastic network model. Biophys. J. 2001 80:505-515. |

## Parse structure¶

We start by importing everything from the ProDy package:

```
In [1]: from prody import *
In [2]: from pylab import *
In [3]: ion()
```

We start with parsing a PDB file by passing an identifier. Note that if a file is not found in the current working directory, it will be downloaded.

```
In [4]: p38 = parsePDB('1p38')
In [5]: p38
Out[5]: <AtomGroup: 1p38 (2962 atoms)>
```

We want to use only Cα atoms, so we select them:

```
In [6]: calphas = p38.select('protein and name CA')
In [7]: calphas
Out[7]: <Selection: 'protein and name CA' from 1p38 (351 atoms)>
```

We can also make the same selection like this:

```
In [8]: calphas2 = p38.select('calpha')
In [9]: calphas2
Out[9]: <Selection: 'calpha' from 1p38 (351 atoms)>
```

To check whether the selections are the same, we can try:

```
In [10]: calphas == calphas2
Out[10]: True
```

Note that, ProDy atom selector gives the flexibility to select any set of atoms to be used in ANM calculations.

## Build Hessian¶

We instantiate an `ANM`

instance:

```
In [11]: anm = ANM('p38 ANM analysis')
```

Then, build the Hessian matrix by passing selected atoms (351 Cα’s)
to `ANM.buildHessian()`

method:

```
In [12]: anm.buildHessian(calphas)
```

We can get a copy of the Hessian matrix using `ANM.getHessian()`

method:

```
In [13]: anm.getHessian().round(3)
Out[13]:
array([[ 9.959, -3.788, 0.624, ..., 0. , 0. , 0. ],
[-3.788, 7.581, 1.051, ..., 0. , 0. , 0. ],
[ 0.624, 1.051, 5.46 , ..., 0. , 0. , 0. ],
...,
[ 0. , 0. , 0. , ..., 1.002, -0.282, 0.607],
[ 0. , 0. , 0. , ..., -0.282, 3.785, -2.504],
[ 0. , 0. , 0. , ..., 0.607, -2.504, 4.214]])
```

## Parameters¶

We didn’t pass any parameters to `ANM.buildHessian()`

method, but it
accepts *cutoff* and *gamma* parameters, for which default values are
`cutoff=15.0`

and `gamma=1.0`

.

```
In [14]: anm.getCutoff()
Out[14]: 15.0
In [15]: anm.getGamma()
Out[15]: 1.0
```

Note that it is also possible to use an externally calculated Hessian
matrix. Just pass it to the ANM instance using `ANM.setHessian()`

method.

## Calculate normal modes¶

Calculate modes using `ANM.calcModes()`

method:

```
In [16]: anm.calcModes()
```

Note that by default 20 non-zero (or non-trivial) and 6 trivial modes are
calculated. Trivial modes are not retained. To calculate a different number
of non-zero modes or to keep zero modes, try `anm.calcModes(50, zeros=True)`

.

## Normal modes data¶

```
In [17]: anm.getEigvals().round(3)
Out[17]:
array([0.179, 0.334, 0.346, 0.791, 0.942, 1.012, 1.188, 1.304, 1.469,
1.546, 1.608, 1.811, 1.925, 1.983, 2.14 , 2.298, 2.33 , 2.364,
2.69 , 2.794])
In [18]: anm.getEigvecs().round(3)
Out[18]:
array([[ 0.039, -0.045, 0.007, ..., 0.105, 0.032, -0.038],
[ 0.009, -0.096, -0.044, ..., 0.091, 0.036, -0.037],
[ 0.058, -0.009, 0.08 , ..., -0.188, -0.08 , -0.063],
...,
[ 0.046, -0.093, -0.131, ..., 0.018, -0.008, 0.006],
[ 0.042, -0.018, -0.023, ..., 0.014, -0.043, 0.037],
[ 0.08 , -0.002, -0.023, ..., 0.024, -0.023, -0.009]])
```

You can get the covariance matrix as follows:

```
In [19]: anm.getCovariance().round(2)
Out[19]:
array([[ 0.03, 0.03, -0. , ..., 0. , 0. , 0.01],
[ 0.03, 0.06, -0.03, ..., 0.01, -0. , 0.01],
[-0. , -0.03, 0.09, ..., -0.01, -0. , 0.01],
...,
[ 0. , 0.01, -0.01, ..., 1.21, 0. , -0.17],
[ 0. , -0. , -0. , ..., 0. , 0.41, 0.38],
[ 0.01, 0.01, 0.01, ..., -0.17, 0.38, 0.4 ]])
```

Covariance matrices are calculated using the available modes (slowest 20 modes in this case). If the user calculates M slowest modes, only they will be used in the calculation of covariances.

## Individual modes¶

Normal mode indices in Python start from 0, so the slowest mode has index 0. By default, modes with zero eigenvalues are excluded. If they were retained, the slowest non-trivial mode would have index 6.

Get the slowest mode by indexing `ANM`

instance as follows:

```
In [20]: slowest_mode = anm[0]
In [21]: slowest_mode.getEigval().round(3)
Out[21]: 0.179
In [22]: slowest_mode.getEigvec().round(3)
Out[22]: array([0.039, 0.009, 0.058, ..., 0.046, 0.042, 0.08 ])
```

## Write NMD file¶

ANM results in NMD format can be visualized using Normal Mode Wizard VMD plugin. The following statement writes the slowest 3 ANM modes into an NMD file:

```
In [23]: writeNMD('p38_anm_modes.nmd', anm[:3], calphas)
Out[23]: 'p38_anm_modes.nmd'
```

Note that slicing an `ANM`

objects returns a list of modes.
In this case, slowest 3 ANM modes were written into NMD file.

## View modes in VMD¶

First make sure that the VMD path is correct

```
In [24]: pathVMD()
Out[24]: '/usr/local/bin/vmd'
```

```
# if this is incorrect use setVMDpath to correct it
In [25]: viewNMDinVMD('p38_anm_modes.nmd')
```

This will show the slowest 3 modes in VMD using NMWiz. This concludes the ANM
example. Many of the methods demonstrated here apply to other NMA models, such
as `GNM`

and `EDA`

.

## Advanced visualization in jupyter notebooks¶

You can visualize structures and modes determined from ANM or GNM calculations in jupyter notebooks using another python module, py3Dmol. It is a java-script library that can visualize structural elements with light weight customization.

You can find an example notebook.