error_estimation_T_boosting.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Error estimate for T-boosted Langevin simulations"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "* Import external libraries:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "%matplotlib inline\n",
    "import math\n",
    "import matplotlib.pyplot as plt\n",
    "from matplotlib import animation, rc\n",
    "from IPython.display import HTML\n",
    "from decimal import *\n",
    "\n",
    "import copy\n",
    "import numpy as np\n",
    "import scipy as scy"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Example trypsin**"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "number_points=13\n",
    "temperature=np.zeros(number_points)\n",
    "temperature_square=np.zeros(number_points)\n",
    "number_events = np.zeros(number_points)\n",
    "number_events_divided_by_T = np.zeros(number_points)\n",
    "number_events_divided_by_T_square = np.zeros(number_points)\n",
    "\n",
    "# Insert number of observed transitions\n",
    "\n",
    "number_events[0]=92\n",
    "number_events[1]=148\n",
    "number_events[2]=306\n",
    "number_events[3]=712\n",
    "number_events[4]=2153\n",
    "number_events[5]=5141\n",
    "number_events[6]=10534\n",
    "number_events[7]=18818\n",
    "number_events[8]=30163\n",
    "number_events[9]=44864\n",
    "number_events[10]=62538\n",
    "number_events[11]=82528\n",
    "number_events[12]=105127\n",
    "\n",
    "# Insert used temperatures\n",
    "\n",
    "temperature[0]=380\n",
    "temperature[1]=400\n",
    "temperature[2]=420\n",
    "temperature[3]=450\n",
    "temperature[4]=500\n",
    "temperature[5]=550\n",
    "temperature[6]=600\n",
    "temperature[7]=650\n",
    "temperature[8]=700\n",
    "temperature[9]=750\n",
    "temperature[10]=800\n",
    "temperature[11]=850\n",
    "temperature[12]=900\n",
    "\n",
    "for i in range(number_points):\n",
    "    temperature_square[i]=temperature[i]*temperature[i]\n",
    "    number_events_divided_by_T[i]=number_events[i]/temperature[i]\n",
    "    number_events_divided_by_T_square[i]=number_events[i]/temperature_square[i]\n",
    "\n",
    "# Calculate sums needed for estimation of error\n",
    "    \n",
    "sum1=np.sum(number_events)\n",
    "sum2=np.sum(number_events_divided_by_T)\n",
    "sum3=np.sum(number_events_divided_by_T_square)\n",
    "Delta=(sum1*sum3)-(sum2*sum2)\n",
    "# Calculate error for linear regression\n",
    "\n",
    "sigma_a_square=np.sum(number_events)/Delta\n",
    "sigma_b_square=np.sum(number_events_divided_by_T_square)/Delta\n",
    "sigma_a=np.sqrt(sigma_a_square)\n",
    "sigma_b=np.sqrt(sigma_b_square)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "* error estimate"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "3.336828187030928\n",
      "3.1169679662113894\n"
     ]
    }
   ],
   "source": [
    "# Use error propagation to derive error at 300 K.\n",
    "\n",
    "error=np.sqrt((sigma_a_square/(300*300))+sigma_b_square)\n",
    "error_just_a=np.sqrt((sigma_a_square/(300*300)))\n",
    "\n",
    "# Complete error in %\n",
    "print(error*100)\n",
    "\n",
    "# Error only from slope in %, we see that it dominates the error\n",
    "print(error_just_a*100)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "This means that we observe an error of 3.3 % for Trypsin"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Example Hsp90**"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {
    "scrolled": true
   },
   "outputs": [],
   "source": [
    "number_points=14\n",
    "temperature=np.zeros(number_points)\n",
    "temperature_square=np.zeros(number_points)\n",
    "number_events = np.zeros(number_points)\n",
    "number_events_divided_by_T = np.zeros(number_points)\n",
    "number_events_divided_by_T_square = np.zeros(number_points)\n",
    "number_events[0]=77\n",
    "number_events[1]=204\n",
    "number_events[2]=422\n",
    "number_events[3]=800\n",
    "number_events[4]=1376\n",
    "number_events[5]=2147\n",
    "number_events[6]=3265\n",
    "number_events[7]=4796\n",
    "number_events[8]=6861\n",
    "number_events[9]=9229\n",
    "number_events[10]=12233\n",
    "number_events[11]=15753\n",
    "number_events[12]=19832\n",
    "number_events[13]=24496\n",
    "\n",
    "temperature[0]=700\n",
    "temperature[1]=750\n",
    "temperature[2]=800\n",
    "temperature[3]=850\n",
    "temperature[4]=900\n",
    "temperature[5]=950\n",
    "temperature[6]=1000\n",
    "temperature[7]=1050\n",
    "temperature[8]=1100\n",
    "temperature[9]=1150\n",
    "temperature[10]=1200\n",
    "temperature[11]=1250\n",
    "temperature[12]=1300\n",
    "temperature[13]=1350\n",
    "\n",
    "for i in range(number_points):\n",
    "    temperature_square[i]=temperature[i]*temperature[i]\n",
    "    number_events_divided_by_T[i]=number_events[i]/temperature[i]\n",
    "    number_events_divided_by_T_square[i]=number_events[i]/temperature_square[i]\n",
    "    \n",
    "sum1=np.sum(number_events)\n",
    "sum2=np.sum(number_events_divided_by_T)\n",
    "sum3=np.sum(number_events_divided_by_T_square)\n",
    "Delta=(sum1*sum3)-(sum2*sum2)\n",
    "sigma_a_square=np.sum(number_events)/Delta\n",
    "sigma_b_square=np.sum(number_events_divided_by_T_square)/Delta\n",
    "sigma_a=np.sqrt(sigma_a_square)\n",
    "sigma_b=np.sqrt(sigma_b_square)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "* error estimate"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "11.030225384847082\n",
      "10.701332719990338\n"
     ]
    }
   ],
   "source": [
    "error=np.sqrt((sigma_a_square/(300*300))+sigma_b_square)\n",
    "error_just_a=np.sqrt((sigma_a_square/(300*300)))\n",
    "print(error*100)\n",
    "print(error_just_a*100)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "This means that we have 11% uncertaintiy for Hsp90."
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "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.6.5"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}