Untitled_amber02.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 2,
   "id": "a045bbb8",
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import matplotlib.pyplot as plt\n",
    "import MDAnalysis as mda\n",
    "from MDAnalysis.analysis import align\n",
    "from MDAnalysis.analysis.rms import RMSD"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "40b0636e",
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "id": "64d385d6",
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "df54129d",
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "8105639b",
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "447bacd1",
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "8fa2f6cc",
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "id": "80e01dc3",
   "metadata": {},
   "outputs": [],
   "source": [
    "\n",
    "\n",
    "# Load the simulation data\n",
    "u = mda.Universe('cpx_solvated.prmtop', 'prod.nc')\n",
    "\n",
    "# Select the protein atoms\n",
    "protein = u.select_atoms('protein')\n",
    "\n",
    "# Calculate the reference structure (the first frame of the trajectory)\n",
    "ref_structure = protein.positions.copy()\n",
    "\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "id": "9d425921",
   "metadata": {},
   "outputs": [],
   "source": [
    "# Initialize the RMSD analysis object\n",
    "R = RMSD(u, select='protein', ref_frame=0)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "id": "1ca4c18b",
   "metadata": {},
   "outputs": [
    {
     "ename": "TypeError",
     "evalue": "'NoneType' object is not subscriptable",
     "output_type": "error",
     "traceback": [
      "\u001b[0;31m---------------------------------------------------------------------------\u001b[0m",
      "\u001b[0;31mTypeError\u001b[0m                                 Traceback (most recent call last)",
      "Cell \u001b[0;32mIn[15], line 1\u001b[0m\n\u001b[0;32m----> 1\u001b[0m \u001b[43mR\u001b[49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43mrun\u001b[49m\u001b[43m(\u001b[49m\u001b[43m)\u001b[49m\n",
      "File \u001b[0;32m~/anaconda3/lib/python3.10/site-packages/MDAnalysis/analysis/base.py:440\u001b[0m, in \u001b[0;36mAnalysisBase.run\u001b[0;34m(self, start, stop, step, frames, verbose, progressbar_kwargs)\u001b[0m\n\u001b[1;32m    436\u001b[0m \u001b[38;5;28mself\u001b[39m\u001b[38;5;241m.\u001b[39m_prepare()\n\u001b[1;32m    437\u001b[0m logger\u001b[38;5;241m.\u001b[39minfo(\u001b[38;5;124m\"\u001b[39m\u001b[38;5;124mStarting analysis loop over \u001b[39m\u001b[38;5;132;01m%d\u001b[39;00m\u001b[38;5;124m trajectory frames\u001b[39m\u001b[38;5;124m\"\u001b[39m,\n\u001b[1;32m    438\u001b[0m             \u001b[38;5;28mself\u001b[39m\u001b[38;5;241m.\u001b[39mn_frames)\n\u001b[0;32m--> 440\u001b[0m \u001b[38;5;28;01mfor\u001b[39;00m i, ts \u001b[38;5;129;01min\u001b[39;00m \u001b[38;5;28menumerate\u001b[39m(ProgressBar(\n\u001b[1;32m    441\u001b[0m         \u001b[38;5;28mself\u001b[39m\u001b[38;5;241m.\u001b[39m_sliced_trajectory,\n\u001b[1;32m    442\u001b[0m         verbose\u001b[38;5;241m=\u001b[39mverbose,\n\u001b[1;32m    443\u001b[0m         \u001b[38;5;241m*\u001b[39m\u001b[38;5;241m*\u001b[39mprogressbar_kwargs)):\n\u001b[1;32m    444\u001b[0m     \u001b[38;5;28mself\u001b[39m\u001b[38;5;241m.\u001b[39m_frame_index \u001b[38;5;241m=\u001b[39m i\n\u001b[1;32m    445\u001b[0m     \u001b[38;5;28mself\u001b[39m\u001b[38;5;241m.\u001b[39m_ts \u001b[38;5;241m=\u001b[39m ts\n",
      "File \u001b[0;32m~/anaconda3/lib/python3.10/site-packages/tqdm/notebook.py:259\u001b[0m, in \u001b[0;36mtqdm_notebook.__iter__\u001b[0;34m(self)\u001b[0m\n\u001b[1;32m    257\u001b[0m \u001b[38;5;28;01mtry\u001b[39;00m:\n\u001b[1;32m    258\u001b[0m     it \u001b[38;5;241m=\u001b[39m \u001b[38;5;28msuper\u001b[39m(tqdm_notebook, \u001b[38;5;28mself\u001b[39m)\u001b[38;5;241m.\u001b[39m\u001b[38;5;21m__iter__\u001b[39m()\n\u001b[0;32m--> 259\u001b[0m     \u001b[38;5;28;01mfor\u001b[39;00m obj \u001b[38;5;129;01min\u001b[39;00m it:\n\u001b[1;32m    260\u001b[0m         \u001b[38;5;66;03m# return super(tqdm...) will not catch exception\u001b[39;00m\n\u001b[1;32m    261\u001b[0m         \u001b[38;5;28;01myield\u001b[39;00m obj\n\u001b[1;32m    262\u001b[0m \u001b[38;5;66;03m# NB: except ... [ as ...] breaks IPython async KeyboardInterrupt\u001b[39;00m\n",
      "File \u001b[0;32m~/anaconda3/lib/python3.10/site-packages/tqdm/std.py:1183\u001b[0m, in \u001b[0;36mtqdm.__iter__\u001b[0;34m(self)\u001b[0m\n\u001b[1;32m   1180\u001b[0m \u001b[38;5;66;03m# If the bar is disabled, then just walk the iterable\u001b[39;00m\n\u001b[1;32m   1181\u001b[0m \u001b[38;5;66;03m# (note: keep this check outside the loop for performance)\u001b[39;00m\n\u001b[1;32m   1182\u001b[0m \u001b[38;5;28;01mif\u001b[39;00m \u001b[38;5;28mself\u001b[39m\u001b[38;5;241m.\u001b[39mdisable:\n\u001b[0;32m-> 1183\u001b[0m     \u001b[38;5;28;01mfor\u001b[39;00m obj \u001b[38;5;129;01min\u001b[39;00m iterable:\n\u001b[1;32m   1184\u001b[0m         \u001b[38;5;28;01myield\u001b[39;00m obj\n\u001b[1;32m   1185\u001b[0m     \u001b[38;5;28;01mreturn\u001b[39;00m\n",
      "File \u001b[0;32m~/anaconda3/lib/python3.10/site-packages/MDAnalysis/coordinates/base.py:714\u001b[0m, in \u001b[0;36mProtoReader.__next__\u001b[0;34m(self)\u001b[0m\n\u001b[1;32m    712\u001b[0m \u001b[38;5;28;01mdef\u001b[39;00m \u001b[38;5;21m__next__\u001b[39m(\u001b[38;5;28mself\u001b[39m) \u001b[38;5;241m-\u001b[39m\u001b[38;5;241m>\u001b[39m Timestep:\n\u001b[1;32m    713\u001b[0m     \u001b[38;5;124;03m\"\"\"Forward one step to next frame when using the `next` builtin.\"\"\"\u001b[39;00m\n\u001b[0;32m--> 714\u001b[0m     \u001b[38;5;28;01mreturn\u001b[39;00m \u001b[38;5;28;43mself\u001b[39;49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43mnext\u001b[49m\u001b[43m(\u001b[49m\u001b[43m)\u001b[49m\n",
      "File \u001b[0;32m~/anaconda3/lib/python3.10/site-packages/MDAnalysis/coordinates/base.py:700\u001b[0m, in \u001b[0;36mProtoReader.next\u001b[0;34m(self)\u001b[0m\n\u001b[1;32m    698\u001b[0m \u001b[38;5;124;03m\"\"\"Forward one step to next frame.\"\"\"\u001b[39;00m\n\u001b[1;32m    699\u001b[0m \u001b[38;5;28;01mtry\u001b[39;00m:\n\u001b[0;32m--> 700\u001b[0m     ts \u001b[38;5;241m=\u001b[39m \u001b[38;5;28;43mself\u001b[39;49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43m_read_next_timestep\u001b[49m\u001b[43m(\u001b[49m\u001b[43m)\u001b[49m\n\u001b[1;32m    701\u001b[0m \u001b[38;5;28;01mexcept\u001b[39;00m (\u001b[38;5;167;01mEOFError\u001b[39;00m, \u001b[38;5;167;01mIOError\u001b[39;00m):\n\u001b[1;32m    702\u001b[0m     \u001b[38;5;28mself\u001b[39m\u001b[38;5;241m.\u001b[39mrewind()\n",
      "File \u001b[0;32m~/anaconda3/lib/python3.10/site-packages/MDAnalysis/coordinates/TRJ.py:687\u001b[0m, in \u001b[0;36mNCDFReader._read_next_timestep\u001b[0;34m(self, ts)\u001b[0m\n\u001b[1;32m    685\u001b[0m     ts \u001b[38;5;241m=\u001b[39m \u001b[38;5;28mself\u001b[39m\u001b[38;5;241m.\u001b[39mts\n\u001b[1;32m    686\u001b[0m \u001b[38;5;28;01mtry\u001b[39;00m:\n\u001b[0;32m--> 687\u001b[0m     \u001b[38;5;28;01mreturn\u001b[39;00m \u001b[38;5;28;43mself\u001b[39;49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43m_read_frame\u001b[49m\u001b[43m(\u001b[49m\u001b[38;5;28;43mself\u001b[39;49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43m_current_frame\u001b[49m\u001b[43m \u001b[49m\u001b[38;5;241;43m+\u001b[39;49m\u001b[43m \u001b[49m\u001b[38;5;241;43m1\u001b[39;49m\u001b[43m)\u001b[49m\n\u001b[1;32m    688\u001b[0m \u001b[38;5;28;01mexcept\u001b[39;00m \u001b[38;5;167;01mIndexError\u001b[39;00m:\n\u001b[1;32m    689\u001b[0m     \u001b[38;5;28;01mraise\u001b[39;00m \u001b[38;5;167;01mIOError\u001b[39;00m \u001b[38;5;28;01mfrom\u001b[39;00m \u001b[38;5;28mNone\u001b[39m\n",
      "File \u001b[0;32m~/anaconda3/lib/python3.10/site-packages/MDAnalysis/coordinates/TRJ.py:675\u001b[0m, in \u001b[0;36mNCDFReader._read_frame\u001b[0;34m(self, frame)\u001b[0m\n\u001b[1;32m    672\u001b[0m         \u001b[38;5;28mself\u001b[39m\u001b[38;5;241m.\u001b[39mconvert_forces_from_native(ts\u001b[38;5;241m.\u001b[39m_forces, inplace\u001b[38;5;241m=\u001b[39m\u001b[38;5;28;01mTrue\u001b[39;00m)\n\u001b[1;32m    673\u001b[0m     \u001b[38;5;28;01mif\u001b[39;00m \u001b[38;5;28mself\u001b[39m\u001b[38;5;241m.\u001b[39mperiodic:\n\u001b[1;32m    674\u001b[0m         \u001b[38;5;66;03m# in-place ! (only lengths)\u001b[39;00m\n\u001b[0;32m--> 675\u001b[0m         \u001b[38;5;28mself\u001b[39m\u001b[38;5;241m.\u001b[39mconvert_pos_from_native(\u001b[43mts\u001b[49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43mdimensions\u001b[49m\u001b[43m[\u001b[49m\u001b[43m:\u001b[49m\u001b[38;5;241;43m3\u001b[39;49m\u001b[43m]\u001b[49m)\n\u001b[1;32m    676\u001b[0m ts\u001b[38;5;241m.\u001b[39mframe \u001b[38;5;241m=\u001b[39m frame  \u001b[38;5;66;03m# frame labels are 0-based\u001b[39;00m\n\u001b[1;32m    677\u001b[0m \u001b[38;5;28mself\u001b[39m\u001b[38;5;241m.\u001b[39m_current_frame \u001b[38;5;241m=\u001b[39m frame\n",
      "\u001b[0;31mTypeError\u001b[0m: 'NoneType' object is not subscriptable"
     ]
    }
   ],
   "source": [
    "\n",
    "R.run()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "4e2d4563",
   "metadata": {},
   "outputs": [],
   "source": [
    "\n",
    "\n",
    "# Extract the RMSD values\n",
    "rmsd_values = R.rmsd[:, 2]\n",
    "frame_numbers = np.arange(len(rmsd_values)) / 100\n",
    "\n",
    "# Plot the RMSD values\n",
    "plt.figure(figsize=(15,5))\n",
    "\n",
    "plt.plot(frame_numbers,rmsd_values)\n",
    "plt.xlabel('Time(ns)')\n",
    "plt.ylabel('RMSD (Å)')\n",
    "plt.title('LTS0041891 and FDFT1 RMSD Plot')\n",
    "plt.savefig('rmsd_plot.png')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "id": "7bb1f677",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([1.24722680e-06, 9.24300016e-01, 1.03212485e+00, ...,\n",
       "       3.60878096e+00, 3.62158403e+00, 3.47682958e+00])"
      ]
     },
     "execution_count": 7,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "rmsd_values"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "id": "1ef21160",
   "metadata": {},
   "outputs": [
    {
     "ename": "TypeError",
     "evalue": "'NoneType' object is not subscriptable",
     "output_type": "error",
     "traceback": [
      "\u001b[0;31m---------------------------------------------------------------------------\u001b[0m",
      "\u001b[0;31mTypeError\u001b[0m                                 Traceback (most recent call last)",
      "Cell \u001b[0;32mIn[6], line 8\u001b[0m\n\u001b[1;32m      6\u001b[0m \u001b[38;5;66;03m# Initialize the RMSF analysis object\u001b[39;00m\n\u001b[1;32m      7\u001b[0m R \u001b[38;5;241m=\u001b[39m RMSF(calphas)\n\u001b[0;32m----> 8\u001b[0m \u001b[43mR\u001b[49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43mrun\u001b[49m\u001b[43m(\u001b[49m\u001b[43m)\u001b[49m\n\u001b[1;32m     10\u001b[0m \u001b[38;5;66;03m# Extract the RMSF values\u001b[39;00m\n\u001b[1;32m     11\u001b[0m rmsf_values \u001b[38;5;241m=\u001b[39m R\u001b[38;5;241m.\u001b[39mrmsf\n",
      "File \u001b[0;32m~/anaconda3/lib/python3.10/site-packages/MDAnalysis/analysis/base.py:440\u001b[0m, in \u001b[0;36mAnalysisBase.run\u001b[0;34m(self, start, stop, step, frames, verbose, progressbar_kwargs)\u001b[0m\n\u001b[1;32m    436\u001b[0m \u001b[38;5;28mself\u001b[39m\u001b[38;5;241m.\u001b[39m_prepare()\n\u001b[1;32m    437\u001b[0m logger\u001b[38;5;241m.\u001b[39minfo(\u001b[38;5;124m\"\u001b[39m\u001b[38;5;124mStarting analysis loop over \u001b[39m\u001b[38;5;132;01m%d\u001b[39;00m\u001b[38;5;124m trajectory frames\u001b[39m\u001b[38;5;124m\"\u001b[39m,\n\u001b[1;32m    438\u001b[0m             \u001b[38;5;28mself\u001b[39m\u001b[38;5;241m.\u001b[39mn_frames)\n\u001b[0;32m--> 440\u001b[0m \u001b[38;5;28;01mfor\u001b[39;00m i, ts \u001b[38;5;129;01min\u001b[39;00m \u001b[38;5;28menumerate\u001b[39m(ProgressBar(\n\u001b[1;32m    441\u001b[0m         \u001b[38;5;28mself\u001b[39m\u001b[38;5;241m.\u001b[39m_sliced_trajectory,\n\u001b[1;32m    442\u001b[0m         verbose\u001b[38;5;241m=\u001b[39mverbose,\n\u001b[1;32m    443\u001b[0m         \u001b[38;5;241m*\u001b[39m\u001b[38;5;241m*\u001b[39mprogressbar_kwargs)):\n\u001b[1;32m    444\u001b[0m     \u001b[38;5;28mself\u001b[39m\u001b[38;5;241m.\u001b[39m_frame_index \u001b[38;5;241m=\u001b[39m i\n\u001b[1;32m    445\u001b[0m     \u001b[38;5;28mself\u001b[39m\u001b[38;5;241m.\u001b[39m_ts \u001b[38;5;241m=\u001b[39m ts\n",
      "File \u001b[0;32m~/anaconda3/lib/python3.10/site-packages/tqdm/notebook.py:259\u001b[0m, in \u001b[0;36mtqdm_notebook.__iter__\u001b[0;34m(self)\u001b[0m\n\u001b[1;32m    257\u001b[0m \u001b[38;5;28;01mtry\u001b[39;00m:\n\u001b[1;32m    258\u001b[0m     it \u001b[38;5;241m=\u001b[39m \u001b[38;5;28msuper\u001b[39m(tqdm_notebook, \u001b[38;5;28mself\u001b[39m)\u001b[38;5;241m.\u001b[39m\u001b[38;5;21m__iter__\u001b[39m()\n\u001b[0;32m--> 259\u001b[0m     \u001b[38;5;28;01mfor\u001b[39;00m obj \u001b[38;5;129;01min\u001b[39;00m it:\n\u001b[1;32m    260\u001b[0m         \u001b[38;5;66;03m# return super(tqdm...) will not catch exception\u001b[39;00m\n\u001b[1;32m    261\u001b[0m         \u001b[38;5;28;01myield\u001b[39;00m obj\n\u001b[1;32m    262\u001b[0m \u001b[38;5;66;03m# NB: except ... [ as ...] breaks IPython async KeyboardInterrupt\u001b[39;00m\n",
      "File \u001b[0;32m~/anaconda3/lib/python3.10/site-packages/tqdm/std.py:1183\u001b[0m, in \u001b[0;36mtqdm.__iter__\u001b[0;34m(self)\u001b[0m\n\u001b[1;32m   1180\u001b[0m \u001b[38;5;66;03m# If the bar is disabled, then just walk the iterable\u001b[39;00m\n\u001b[1;32m   1181\u001b[0m \u001b[38;5;66;03m# (note: keep this check outside the loop for performance)\u001b[39;00m\n\u001b[1;32m   1182\u001b[0m \u001b[38;5;28;01mif\u001b[39;00m \u001b[38;5;28mself\u001b[39m\u001b[38;5;241m.\u001b[39mdisable:\n\u001b[0;32m-> 1183\u001b[0m     \u001b[38;5;28;01mfor\u001b[39;00m obj \u001b[38;5;129;01min\u001b[39;00m iterable:\n\u001b[1;32m   1184\u001b[0m         \u001b[38;5;28;01myield\u001b[39;00m obj\n\u001b[1;32m   1185\u001b[0m     \u001b[38;5;28;01mreturn\u001b[39;00m\n",
      "File \u001b[0;32m~/anaconda3/lib/python3.10/site-packages/MDAnalysis/coordinates/base.py:714\u001b[0m, in \u001b[0;36mProtoReader.__next__\u001b[0;34m(self)\u001b[0m\n\u001b[1;32m    712\u001b[0m \u001b[38;5;28;01mdef\u001b[39;00m \u001b[38;5;21m__next__\u001b[39m(\u001b[38;5;28mself\u001b[39m) \u001b[38;5;241m-\u001b[39m\u001b[38;5;241m>\u001b[39m Timestep:\n\u001b[1;32m    713\u001b[0m     \u001b[38;5;124;03m\"\"\"Forward one step to next frame when using the `next` builtin.\"\"\"\u001b[39;00m\n\u001b[0;32m--> 714\u001b[0m     \u001b[38;5;28;01mreturn\u001b[39;00m \u001b[38;5;28;43mself\u001b[39;49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43mnext\u001b[49m\u001b[43m(\u001b[49m\u001b[43m)\u001b[49m\n",
      "File \u001b[0;32m~/anaconda3/lib/python3.10/site-packages/MDAnalysis/coordinates/base.py:700\u001b[0m, in \u001b[0;36mProtoReader.next\u001b[0;34m(self)\u001b[0m\n\u001b[1;32m    698\u001b[0m \u001b[38;5;124;03m\"\"\"Forward one step to next frame.\"\"\"\u001b[39;00m\n\u001b[1;32m    699\u001b[0m \u001b[38;5;28;01mtry\u001b[39;00m:\n\u001b[0;32m--> 700\u001b[0m     ts \u001b[38;5;241m=\u001b[39m \u001b[38;5;28;43mself\u001b[39;49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43m_read_next_timestep\u001b[49m\u001b[43m(\u001b[49m\u001b[43m)\u001b[49m\n\u001b[1;32m    701\u001b[0m \u001b[38;5;28;01mexcept\u001b[39;00m (\u001b[38;5;167;01mEOFError\u001b[39;00m, \u001b[38;5;167;01mIOError\u001b[39;00m):\n\u001b[1;32m    702\u001b[0m     \u001b[38;5;28mself\u001b[39m\u001b[38;5;241m.\u001b[39mrewind()\n",
      "File \u001b[0;32m~/anaconda3/lib/python3.10/site-packages/MDAnalysis/coordinates/TRJ.py:687\u001b[0m, in \u001b[0;36mNCDFReader._read_next_timestep\u001b[0;34m(self, ts)\u001b[0m\n\u001b[1;32m    685\u001b[0m     ts \u001b[38;5;241m=\u001b[39m \u001b[38;5;28mself\u001b[39m\u001b[38;5;241m.\u001b[39mts\n\u001b[1;32m    686\u001b[0m \u001b[38;5;28;01mtry\u001b[39;00m:\n\u001b[0;32m--> 687\u001b[0m     \u001b[38;5;28;01mreturn\u001b[39;00m \u001b[38;5;28;43mself\u001b[39;49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43m_read_frame\u001b[49m\u001b[43m(\u001b[49m\u001b[38;5;28;43mself\u001b[39;49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43m_current_frame\u001b[49m\u001b[43m \u001b[49m\u001b[38;5;241;43m+\u001b[39;49m\u001b[43m \u001b[49m\u001b[38;5;241;43m1\u001b[39;49m\u001b[43m)\u001b[49m\n\u001b[1;32m    688\u001b[0m \u001b[38;5;28;01mexcept\u001b[39;00m \u001b[38;5;167;01mIndexError\u001b[39;00m:\n\u001b[1;32m    689\u001b[0m     \u001b[38;5;28;01mraise\u001b[39;00m \u001b[38;5;167;01mIOError\u001b[39;00m \u001b[38;5;28;01mfrom\u001b[39;00m \u001b[38;5;28mNone\u001b[39m\n",
      "File \u001b[0;32m~/anaconda3/lib/python3.10/site-packages/MDAnalysis/coordinates/TRJ.py:675\u001b[0m, in \u001b[0;36mNCDFReader._read_frame\u001b[0;34m(self, frame)\u001b[0m\n\u001b[1;32m    672\u001b[0m         \u001b[38;5;28mself\u001b[39m\u001b[38;5;241m.\u001b[39mconvert_forces_from_native(ts\u001b[38;5;241m.\u001b[39m_forces, inplace\u001b[38;5;241m=\u001b[39m\u001b[38;5;28;01mTrue\u001b[39;00m)\n\u001b[1;32m    673\u001b[0m     \u001b[38;5;28;01mif\u001b[39;00m \u001b[38;5;28mself\u001b[39m\u001b[38;5;241m.\u001b[39mperiodic:\n\u001b[1;32m    674\u001b[0m         \u001b[38;5;66;03m# in-place ! (only lengths)\u001b[39;00m\n\u001b[0;32m--> 675\u001b[0m         \u001b[38;5;28mself\u001b[39m\u001b[38;5;241m.\u001b[39mconvert_pos_from_native(\u001b[43mts\u001b[49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43mdimensions\u001b[49m\u001b[43m[\u001b[49m\u001b[43m:\u001b[49m\u001b[38;5;241;43m3\u001b[39;49m\u001b[43m]\u001b[49m)\n\u001b[1;32m    676\u001b[0m ts\u001b[38;5;241m.\u001b[39mframe \u001b[38;5;241m=\u001b[39m frame  \u001b[38;5;66;03m# frame labels are 0-based\u001b[39;00m\n\u001b[1;32m    677\u001b[0m \u001b[38;5;28mself\u001b[39m\u001b[38;5;241m.\u001b[39m_current_frame \u001b[38;5;241m=\u001b[39m frame\n",
      "\u001b[0;31mTypeError\u001b[0m: 'NoneType' object is not subscriptable"
     ]
    }
   ],
   "source": [
    "from MDAnalysis.analysis.rms import RMSF\n",
    "\n",
    "# Select the protein C-alpha atoms\n",
    "calphas = u.select_atoms('protein and name CA')\n",
    "\n",
    "# Initialize the RMSF analysis object\n",
    "R = RMSF(calphas)\n",
    "R.run()\n",
    "\n",
    "# Extract the RMSF values\n",
    "rmsf_values = R.rmsf\n",
    "\n",
    "# Plot the RMSF values\n",
    "plt.figure(figsize=(5,5))\n",
    "plt.plot(calphas.resids, rmsf_values, 'o-')\n",
    "plt.xlabel('Residue number')\n",
    "plt.ylabel('RMSF (Å)')\n",
    "plt.title('RMSF Plot')\n",
    "plt.grid(True)\n",
    "plt.savefig('rmsf_plot.png')\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "id": "94958806",
   "metadata": {},
   "outputs": [
    {
     "ename": "NameError",
     "evalue": "name 'rmsf_values' is not defined",
     "output_type": "error",
     "traceback": [
      "\u001b[0;31m---------------------------------------------------------------------------\u001b[0m",
      "\u001b[0;31mNameError\u001b[0m                                 Traceback (most recent call last)",
      "Cell \u001b[0;32mIn[16], line 2\u001b[0m\n\u001b[1;32m      1\u001b[0m plt\u001b[38;5;241m.\u001b[39mfigure(figsize\u001b[38;5;241m=\u001b[39m(\u001b[38;5;241m10\u001b[39m,\u001b[38;5;241m5\u001b[39m))\n\u001b[0;32m----> 2\u001b[0m plt\u001b[38;5;241m.\u001b[39mplot(calphas\u001b[38;5;241m.\u001b[39mresids, \u001b[43mrmsf_values\u001b[49m, \u001b[38;5;124m'\u001b[39m\u001b[38;5;124mo-\u001b[39m\u001b[38;5;124m'\u001b[39m)\n\u001b[1;32m      3\u001b[0m plt\u001b[38;5;241m.\u001b[39mxlabel(\u001b[38;5;124m'\u001b[39m\u001b[38;5;124mResidue number\u001b[39m\u001b[38;5;124m'\u001b[39m)\n\u001b[1;32m      4\u001b[0m plt\u001b[38;5;241m.\u001b[39mylabel(\u001b[38;5;124m'\u001b[39m\u001b[38;5;124mRMSD (Å)\u001b[39m\u001b[38;5;124m'\u001b[39m)\n",
      "\u001b[0;31mNameError\u001b[0m: name 'rmsf_values' is not defined"
     ]
    },
    {
     "data": {
      "text/plain": [
       "<Figure size 1000x500 with 0 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "plt.figure(figsize=(10,5))\n",
    "plt.plot(calphas.resids, rmsf_values, 'o-')\n",
    "plt.xlabel('Residue number')\n",
    "plt.ylabel('RMSD (Å)')\n",
    "plt.title('RMSF Plot')\n",
    "plt.grid(True)\n",
    "plt.savefig('rmsf_plot.png')\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "c1829c2f",
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3 (ipykernel)",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.10.9"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}