FRETpredict/tutorials/Tutorial_FRETpredict_Hsp90.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "d05f6403",
   "metadata": {},
   "source": [
    "# FRETpredict tutorial (Hsp90)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "id": "8f0eaeed",
   "metadata": {},
   "outputs": [],
   "source": [
    "import MDAnalysis\n",
    "import numpy as np\n",
    "%matplotlib inline\n",
    "import matplotlib.pyplot as plt\n",
    "import matplotlib as mpl\n",
    "import os\n",
    "import pandas as pd\n",
    "import seaborn as sns\n",
    "import yaml\n",
    "from FRETpredict import FRETpredict\n",
    "import warnings\n",
    "\n",
    "warnings.filterwarnings('ignore')"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "89d3e848",
   "metadata": {},
   "source": [
    "### Quick biological background"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "c773bea2",
   "metadata": {},
   "source": [
    "Cartoon representation of the open structure of Hsp90 (apoprotein). Chains A and B are colored in red and blue, respectively. On each chain, the residues that will be used for FRET Efficiency calculations are represented as spheres."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "d8cbcb0b",
   "metadata": {},
   "source": [
    "![title](Hsp90_chains.png)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "bf5ecdea",
   "metadata": {},
   "source": [
    "Closed structure of Hsp90 (Holo protein, with ATP)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "9a87002e",
   "metadata": {},
   "source": [
    "![title](Hsp90_closed_chains.png)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "b71c73ad",
   "metadata": {},
   "source": [
    "# Single protein structure"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "ce1fd0a4",
   "metadata": {},
   "source": [
    "In this part of the tutorial, we will perform FRET Efficiency calculations using only a single residue pair to put the rotamer libraries on.<br>\n",
    "<br>First, let's have a look at the possible rotamer libraries we can use and how they're called:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "id": "1160a971",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "{'ATTO 390 C2R': {'author': 'D Montepietra, G Tesei, JM Martins, MBA Kunze, RB Best, K Lindorff-Larsen',\n",
       "  'citation': 'TBD',\n",
       "  'filename': 'T39_C2R_cutoff30',\n",
       "  'licence': 'GPLv3',\n",
       "  'mu': ['C1', 'C10 and resname T39'],\n",
       "  'negative': [],\n",
       "  'positive': [],\n",
       "  'r': ['C7 and resname T39']},\n",
       " 'ATTO 390 L1R': {'author': 'D Montepietra, G Tesei, JM Martins, MBA Kunze, RB Best, K Lindorff-Larsen',\n",
       "  'citation': 'TBD',\n",
       "  'filename': 'T39_L1R_cutoff30',\n",
       "  'licence': 'GPLv3',\n",
       "  'mu': ['C1', 'C10 and resname T39'],\n",
       "  'negative': [],\n",
       "  'positive': [],\n",
       "  'r': ['C7 and resname T39']},\n",
       " 'ATTO 425 C2R': {'author': 'D Montepietra, G Tesei, JM Martins, MBA Kunze, RB Best, K Lindorff-Larsen',\n",
       "  'citation': 'TBD',\n",
       "  'filename': 'T42_C2R_cutoff30',\n",
       "  'licence': 'GPLv3',\n",
       "  'mu': ['C1', 'C10 and resname T42'],\n",
       "  'negative': [],\n",
       "  'positive': [],\n",
       "  'r': ['C7 and resname T42']},\n",
       " 'ATTO 425 L1R': {'author': 'D Montepietra, G Tesei, JM Martins, MBA Kunze, RB Best, K Lindorff-Larsen',\n",
       "  'citation': 'TBD',\n",
       "  'filename': 'T42_L1R_cutoff30',\n",
       "  'licence': 'GPLv3',\n",
       "  'mu': ['C1', 'C10 and resname T42'],\n",
       "  'negative': [],\n",
       "  'positive': [],\n",
       "  'r': ['C7 and resname T42']},\n",
       " 'ATTO 465 C2R': {'author': 'D Montepietra, G Tesei, JM Martins, MBA Kunze, RB Best, K Lindorff-Larsen',\n",
       "  'citation': 'TBD',\n",
       "  'filename': 'T46_C2R_cutoff30',\n",
       "  'licence': 'GPLv3',\n",
       "  'mu': ['C1', 'C10 and resname T46'],\n",
       "  'negative': [],\n",
       "  'positive': ['N1 and resname T46',\n",
       "   'N2 and resname T46',\n",
       "   'N3 and resname T46'],\n",
       "  'r': ['N1 and resname T46']},\n",
       " 'ATTO 465 L1R': {'author': 'D Montepietra, G Tesei, JM Martins, MBA Kunze, RB Best, K Lindorff-Larsen',\n",
       "  'citation': 'TBD',\n",
       "  'filename': 'T46_L1R_cutoff30',\n",
       "  'licence': 'GPLv3',\n",
       "  'mu': ['C1', 'C10 and resname T46'],\n",
       "  'negative': [],\n",
       "  'positive': ['N1 and resname T46',\n",
       "   'N2 and resname T46',\n",
       "   'N3 and resname T46'],\n",
       "  'r': ['N1 and resname T46']},\n",
       " 'ATTO 488 C3R': {'author': 'D Montepietra, G Tesei, JM Martins, MBA Kunze, RB Best, K Lindorff-Larsen',\n",
       "  'citation': 'TBD',\n",
       "  'filename': 'T48_C3R_cutoff30',\n",
       "  'licence': 'GPLv3',\n",
       "  'mu': ['C1', 'C9 and resname T48'],\n",
       "  'negative': ['S1 and resname T48', 'S2 and resname T48'],\n",
       "  'positive': ['N2 and resname T48', 'N3 and resname T48'],\n",
       "  'r': ['C15 and resname T48']},\n",
       " 'ATTO 495 C2R': {'author': 'D Montepietra, G Tesei, JM Martins, MBA Kunze, RB Best, K Lindorff-Larsen',\n",
       "  'citation': 'TBD',\n",
       "  'filename': 'T49_C2R_cutoff30',\n",
       "  'licence': 'GPLv3',\n",
       "  'mu': ['C1', 'C10 and resname T49'],\n",
       "  'negative': [],\n",
       "  'positive': ['N1 and resname T49'],\n",
       "  'r': ['N1 and resname T49']},\n",
       " 'ATTO 495 L1R': {'author': 'D Montepietra, G Tesei, JM Martins, MBA Kunze, RB Best, K Lindorff-Larsen',\n",
       "  'citation': 'TBD',\n",
       "  'filename': 'T49_L1R_cutoff30',\n",
       "  'licence': 'GPLv3',\n",
       "  'mu': ['C1', 'C10 and resname T49'],\n",
       "  'negative': [],\n",
       "  'positive': ['N1 and resname T49'],\n",
       "  'r': ['N1 and resname T49']},\n",
       " 'ATTO 520 C2R': {'author': 'D Montepietra, G Tesei, JM Martins, MBA Kunze, RB Best, K Lindorff-Larsen',\n",
       "  'citation': 'TBD',\n",
       "  'filename': 'T52_C2R_cutoff30',\n",
       "  'licence': 'GPLv3',\n",
       "  'mu': ['C4', 'C11 and resname T52'],\n",
       "  'negative': [],\n",
       "  'positive': ['N1 and resname T52', 'N2 and resname T52'],\n",
       "  'r': ['C7 and resname T52']},\n",
       " 'ATTO 520 L1R': {'author': 'D Montepietra, G Tesei, JM Martins, MBA Kunze, RB Best, K Lindorff-Larsen',\n",
       "  'citation': 'TBD',\n",
       "  'filename': 'T52_L1R_cutoff30',\n",
       "  'licence': 'GPLv3',\n",
       "  'mu': ['C4', 'C11 and resname T52'],\n",
       "  'negative': [],\n",
       "  'positive': ['N1 and resname T52', 'N2 and resname T52'],\n",
       "  'r': ['C7 and resname T52']},\n",
       " 'ATTO 610 C2R': {'author': 'D Montepietra, G Tesei, JM Martins, MBA Kunze, RB Best, K Lindorff-Larsen',\n",
       "  'citation': 'TBD',\n",
       "  'filename': 'T61_C2R_cutoff30',\n",
       "  'licence': 'GPLv3',\n",
       "  'mu': ['C2', 'C17 and resname T61'],\n",
       "  'negative': [],\n",
       "  'positive': ['N1 and resname T61', 'N2 and resname T61'],\n",
       "  'r': ['C8 and resname T61']},\n",
       " 'ATTO 610 L1R': {'author': 'D Montepietra, G Tesei, JM Martins, MBA Kunze, RB Best, K Lindorff-Larsen',\n",
       "  'citation': 'TBD',\n",
       "  'filename': 'T61_L1R_cutoff30',\n",
       "  'licence': 'GPLv3',\n",
       "  'mu': ['C2', 'C17 and resname T61'],\n",
       "  'negative': [],\n",
       "  'positive': ['N1 and resname T61', 'N2 and resname T61'],\n",
       "  'r': ['C8 and resname T61']},\n",
       " 'ATTO Thio12 C3R': {'author': 'D Montepietra, G Tesei, JM Martins, MBA Kunze, RB Best, K Lindorff-Larsen',\n",
       "  'citation': 'TBD',\n",
       "  'filename': 'Tth_C3R_cutoff30',\n",
       "  'licence': 'GPLv3',\n",
       "  'mu': ['C5', 'C11 and resname Tth'],\n",
       "  'negative': [],\n",
       "  'positive': ['N2 and resname Tth', 'N3 and resname Tth'],\n",
       "  'r': ['C9 and resname Tth']},\n",
       " 'AlexaFluor 350 C1R': {'author': 'D Montepietra, G Tesei, JM Martins, MBA Kunze, RB Best, K Lindorff-Larsen',\n",
       "  'citation': 'TBD',\n",
       "  'filename': 'A35_C1R_cutoff30',\n",
       "  'licence': 'GPLv3',\n",
       "  'mu': ['C2', 'C9 and resname A35'],\n",
       "  'negative': ['S1 and resname A35'],\n",
       "  'positive': ['N1 and resname A35'],\n",
       "  'r': ['C5 and resname A35']},\n",
       " 'AlexaFluor 350 L1R': {'author': 'D Montepietra, G Tesei, JM Martins, MBA Kunze, RB Best, K Lindorff-Larsen',\n",
       "  'citation': 'TBD',\n",
       "  'filename': 'A35_L1R_cutoff30',\n",
       "  'licence': 'GPLv3',\n",
       "  'mu': ['C2', 'C9 and resname A35'],\n",
       "  'negative': ['S1 and resname A35'],\n",
       "  'positive': ['N1 and resname A35'],\n",
       "  'r': ['C5 and resname A35']},\n",
       " 'AlexaFluor 488 B1R': {'author': 'D Montepietra, G Tesei, JM Martins, MBA Kunze, RB Best, K Lindorff-Larsen',\n",
       "  'citation': 'TBD',\n",
       "  'filename': 'A48_B1R_cutoff30',\n",
       "  'licence': 'GPLv3',\n",
       "  'mu': ['C2 and resname A48', 'C13 and resname A48'],\n",
       "  'negative': ['S1 and resname A48',\n",
       "   'S2 and resname A48',\n",
       "   'C20 and resname A48'],\n",
       "  'positive': ['N1 and resname A48', 'N2 and resname A48'],\n",
       "  'r': ['C7 and resname A48']},\n",
       " 'AlexaFluor 488 C1R': {'author': 'D Montepietra, G Tesei, JM Martins, MBA Kunze, RB Best, K Lindorff-Larsen',\n",
       "  'citation': 'TBD',\n",
       "  'filename': 'A48_C1R_cutoff30',\n",
       "  'licence': 'GPLv3',\n",
       "  'mu': ['C2', 'C13 and resname A48'],\n",
       "  'negative': ['S1 and resname A48',\n",
       "   'S2 and resname A48',\n",
       "   'C20 and resname A48'],\n",
       "  'positive': ['N1 and resname A48', 'N2 and resname A48'],\n",
       "  'r': ['C7 and resname A48']},\n",
       " 'AlexaFluor 488 L1R': {'author': 'D Montepietra, G Tesei, JM Martins, MBA Kunze, RB Best, K Lindorff-Larsen',\n",
       "  'citation': 'TBD',\n",
       "  'filename': 'A48_L1R_cutoff30',\n",
       "  'licence': 'GPLv3',\n",
       "  'mu': ['C2', 'C13 and resname A48'],\n",
       "  'negative': ['S1 and resname A48',\n",
       "   'S2 and resname A48',\n",
       "   'C20 and resname A48'],\n",
       "  'positive': ['N1 and resname A48', 'N2 and resname A48'],\n",
       "  'r': ['C7 and resname A48']},\n",
       " 'AlexaFluor 532 C1R': {'author': 'D Montepietra, G Tesei, JM Martins, MBA Kunze, RB Best, K Lindorff-Larsen',\n",
       "  'citation': 'TBD',\n",
       "  'filename': 'A53_C1R_cutoff30',\n",
       "  'licence': 'GPLv3',\n",
       "  'mu': ['C2', 'C13 and resname A53'],\n",
       "  'negative': ['S1 and resname A53', 'S2 and resname A53'],\n",
       "  'positive': ['N1 and resname A53', 'N2 and resname A53'],\n",
       "  'r': ['C7 and resname A53']},\n",
       " 'AlexaFluor 532 L1R': {'author': 'D Montepietra, G Tesei, JM Martins, MBA Kunze, RB Best, K Lindorff-Larsen',\n",
       "  'citation': 'TBD',\n",
       "  'filename': 'A53_L1R_cutoff30',\n",
       "  'licence': 'GPLv3',\n",
       "  'mu': ['C2', 'C13 and resname A53'],\n",
       "  'negative': ['S1 and resname A53', 'S2 and resname A53'],\n",
       "  'positive': ['N1 and resname A53', 'N2 and resname A53'],\n",
       "  'r': ['C7 and resname A53']},\n",
       " 'AlexaFluor 568 C1R': {'author': 'D Montepietra, G Tesei, JM Martins, MBA Kunze, RB Best, K Lindorff-Larsen',\n",
       "  'citation': 'TBD',\n",
       "  'filename': 'A56_C1R_cutoff30',\n",
       "  'licence': 'GPLv3',\n",
       "  'mu': ['C5', 'C16 and resname A56'],\n",
       "  'negative': ['S1 and resname A56',\n",
       "   'S2 and resname A56',\n",
       "   'C32 and resname A56'],\n",
       "  'positive': ['N1 and resname A56', 'N2 and resname A56'],\n",
       "  'r': ['C10 and resname A56']},\n",
       " 'AlexaFluor 568 L1R': {'author': 'D Montepietra, G Tesei, JM Martins, MBA Kunze, RB Best, K Lindorff-Larsen',\n",
       "  'citation': 'TBD',\n",
       "  'filename': 'A56_L1R_cutoff30',\n",
       "  'licence': 'GPLv3',\n",
       "  'mu': ['C5', 'C16 and resname A56'],\n",
       "  'negative': ['S1 and resname A56',\n",
       "   'S2 and resname A56',\n",
       "   'C32 and resname A56'],\n",
       "  'positive': ['N1 and resname A56', 'N2 and resname A56'],\n",
       "  'r': ['C10 and resname A56']},\n",
       " 'AlexaFluor 594 C1R': {'author': 'D Montepietra, G Tesei, JM Martins, MBA Kunze, RB Best, K Lindorff-Larsen',\n",
       "  'citation': 'TBD',\n",
       "  'filename': 'A59_C1R_cutoff30',\n",
       "  'licence': 'GPLv3',\n",
       "  'mu': ['C5', 'C16 and resname A59'],\n",
       "  'negative': ['S1 and resname A59',\n",
       "   'S2 and resname A59',\n",
       "   'C26 and resname A56'],\n",
       "  'positive': ['N1 and resname A59', 'N2 and resname A59'],\n",
       "  'r': ['C10 and resname A59']},\n",
       " 'AlexaFluor 594 L1R': {'author': 'D Montepietra, G Tesei, JM Martins, MBA Kunze, RB Best, K Lindorff-Larsen',\n",
       "  'citation': 'TBD',\n",
       "  'filename': 'A59_L1R_cutoff30',\n",
       "  'licence': 'GPLv3',\n",
       "  'mu': ['C5', 'C16 and resname A59'],\n",
       "  'negative': ['S1 and resname A59',\n",
       "   'S2 and resname A59',\n",
       "   'C26 and resname A56'],\n",
       "  'positive': ['N1 and resname A59', 'N2 and resname A59'],\n",
       "  'r': ['C10 and resname A59']},\n",
       " 'AlexaFluor 647 C2R': {'author': 'D Montepietra, G Tesei, JM Martins, MBA Kunze, RB Best, K Lindorff-Larsen',\n",
       "  'citation': 'TBD',\n",
       "  'filename': 'A64_C2R_cutoff30',\n",
       "  'licence': 'GPLv3',\n",
       "  'mu': ['C4', 'C15 and resname A64'],\n",
       "  'negative': ['S1 and resname A64',\n",
       "   'S2 and resname A64',\n",
       "   'S3 and resname A64',\n",
       "   'S4 and resname A64'],\n",
       "  'positive': ['N1 and resname A64', 'N2 and resname A64'],\n",
       "  'r': ['C29 and resname A64']},\n",
       " 'AlexaFluor 647 L1R': {'author': 'D Montepietra, G Tesei, JM Martins, MBA Kunze, RB Best, K Lindorff-Larsen',\n",
       "  'citation': 'TBD',\n",
       "  'filename': 'A64_L1R_cutoff30',\n",
       "  'licence': 'GPLv3',\n",
       "  'mu': ['C4', 'C15 and resname A64'],\n",
       "  'negative': ['S1 and resname A64',\n",
       "   'S2 and resname A64',\n",
       "   'S3 and resname A64',\n",
       "   'S4 and resname A64'],\n",
       "  'positive': ['N1 and resname A64', 'N2 and resname A64'],\n",
       "  'r': ['C29 and resname A64']},\n",
       " 'Lumiprobe Cy3 C2R': {'author': 'D Montepietra, G Tesei, JM Martins, MBA Kunze, RB Best, K Lindorff-Larsen',\n",
       "  'citation': 'TBD',\n",
       "  'filename': 'C3N_C2R_cutoff30',\n",
       "  'licence': 'GPLv3',\n",
       "  'mu': ['C5', 'C16 and resname C3N'],\n",
       "  'negative': [],\n",
       "  'positive': ['N1 and resname C3N', 'N2 and resname C3N'],\n",
       "  'r': ['C10 and resname C3N']},\n",
       " 'Lumiprobe Cy5.5 C2R': {'author': 'D Montepietra, G Tesei, JM Martins, MBA Kunze, RB Best, K Lindorff-Larsen',\n",
       "  'citation': 'TBD',\n",
       "  'filename': 'C55_C2R_cutoff30',\n",
       "  'licence': 'GPLv3',\n",
       "  'mu': ['C5', 'C16 and resname C55'],\n",
       "  'negative': [],\n",
       "  'positive': ['N1 and resname C55', 'N2 and resname C55'],\n",
       "  'r': ['C39 and resname C55']},\n",
       " 'Lumiprobe Cy5.5 L1R': {'author': 'D Montepietra, G Tesei, JM Martins, MBA Kunze, RB Best, K Lindorff-Larsen',\n",
       "  'citation': 'TBD',\n",
       "  'filename': 'C55_L1R_cutoff30',\n",
       "  'licence': 'GPLv3',\n",
       "  'mu': ['C5', 'C16 and resname C55'],\n",
       "  'negative': [],\n",
       "  'positive': ['N1 and resname C55', 'N2 and resname C55'],\n",
       "  'r': ['C39 and resname C55']},\n",
       " 'Lumiprobe Cy7.5 L1R': {'author': 'D Montepietra, G Tesei, JM Martins, MBA Kunze, RB Best, K Lindorff-Larsen',\n",
       "  'citation': 'TBD',\n",
       "  'filename': 'C75_L1R_cutoff30',\n",
       "  'licence': 'GPLv3',\n",
       "  'mu': ['C5', 'C16 and resname C75'],\n",
       "  'negative': [],\n",
       "  'positive': ['N1 and resname C75', 'N2 and resname C75'],\n",
       "  'r': ['C43 and resname C75']}}"
      ]
     },
     "execution_count": 2,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "with open('../lib/libraries.yml') as f:\n",
    "    libraries = yaml.load(f)\n",
    "libraries"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "6022f912",
   "metadata": {},
   "source": [
    "Every Rotamer Library name is composed of three parts: the producer (AlexaFluor, ATTO, Lumiprobe), the peak wavelength (e.g. 488, 550, 647), and the linker that connects the dye to the protein (C1R, C2R, C3R, L1R, L2R, B1R).<br>\n",
    "<br>To learn more about the rotamer libraries, see `Tutorial_generate_new_rotamer libraries`."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "ff2d7ef0",
   "metadata": {},
   "source": [
    "Now, we will select the parameters for the FRET Efficiency calculation: \n",
    "- The `residues` to put the rotamer libraries on, their protein `chain`, and the corresponding experimental data file.\n",
    "For this tutorial we're going to use a residue pair which FRET Efficiency was obtained on an open (APO) structure of Hsp90. \n",
    "- The rotamer libraries we will use: AlexaFluor dyes 594 and 568, with C1R linkers.\n",
    "`donor` and `acceptor` are used for R0 calculations, while `libname_1` and `libname_2` for the FRET Efficiency calculation.\n",
    "- The clustering `cutoff` for the rotamer libraries.\n",
    "- The Universe object for the `protein` structure.\n",
    "- `electrostatic` calculations will be enabled, and the `temperature` will be set at 293K."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "id": "6fc11f8b",
   "metadata": {},
   "outputs": [],
   "source": [
    "res1 = 452\n",
    "chain_res1 = 'A'\n",
    "dye_1 = 'AlexaFluor 594'\n",
    "linker_1 = 'C1R'\n",
    "\n",
    "res2 = 637\n",
    "chain_res2 = 'B'\n",
    "dye_2 = 'AlexaFluor 568'\n",
    "linker_2 = 'C1R'\n",
    "\n",
    "# Clustering cutoff used for the creation of the rotamer library. It can be 10, 20, or 30.\n",
    "cutoff = 10\n",
    "\n",
    "# Experimental data to compare our predictions.\n",
    "exp_data_file = '452_532_637_647N_APO.txt'\n",
    "\n",
    "# Create MDAnalysis.Universe object for the protein (open Hsp90)\n",
    "u = MDAnalysis.Universe('../test_systems/Hsp90/openHsp90.pdb')"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "94cb7928",
   "metadata": {},
   "source": [
    "Let's create an instance of the FRETpredict class"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "id": "b2b5f89d",
   "metadata": {},
   "outputs": [],
   "source": [
    "FRET = FRETpredict(protein=u, residues=[res1, res2], temperature=293, \n",
    "                   chains=[chain_res1, chain_res2], \n",
    "                   donor=dye_1, acceptor=dye_2, electrostatic=True,\n",
    "                   r0lib='../lib/R0/',\n",
    "                   libname_1=f'{dye_1} {linker_1} cutoff{cutoff}',\n",
    "                   libname_2=f'{dye_2} {linker_2} cutoff{cutoff}', \n",
    "                   output_prefix=f'test/E{cutoff}')"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "744de18d",
   "metadata": {},
   "source": [
    "Run FRET efficiency calculations."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "id": "7f7c63ed",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\n",
      "Frame 1/1\n",
      "\n",
      "Done.\n"
     ]
    }
   ],
   "source": [
    "FRET.run()"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "2774d116",
   "metadata": {},
   "source": [
    "Retrieve experimental data from file for the comparison."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "id": "f39a8a05",
   "metadata": {},
   "outputs": [],
   "source": [
    "filename = '../test_systems/Hsp90/MDA/histograms/{}'.format(exp_data_file)\n",
    "\n",
    "if os.path.isfile(filename):\n",
    "    \n",
    "        df = pd.read_csv(filename, skiprows=6, nrows=40, header=None, sep='\\s+')\n",
    "        Ex = np.average(df[0], weights=df[1])"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "9067bd39",
   "metadata": {},
   "source": [
    "Create DataFrame of the data (experimental and predicted) and show the results."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "id": "3845e162",
   "metadata": {},
   "outputs": [],
   "source": [
    "results = []\n",
    "\n",
    "data = pd.read_pickle(r'test/E{:d}-data-{:d}-{:d}.pkl'.format(cutoff, res1, res2))\n",
    "results.append({'res': '{:d}-{:d}'.format(res1, res2),\n",
    "                'k2': float(data.loc['k2']),\n",
    "                'Ex': Ex,\n",
    "                'Es': float(data.loc['Estatic']),\n",
    "                'Ed': float(data.loc['Edynamic1']),\n",
    "                'Ed2': float(data.loc['Edynamic2']),\n",
    "                'diffS': np.abs(Ex-data.loc['Estatic']),\n",
    "                'diffD': np.abs(Ex-data.loc['Edynamic1']), \n",
    "                'diffD2': np.abs(Ex-data.loc['Edynamic2'])})\n",
    "        \n",
    "results_df = pd.DataFrame(results).set_index('res')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "id": "34b018fc",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>k2</th>\n",
       "      <th>Ex</th>\n",
       "      <th>Es</th>\n",
       "      <th>Ed</th>\n",
       "      <th>Ed2</th>\n",
       "      <th>diffS</th>\n",
       "      <th>diffD</th>\n",
       "      <th>diffD2</th>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>res</th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>452-637</th>\n",
       "      <td>0.584892</td>\n",
       "      <td>0.592377</td>\n",
       "      <td>0.016119</td>\n",
       "      <td>0.017959</td>\n",
       "      <td>0.987855</td>\n",
       "      <td>0.576258</td>\n",
       "      <td>0.574418</td>\n",
       "      <td>0.395478</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "               k2        Ex        Es        Ed       Ed2     diffS     diffD  \\\n",
       "res                                                                             \n",
       "452-637  0.584892  0.592377  0.016119  0.017959  0.987855  0.576258  0.574418   \n",
       "\n",
       "           diffD2  \n",
       "res                \n",
       "452-637  0.395478  "
      ]
     },
     "execution_count": 8,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Show results\n",
    "results_df"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "0b5b2ee6",
   "metadata": {},
   "source": [
    ">For single frame or structures predictions, the best averaging method is \"Static\" (`Es`)."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "00f31a17",
   "metadata": {},
   "source": [
    "## R0 calculation"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "c5e420a0",
   "metadata": {},
   "source": [
    "It is possible to select the R0 calculation method for the FRET Efficiency calculation.<br>\n",
    "<br>By default, R0 is calculated for the provided dyes pair (`donor` and `acceptor`) for every rotamer conformation, taking the relative orientation into account.<br>\n",
    "<br>However, it is possible to select fixed R0 value. The only change is in the parameters passed to FRETpredict.<br>\n",
    "<br>Let's see how it's done."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "id": "5d5e9580",
   "metadata": {},
   "outputs": [],
   "source": [
    "FRET_fixedR0 = FRETpredict(protein=u, residues=[res1, res2], temperature=293, \n",
    "                           chains=[chain_res1, chain_res2], \n",
    "                           fixed_R0=True, r0=5.5, electrostatic=True,\n",
    "                           libname_1=f'{dye_1} {linker_1} cutoff{cutoff}',\n",
    "                           libname_2=f'{dye_2} {linker_2} cutoff{cutoff}', \n",
    "                           output_prefix=f'test/E{cutoff}_fixedR0')"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "79e04d0e",
   "metadata": {},
   "source": [
    "We set `fixed_R0=True` and passed the `r0` parameter the R0 value we want our dye pair to have.<br>In __[`lib/R0/R0_pairs.csv`](https://github.com/Monte95/FRETpredict/blob/ce327c60cf5f86a33251b2e0c60cf1f668bc46e2/FRETpredict/lib/R0/R0_pairs.csv)__ are reported the R0 values for many dye pairs.<br>\n",
    "<br>Otherwise, if the R0 value is not present in the file, online services as __[FPBase](https://www.fpbase.org/fret/)__ can be used to obtain it.\n",
    "<br>Now let's run FRET Efficiency calculations."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "id": "847dcdc1",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\n",
      "Frame 1/1\n",
      "\n",
      "Done.\n"
     ]
    }
   ],
   "source": [
    "FRET_fixedR0.run()"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "46909d67",
   "metadata": {},
   "source": [
    "Retrieve experimental data from file for the comparison."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "id": "e2f9cb96",
   "metadata": {},
   "outputs": [],
   "source": [
    "filename = '../test_systems/Hsp90/MDA/histograms/{}'.format(exp_data_file)\n",
    "\n",
    "if os.path.isfile(filename):\n",
    "    \n",
    "        df = pd.read_csv(filename, skiprows=6, nrows=40, header=None, sep='\\s+')\n",
    "        Ex = np.average(df[0], weights=df[1])"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "7647f5c2",
   "metadata": {},
   "source": [
    "Create DataFrame of the data (experimental and predicted) and show the results."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "id": "e4f7a63b",
   "metadata": {},
   "outputs": [],
   "source": [
    "results_fixedR0 = []\n",
    "\n",
    "data = pd.read_pickle(r'test/E{:d}_fixedR0-data-{:d}-{:d}.pkl'.format(cutoff, res1, res2))\n",
    "results_fixedR0.append({'res': '{:d}-{:d}'.format(res1, res2),\n",
    "                'k2': float(data.loc['k2']),\n",
    "                'Ex': Ex,\n",
    "                'Es': float(data.loc['Estatic']),\n",
    "                'Ed': float(data.loc['Edynamic1']),\n",
    "                'Ed2': float(data.loc['Edynamic2']),\n",
    "                'diffS': np.abs(Ex-data.loc['Estatic']),\n",
    "                'diffD': np.abs(Ex-data.loc['Edynamic1']), \n",
    "                'diffD2': np.abs(Ex-data.loc['Edynamic2'])})\n",
    "        \n",
    "results_fixedR0_df = pd.DataFrame(results_fixedR0).set_index('res')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "id": "ccb3b44f",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>k2</th>\n",
       "      <th>Ex</th>\n",
       "      <th>Es</th>\n",
       "      <th>Ed</th>\n",
       "      <th>Ed2</th>\n",
       "      <th>diffS</th>\n",
       "      <th>diffD</th>\n",
       "      <th>diffD2</th>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>res</th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>452-637</th>\n",
       "      <td>0.584892</td>\n",
       "      <td>0.592377</td>\n",
       "      <td>0.025657</td>\n",
       "      <td>0.028789</td>\n",
       "      <td>0.995613</td>\n",
       "      <td>0.56672</td>\n",
       "      <td>0.563588</td>\n",
       "      <td>0.403236</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "               k2        Ex        Es        Ed       Ed2    diffS     diffD  \\\n",
       "res                                                                            \n",
       "452-637  0.584892  0.592377  0.025657  0.028789  0.995613  0.56672  0.563588   \n",
       "\n",
       "           diffD2  \n",
       "res                \n",
       "452-637  0.403236  "
      ]
     },
     "execution_count": 13,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "results_fixedR0_df"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "d08d7fa6",
   "metadata": {},
   "source": [
    "# Protein trajectory (open conformations)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "6cffddf4",
   "metadata": {},
   "source": [
    "Now, we will see how to predict FRET Efficiency values for protein trajectories. The procedure is exactly the same, as the only difference is the Universe objects now containing multiple structures.<br>\n",
    "<br>We will use the same residues and chromophores as in the previous example (here reported for clarity)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "id": "3d923b27",
   "metadata": {},
   "outputs": [],
   "source": [
    "res1 = 452\n",
    "chain_res1 = 'A'\n",
    "dye_1 = 'AlexaFluor 594'\n",
    "linker_1 = 'C1R'\n",
    "\n",
    "res2 = 637\n",
    "chain_res2 = 'B'\n",
    "dye_2 = 'AlexaFluor 568'\n",
    "linker_2 = 'C1R'\n",
    "\n",
    "# Clustering cutoff used for the creation of the rotamer library. It can be 10, 20, or 30.\n",
    "cutoff = 10\n",
    "\n",
    "# Experimental data to compare our predictions.\n",
    "exp_data_file = '452_532_637_647N_APO.txt'\n",
    "\n",
    "# Create universe for protein trajectory (3 different open conformations of Hsp90).\n",
    "u_traj = MDAnalysis.Universe('../test_systems/Hsp90/conf.pdb', \n",
    "                             '../test_systems/Hsp90/Hsp90_open_all.xtc')"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "52d28749",
   "metadata": {},
   "source": [
    "Create instance of the FRETpredict class"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "id": "ae49cd5e",
   "metadata": {},
   "outputs": [],
   "source": [
    "FRET_traj = FRETpredict(protein=u_traj, residues=[res1, res2], temperature=293, \n",
    "                        chains=[chain_res1, chain_res2], \n",
    "                        donor=dye_1, acceptor=dye_2, electrostatic=True,\n",
    "                        r0lib='../lib/R0/',\n",
    "                        libname_1=f'{dye_1} {linker_1} cutoff{cutoff}',\n",
    "                        libname_2=f'{dye_2} {linker_2} cutoff{cutoff}', \n",
    "                        output_prefix=f'test/E{cutoff}_traj')"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "a3c2f350",
   "metadata": {},
   "source": [
    "Run FRET Efficiency calculations"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "id": "3c57702b",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\n",
      "Frame 1/3\n",
      "\n",
      "Frame 2/3\n",
      "\n",
      "Frame 3/3\n",
      "\n",
      "Done.\n"
     ]
    }
   ],
   "source": [
    "FRET_traj.run()"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "87920b18",
   "metadata": {},
   "source": [
    "Retrieve experimental data from file for the comparison"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "id": "13ebf35f",
   "metadata": {},
   "outputs": [],
   "source": [
    "filename = '../test_systems/Hsp90/MDA/histograms/{}'.format(exp_data_file)\n",
    "\n",
    "if os.path.isfile(filename):\n",
    "\n",
    "        df = pd.read_csv(filename, skiprows=6, nrows=40, header=None, sep='\\s+')\n",
    "        Ex = np.average(df[0], weights=df[1])"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "be2c4db8",
   "metadata": {},
   "source": [
    "Create DataFrame for the data "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "id": "298263d9",
   "metadata": {},
   "outputs": [],
   "source": [
    "results_traj = []\n",
    "\n",
    "data_traj = pd.read_pickle(r'test/E{:d}_traj-data-{:d}-{:d}.pkl'.format(cutoff, \n",
    "                                                                                  res1, res2))\n",
    "results_traj.append({'res': '{:d}-{:d}'.format(res1, res2),\n",
    "                     'k2': data_traj.loc['k2', 'Average'],\n",
    "                     'Ex': Ex,\n",
    "                     'Es': data_traj.loc['Estatic', 'Average'],\n",
    "                     'Ed': data_traj.loc['Edynamic1', 'Average'],\n",
    "                     'Ed2': data_traj.loc['Edynamic2', 'Average'],\n",
    "                     'diffS': np.abs(Ex-data_traj.loc['Estatic', 'Average']),\n",
    "                     'diffD': np.abs(Ex-data_traj.loc['Edynamic1', 'Average']), \n",
    "                     'diffD2': np.abs(Ex-data_traj.loc['Edynamic2', 'Average'])})\n",
    "        \n",
    "results_traj_df = pd.DataFrame(results_traj).set_index('res')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "id": "d69797cc",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>k2</th>\n",
       "      <th>Ex</th>\n",
       "      <th>Es</th>\n",
       "      <th>Ed</th>\n",
       "      <th>Ed2</th>\n",
       "      <th>diffS</th>\n",
       "      <th>diffD</th>\n",
       "      <th>diffD2</th>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>res</th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>452-637</th>\n",
       "      <td>0.603643</td>\n",
       "      <td>0.592377</td>\n",
       "      <td>0.023215</td>\n",
       "      <td>0.025707</td>\n",
       "      <td>0.950529</td>\n",
       "      <td>0.569163</td>\n",
       "      <td>0.56667</td>\n",
       "      <td>0.358151</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "               k2        Ex        Es        Ed       Ed2     diffS    diffD  \\\n",
       "res                                                                            \n",
       "452-637  0.603643  0.592377  0.023215  0.025707  0.950529  0.569163  0.56667   \n",
       "\n",
       "           diffD2  \n",
       "res                \n",
       "452-637  0.358151  "
      ]
     },
     "execution_count": 19,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Show results\n",
    "results_traj_df"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "9984397d",
   "metadata": {},
   "source": [
    ">For trajectory data or multiple protein structures, the best averaging method is \"Dynamic2\" (`Ed2`)"
   ]
  }
 ],
 "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.9.7"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}