Alanine RMSF analysis¶

Uses MDAnalysis to analyze the alanine trajectory produced by docs/episodes/scripts/03-simulaciones-clasicas_simple.py. This mirrors the scripts in https://github.com/JordiVillaFreixa/Pau_TFG_DAO/tree/main/analysis.

Table of contents¶

  • Alanine RMSF analysis
  • Alanine RMSF analysis
  • Step 1

Step 1¶

In [ ]:
from pathlib import Path
import os
import MDAnalysis as mda
from MDAnalysis.analysis import align
import pandas as pd
import numpy as np
import matplotlib.pyplot as plt

COURSE_DIR = Path(os.environ.get('COURSE_DIR','~/Concepcion26')).expanduser()
topo = COURSE_DIR / 'data' / 'alanine-dipeptide.pdb'
traj = COURSE_DIR / 'results' / '03-simulaciones-clasicas' / 'simple' / 'traj.dcd'
if not topo.exists() or not traj.exists():
    raise FileNotFoundError(f'Missing segments - run 03-simulaciones-clasicas_simple.py first')

u=mda.Universe(topo,traj)
sel=u.select_atoms('protein and name CA')
average=align.AverageStructure(u,u,select='protein and name CA',ref_frame=0).run().results.universe
align.AlignTraj(u,average,select='protein and name CA',in_memory=True).run()
coords=[]
for ts in u.trajectory[100:]:
    coords.append(sel.positions.copy())
coords=np.array(coords)
avg=coords.mean(axis=0)
rmsf=np.sqrt(((coords-avg)**2).sum(axis=2).mean(axis=0))
residues=[res.resid for res in sel.residues]
df=pd.DataFrame({'Residue ID':residues,'RMSF (Å)':rmsf})
fig,ax=plt.subplots()
ax.plot(df['Residue ID'],df['RMSF (Å)'])
ax.set_xlabel('Residue')
ax.set_ylabel('RMSF (Å)')
ax.set_title('Alanine Cα RMSF')
ax.grid(True)
fig.tight_layout()
fig.savefig('index_rmsf.png')
print(df.head())