LOPO_predict_missingFeatureMuts.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false,
    "deletable": true,
    "editable": true
   },
   "outputs": [],
   "source": [
    "import graphlab as gl\n",
    "import pandas as pd\n",
    "import numpy as np\n",
    "import scipy.stats.stats as stats\n",
    "#%load_ext rpy2.ipython\n",
    "import itertools\n",
    "from collections import Counter\n",
    "from sklearn.ensemble import GradientBoostingClassifier  #GBM algorithm\n",
    "from sklearn import cross_validation, metrics, preprocessing, datasets   #Additional scklearn functions\n",
    "from sklearn.grid_search import GridSearchCV   #Performing grid search\n",
    "from sklearn.utils import shuffle\n",
    "from sklearn.metrics import mean_squared_error\n",
    "from sklearn.linear_model import LogisticRegression\n",
    "from sklearn.feature_selection import SelectKBest\n",
    "from sklearn.feature_extraction import DictVectorizer\n",
    "from sklearn.metrics import confusion_matrix\n",
    "from sklearn.cross_validation import LabelKFold\n",
    "from sklearn.cross_validation import LeaveOneLabelOut\n",
    "import matplotlib.pylab as plt\n",
    "%matplotlib inline\n",
    "from matplotlib.pylab import rcParams\n",
    "rcParams['figure.figsize'] = 8, 8\n",
    "import sys\n",
    "import pickle\n",
    "from graphlab.util import cloudpickle\n",
    "import graphlab.aggregate as agg\n",
    "from math import sqrt"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false,
    "deletable": true,
    "editable": true
   },
   "outputs": [],
   "source": [
    "# Load the data\n",
    "data =  gl.SFrame.read_csv('/net/fowler/vol1/home/vegray/metaDMS/dato/scaled/data/dmsTraining_2017-02-20.csv')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true,
    "deletable": true,
    "editable": true
   },
   "outputs": [],
   "source": [
    "## Use predictors3\n",
    "predictors3 = ['aa1', 'aa2', 'wt_mut', 'aa1_polarity', \n",
    "               'aa2_polarity', 'aa1_PI', 'aa2_PI', 'deltaPI', \n",
    "               'Grantham', 'aa1_weight', 'aa2_weight', 'deltaWeight',\n",
    "               'aa1vol', 'aa2vol', 'deltavolume', 'aa1_psic',\n",
    "               'aa2_psic', 'delta_psic',  'accessibility',\n",
    "               'dssp_sec_str', 'phi_psi_reg', 'delta_solvent_accessibility',\n",
    "               'b_factor', 'mut_msa_congruency', 'mut_mut_msa_congruency', \n",
    "               'seq_ind_closest_mut',  'evolutionary_coupling_avg']"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false,
    "deletable": true,
    "editable": true
   },
   "outputs": [],
   "source": [
    "## Global model dataset cleanup\n",
    "data2 = data.filter_by([x for x in data['dms_id'].unique() if x not in ['CcdB','beta-lactamase_2500', 'beta-lactamase_0', \n",
    "                                                                                      'kka2_1:1','kka2_1:4','kka2_1:8',\n",
    "                                                                     'beta-lactamase_39', 'gal4','hemagglutinin', 'np',\n",
    "                                                                                      'dbr1', 'beta-lactamase_156',\n",
    "                                                                                      'beta-lactamase_625', 'ERK2','GFP',\n",
    "                                                                                      'Brca1_HDR', ]], 'dms_id')\n",
    "#data2 = data2[(data2['binary_mut_type_prob'] < 0.4) | (data2['binary_mut_type_prob'] > 0.6)]\n",
    "data2 = data2.filter_by([x for x in data2['mut_type'].unique() if x in ['missense']], 'mut_type')\n",
    "#del data2['id']\n",
    "#data2 = data2.add_row_number('id',0)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false,
    "deletable": true,
    "editable": true
   },
   "outputs": [],
   "source": [
    "data2['dms_id'].unique()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false,
    "deletable": true,
    "editable": true
   },
   "outputs": [],
   "source": [
    "## Make data set w/o structural info\n",
    "\n",
    "data_str = data2\n",
    "\n",
    "data_str['accessibility'] = None\n",
    "data_str['dssp_sec_str'] = None\n",
    "data_str['dssp_sec_str'] = data_str['dssp_sec_str'].astype(str)\n",
    "data_str['phi_psi_reg'] = None\n",
    "data_str['phi_psi_reg'] = data_str['phi_psi_reg'].astype(str)\n",
    "data_str['delta_solvent_accessibility'] = None\n",
    "data_str['b_factor'] = None\n",
    "\n",
    "data_str['mut_mut_msa_congruency'] = data_str['mut_mut_msa_congruency'].astype(str)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true,
    "deletable": true,
    "editable": true
   },
   "outputs": [],
   "source": [
    "## Make data set w/o structural info\n",
    "\n",
    "data_bfact = data2\n",
    "\n",
    "data_bfact['accessibility'] = None\n",
    "data_bfact['dssp_sec_str'] = None\n",
    "data_bfact['dssp_sec_str'] = data_bfact['dssp_sec_str'].astype(str)\n",
    "data_bfact['phi_psi_reg'] = None\n",
    "data_bfact['phi_psi_reg'] = data_bfact['phi_psi_reg'].astype(str)\n",
    "data_bfact['delta_solvent_accessibility'] = None\n",
    "data_bfact['b_factor'] = None\n",
    "\n",
    "data_bfact['mut_mut_msa_congruency'] = data_bfact['mut_mut_msa_congruency'].astype(str)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {
    "collapsed": true,
    "deletable": true,
    "editable": true
   },
   "outputs": [],
   "source": [
    "## Make data set w/o evolutionary info\n",
    "\n",
    "data_evol = data2\n",
    "\n",
    "data_evol['aa1_psic'] = None\n",
    "data_evol['aa2_psic'] = None\n",
    "data_evol['delta_psic'] = None\n",
    "data_evol['mut_msa_congruency'] = None\n",
    "data_evol['mut_mut_msa_congruency'] = None\n",
    "data_evol['seq_ind_closest_mut'] = None\n",
    "data_evol['evolutionary_coupling_avg'] = None\n",
    "\n",
    "data_evol['mut_mut_msa_congruency'] = data_evol['mut_mut_msa_congruency'].astype(str)\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "deletable": true,
    "editable": true
   },
   "source": [
    "## KKA2"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {
    "collapsed": false,
    "deletable": true,
    "editable": true
   },
   "outputs": [],
   "source": [
    "unpickler = gl._gl_pickle.GLUnpickler(filename = '/net/fowler/vol1/home/vegray/metaDMS/dato/scaled/models/LOPO_models/NoKka2_model_regressionModel_noSELCO_limited_2017-02-08')\n",
    "NoKka2_model = unpickler.load()\n",
    "unpickler.close()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false,
    "deletable": true,
    "editable": true
   },
   "outputs": [],
   "source": [
    "## Missing structural features\n",
    "Kka2_str = data_str.filter_by([x for x in data_str['dms_id'].unique() if x in ['kka2_1:2']], 'dms_id')\n",
    "Kka2_str = Kka2_str.filter_by([x for x in Kka2_str['mut_type'].unique() if x not in ['synonymous']], 'mut_type')\n",
    "del Kka2_str['id']\n",
    "Kka2_str = Kka2_str.add_row_number('id',0)\n",
    "\n",
    "pred = NoKka2_model.predict(Kka2_str)\n",
    "\n",
    "stats.pearsonr(pred, Kka2_str['scaled_effect1'])\n",
    "stats.spearmanr(pred, Kka2_str['scaled_effect1'])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {
    "collapsed": false,
    "deletable": true,
    "editable": true
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "SpearmanrResult(correlation=0.49501403312345593, pvalue=0.0)"
      ]
     },
     "execution_count": 10,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "## Missing evolutionary features\n",
    "Kka2_evol = data_evol.filter_by([x for x in data_evol['dms_id'].unique() if x in ['kka2_1:2']], 'dms_id')\n",
    "Kka2_evol = Kka2_evol.filter_by([x for x in Kka2_evol['mut_type'].unique() if x not in ['synonymous']], 'mut_type')\n",
    "del Kka2_evol['id']\n",
    "Kka2_evol = Kka2_evol.add_row_number('id',0)\n",
    "\n",
    "pred = NoKka2_model.predict(Kka2_evol)\n",
    "\n",
    "\n",
    "stats.pearsonr(pred, Kka2_evol['scaled_effect1'])\n",
    "stats.spearmanr(pred, Kka2_evol['scaled_effect1'])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "## Missing bfact feature\n",
    "Kka2_bfact = data_bfact.filter_by([x for x in data_bfact['dms_id'].unique() if x in ['kka2_1:2']], 'dms_id')\n",
    "Kka2_bfact = Kka2_bfact.filter_by([x for x in Kka2_bfact['mut_type'].unique() if x not in ['synonymous']], 'mut_type')\n",
    "del Kka2_bfact['id']\n",
    "Kka2_bfact = Kka2_bfact.add_row_number('id',0)\n",
    "\n",
    "pred = NoKka2_model.predict(Kka2_bfact)\n",
    "\n",
    "\n",
    "stats.pearsonr(pred, Kka2_bfact['scaled_effect1'])\n",
    "stats.spearmanr(pred, Kka2_bfact['scaled_effect1'])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "deletable": true,
    "editable": true
   },
   "source": [
    "## PDZ domain"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {
    "collapsed": true,
    "deletable": true,
    "editable": true
   },
   "outputs": [],
   "source": [
    "unpickler = gl._gl_pickle.GLUnpickler(filename = '/net/fowler/vol1/home/vegray/metaDMS/dato/scaled/models/LOPO_models/NoPSD95pdz3_model_regressionModel_noSELCO_limited_2017-02-08')\n",
    "NoPSD95pdz3_model = unpickler.load()\n",
    "unpickler.close()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false,
    "deletable": true,
    "editable": true
   },
   "outputs": [],
   "source": [
    "## Missing structural features\n",
    "PSD95pdz3_str = data_str.filter_by([x for x in data_str['dms_id'].unique() if x in ['PSD95pdz3']], 'dms_id')\n",
    "PSD95pdz3_str = PSD95pdz3_str.filter_by([x for x in PSD95pdz3_str['mut_type'].unique() if x not in ['synonymous']], 'mut_type')\n",
    "del PSD95pdz3_str['id']\n",
    "PSD95pdz3_str = PSD95pdz3_str.add_row_number('id',0)\n",
    "\n",
    "pred = NoPSD95pdz3_model.predict(PSD95pdz3_str)\n",
    "\n",
    "stats.pearsonr(pred, PSD95pdz3_str['scaled_effect1'])\n",
    "stats.spearmanr(pred, PSD95pdz3_str['scaled_effect1'])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {
    "collapsed": false,
    "deletable": true,
    "editable": true
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "SpearmanrResult(correlation=0.46863990702842956, pvalue=6.6225104489696682e-87)"
      ]
     },
     "execution_count": 11,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "## Missing evolutionary features\n",
    "PSD95pdz3_evol = data_evol.filter_by([x for x in data_evol['dms_id'].unique() if x in ['PSD95pdz3']], 'dms_id')\n",
    "PSD95pdz3_evol = PSD95pdz3_evol.filter_by([x for x in PSD95pdz3_evol['mut_type'].unique() if x not in ['synonymous']], 'mut_type')\n",
    "del PSD95pdz3_evol['id']\n",
    "PSD95pdz3_evol = PSD95pdz3_evol.add_row_number('id',0)\n",
    "\n",
    "pred = NoPSD95pdz3_model.predict(PSD95pdz3_evol)\n",
    "stats.pearsonr(pred, PSD95pdz3_evol['scaled_effect1'])\n",
    "\n",
    "stats.pearsonr(pred, PSD95pdz3_evol['scaled_effect1'])\n",
    "stats.spearmanr(pred, PSD95pdz3_evol['scaled_effect1'])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "deletable": true,
    "editable": true
   },
   "source": [
    "## Pab1"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {
    "collapsed": true,
    "deletable": true,
    "editable": true
   },
   "outputs": [],
   "source": [
    "unpickler = gl._gl_pickle.GLUnpickler(filename = '/net/fowler/vol1/home/vegray/metaDMS/dato/scaled/models/LOPO_models/NoPab1_model_regressionModel_noSELCO_limited_2017-02-08')\n",
    "NoPab1_model = unpickler.load()\n",
    "unpickler.close()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false,
    "deletable": true,
    "editable": true
   },
   "outputs": [],
   "source": [
    "## Missing structural features\n",
    "pab1_str = data_str.filter_by([x for x in data_str['dms_id'].unique() if x in ['Pab1']], 'dms_id')\n",
    "pab1_str = pab1_str.filter_by([x for x in pab1_str['mut_type'].unique() if x not in ['synonymous']], 'mut_type')\n",
    "del pab1_str['id']\n",
    "pab1_str = pab1_str.add_row_number('id',0)\n",
    "\n",
    "pred = NoPab1_model.predict(pab1_str)\n",
    "\n",
    "stats.pearsonr(pred, pab1_str['scaled_effect1'])\n",
    "#stats.spearmanr(pred, pab1_str['scaled_effect1'])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {
    "collapsed": false,
    "deletable": true,
    "editable": true
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(0.49296129877674449, 9.3547008679135554e-74)"
      ]
     },
     "execution_count": 14,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "## Missing evolutionary features\n",
    "pab1_evol = data_evol.filter_by([x for x in data_evol['dms_id'].unique() if x in ['Pab1']], 'dms_id')\n",
    "pab1_evol = pab1_evol.filter_by([x for x in pab1_evol['mut_type'].unique() if x not in ['synonymous']], 'mut_type')\n",
    "del pab1_evol['id']\n",
    "pab1_evol = pab1_evol.add_row_number('id',0)\n",
    "\n",
    "pred = NoPab1_model.predict(pab1_evol)\n",
    "stats.pearsonr(pred, pab1_evol['scaled_effect1'])\n",
    "\n",
    "stats.pearsonr(pred, pab1_evol['scaled_effect1'])\n",
    "#stats.spearmanr(pred, pab1_evol['scaled_effect1'])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "deletable": true,
    "editable": true
   },
   "source": [
    "## Beta - lactamase"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {
    "collapsed": false,
    "deletable": true,
    "editable": true
   },
   "outputs": [],
   "source": [
    "unpickler = gl._gl_pickle.GLUnpickler(filename = '/net/fowler/vol1/home/vegray/metaDMS/dato/scaled/models/LOPO_models/Nobeta_lact_model_regressionModel_noSELCO_limited_2017-02-08')\n",
    "Nobeta_lact_model = unpickler.load()\n",
    "unpickler.close()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false,
    "deletable": true,
    "editable": true
   },
   "outputs": [],
   "source": [
    "## Missing structural features\n",
    "beta_lact_str = data_str.filter_by([x for x in data_str['dms_id'].unique() if x in ['beta-lactamase']], 'dms_id')\n",
    "beta_lact_str = beta_lact_str.filter_by([x for x in beta_lact_str['mut_type'].unique() if x not in ['synonymous']], 'mut_type')\n",
    "del beta_lact_str['id']\n",
    "beta_lact_str = beta_lact_str.add_row_number('id',0)\n",
    "\n",
    "pred = Nobeta_lact_model.predict(beta_lact_str)\n",
    "\n",
    "stats.pearsonr(pred, beta_lact_str['scaled_effect1'])\n",
    "#stats.spearmanr(pred, beta_lact_str['scaled_effect1'])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {
    "collapsed": false,
    "deletable": true,
    "editable": true
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(0.60320969082891507, 0.0)"
      ]
     },
     "execution_count": 17,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "## Missing evolutionary features\n",
    "beta_lact_evol = data_evol.filter_by([x for x in data_evol['dms_id'].unique() if x in ['beta-lactamase']], 'dms_id')\n",
    "beta_lact_evol = beta_lact_evol.filter_by([x for x in beta_lact_evol['mut_type'].unique() if x not in ['synonymous']], 'mut_type')\n",
    "del beta_lact_evol['id']\n",
    "beta_lact_evol = beta_lact_evol.add_row_number('id',0)\n",
    "\n",
    "pred = Nobeta_lact_model.predict(beta_lact_evol)\n",
    "stats.pearsonr(pred, beta_lact_evol['scaled_effect1'])\n",
    "\n",
    "stats.pearsonr(pred, beta_lact_evol['scaled_effect1'])\n",
    "#stats.spearmanr(pred, beta_lact_evol['scaled_effect1'])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "deletable": true,
    "editable": true
   },
   "source": [
    "## Ubiquitin"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "metadata": {
    "collapsed": true,
    "deletable": true,
    "editable": true
   },
   "outputs": [],
   "source": [
    "unpickler = gl._gl_pickle.GLUnpickler(filename = '/net/fowler/vol1/home/vegray/metaDMS/dato/scaled/models/LOPO_models/NoUba1_model_regressionModel_noSELCO_limited_2017-02-08')\n",
    "NoUba1_model = unpickler.load()\n",
    "unpickler.close()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false,
    "deletable": true,
    "editable": true
   },
   "outputs": [],
   "source": [
    "## Missing structural features\n",
    "E1_Ubiquitin_str = data_str.filter_by([x for x in data_str['dms_id'].unique() if x in ['E1_Ubiquitin']], 'dms_id')\n",
    "E1_Ubiquitin_str = E1_Ubiquitin_str.filter_by([x for x in E1_Ubiquitin_str['mut_type'].unique() if x not in ['synonymous']], 'mut_type')\n",
    "del E1_Ubiquitin_str['id']\n",
    "E1_Ubiquitin_str = E1_Ubiquitin_str.add_row_number('id',0)\n",
    "\n",
    "pred = NoUba1_model.predict(E1_Ubiquitin_str)\n",
    "\n",
    "stats.pearsonr(pred, E1_Ubiquitin_str['scaled_effect1'])\n",
    "#stats.spearmanr(pred, E1_Ubiquitin_str['scaled_effect1'])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 22,
   "metadata": {
    "collapsed": false,
    "deletable": true,
    "editable": true
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(0.53198854885273017, 2.6143410920958567e-80)"
      ]
     },
     "execution_count": 22,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "## Missing evolutionary features\n",
    "E1_Ubiquitin_evol = data_evol.filter_by([x for x in data_evol['dms_id'].unique() if x in ['E1_Ubiquitin']], 'dms_id')\n",
    "E1_Ubiquitin_evol = E1_Ubiquitin_evol.filter_by([x for x in E1_Ubiquitin_evol['mut_type'].unique() if x not in ['synonymous']], 'mut_type')\n",
    "del E1_Ubiquitin_evol['id']\n",
    "E1_Ubiquitin_evol = E1_Ubiquitin_evol.add_row_number('id',0)\n",
    "\n",
    "pred = NoUba1_model.predict(E1_Ubiquitin_evol)\n",
    "stats.pearsonr(pred, E1_Ubiquitin_evol['scaled_effect1'])\n",
    "\n",
    "stats.pearsonr(pred, E1_Ubiquitin_evol['scaled_effect1'])\n",
    "#stats.spearmanr(pred, E1_Ubiquitin_evol['scaled_effect1'])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false,
    "deletable": true,
    "editable": true
   },
   "outputs": [],
   "source": [
    "## Missing structural features\n",
    "Ubiquitin_str = data_str.filter_by([x for x in data_str['dms_id'].unique() if x in ['Ubiquitin']], 'dms_id')\n",
    "Ubiquitin_str = Ubiquitin_str.filter_by([x for x in Ubiquitin_str['mut_type'].unique() if x not in ['synonymous']], 'mut_type')\n",
    "del Ubiquitin_str['id']\n",
    "Ubiquitin_str = Ubiquitin_str.add_row_number('id',0)\n",
    "\n",
    "pred = NoUba1_model.predict(Ubiquitin_str)\n",
    "\n",
    "stats.pearsonr(pred, Ubiquitin_str['scaled_effect1'])\n",
    "#stats.spearmanr(pred, Ubiquitin_str['scaled_effect1'])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 23,
   "metadata": {
    "collapsed": false,
    "deletable": true,
    "editable": true
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "SpearmanrResult(correlation=0.40505687537971141, pvalue=1.6490824380950779e-50)"
      ]
     },
     "execution_count": 23,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "## Missing evolutionary features\n",
    "Ubiquitin_evol = data_evol.filter_by([x for x in data_evol['dms_id'].unique() if x in ['Ubiquitin']], 'dms_id')\n",
    "Ubiquitin_evol = Ubiquitin_evol.filter_by([x for x in Ubiquitin_evol['mut_type'].unique() if x not in ['synonymous']], 'mut_type')\n",
    "del Ubiquitin_evol['id']\n",
    "Ubiquitin_evol = Ubiquitin_evol.add_row_number('id',0)\n",
    "\n",
    "pred = NoUba1_model.predict(Ubiquitin_evol)\n",
    "stats.pearsonr(pred, Ubiquitin_evol['scaled_effect1'])\n",
    "\n",
    "stats.pearsonr(pred, Ubiquitin_evol['scaled_effect1'])\n",
    "stats.spearmanr(pred, Ubiquitin_evol['scaled_effect1'])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "deletable": true,
    "editable": true
   },
   "source": [
    "## Yap65"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 24,
   "metadata": {
    "collapsed": true,
    "deletable": true,
    "editable": true
   },
   "outputs": [],
   "source": [
    "unpickler = gl._gl_pickle.GLUnpickler(filename = '/net/fowler/vol1/home/vegray/metaDMS/dato/scaled/models/LOPO_models/NoYap65_model_regressionModel_noSELCO_limited_2017-02-08')\n",
    "NoYap65_model = unpickler.load()\n",
    "unpickler.close()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false,
    "deletable": true,
    "editable": true
   },
   "outputs": [],
   "source": [
    "## Missing structural features\n",
    "yap65_str = data_str.filter_by([x for x in data_str['dms_id'].unique() if x in ['WW_domain']], 'dms_id')\n",
    "yap65_str = yap65_str.filter_by([x for x in yap65_str['mut_type'].unique() if x not in ['synonymous']], 'mut_type')\n",
    "del yap65_str['id']\n",
    "yap65_str = yap65_str.add_row_number('id',0)\n",
    "\n",
    "pred = NoYap65_model.predict(yap65_str)\n",
    "\n",
    "stats.pearsonr(pred, yap65_str['scaled_effect1'])\n",
    "#stats.spearmanr(pred, yap65_str['scaled_effect1'])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 26,
   "metadata": {
    "collapsed": false,
    "deletable": true,
    "editable": true
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "SpearmanrResult(correlation=-0.027704038256774081, pvalue=0.59880947272263896)"
      ]
     },
     "execution_count": 26,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "## Missing evolutionary features\n",
    "yap65_evol = data_evol.filter_by([x for x in data_evol['dms_id'].unique() if x in ['WW_domain']], 'dms_id')\n",
    "yap65_evol = yap65_evol.filter_by([x for x in yap65_evol['mut_type'].unique() if x not in ['synonymous']], 'mut_type')\n",
    "del yap65_evol['id']\n",
    "yap65_evol = yap65_evol.add_row_number('id',0)\n",
    "\n",
    "pred = NoYap65_model.predict(yap65_evol)\n",
    "stats.pearsonr(pred, yap65_evol['scaled_effect1'])\n",
    "\n",
    "stats.pearsonr(pred, yap65_evol['scaled_effect1'])\n",
    "stats.spearmanr(pred, yap65_evol['scaled_effect1'])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "deletable": true,
    "editable": true
   },
   "source": [
    "## Hsp90 "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 27,
   "metadata": {
    "collapsed": true,
    "deletable": true,
    "editable": true
   },
   "outputs": [],
   "source": [
    "unpickler = gl._gl_pickle.GLUnpickler(filename = '/net/fowler/vol1/home/vegray/metaDMS/dato/scaled/models/LOPO_models/Nohsp90_model_regressionModel_noSELCO_limited_2017-02-08')\n",
    "Nohsp90_model = unpickler.load()\n",
    "unpickler.close()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false,
    "deletable": true,
    "editable": true
   },
   "outputs": [],
   "source": [
    "## Missing structural features\n",
    "hsp90_str = data_str.filter_by([x for x in data_str['dms_id'].unique() if x in ['hsp90']], 'dms_id')\n",
    "hsp90_str = hsp90_str.filter_by([x for x in hsp90_str['mut_type'].unique() if x not in ['synonymous']], 'mut_type')\n",
    "del hsp90_str['id']\n",
    "hsp90_str = hsp90_str.add_row_number('id',0)\n",
    "\n",
    "pred = Nohsp90_model.predict(hsp90_str)\n",
    "\n",
    "stats.pearsonr(pred, hsp90_str['scaled_effect1'])\n",
    "#stats.spearmanr(pred, hsp90_str['scaled_effect1'])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 29,
   "metadata": {
    "collapsed": false,
    "deletable": true,
    "editable": true
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(0.26119716418033989, 1.0302937476021656e-63)"
      ]
     },
     "execution_count": 29,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "## Missing evolutionary features\n",
    "hsp90_evol = data_evol.filter_by([x for x in data_evol['dms_id'].unique() if x in ['hsp90']], 'dms_id')\n",
    "hsp90_evol = hsp90_evol.filter_by([x for x in hsp90_evol['mut_type'].unique() if x not in ['synonymous']], 'mut_type')\n",
    "del hsp90_evol['id']\n",
    "hsp90_evol = hsp90_evol.add_row_number('id',0)\n",
    "\n",
    "pred = Nohsp90_model.predict(hsp90_evol)\n",
    "stats.pearsonr(pred, hsp90_evol['scaled_effect1'])\n",
    "\n",
    "stats.pearsonr(pred, hsp90_evol['scaled_effect1'])\n",
    "#stats.spearmanr(pred, hsp90_evol['scaled_effect1'])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "deletable": true,
    "editable": true
   },
   "source": [
    "## Protein G"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 30,
   "metadata": {
    "collapsed": true,
    "deletable": true,
    "editable": true
   },
   "outputs": [],
   "source": [
    "unpickler = gl._gl_pickle.GLUnpickler(filename = '/net/fowler/vol1/home/vegray/metaDMS/dato/scaled/models/LOPO_models/Nogb1_model_regressionModel_noSELCO_limited_2017-02-08')\n",
    "Nogb1_model = unpickler.load()\n",
    "unpickler.close()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false,
    "deletable": true,
    "editable": true
   },
   "outputs": [],
   "source": [
    "## Missing structural features\n",
    "gb1_str = data_str.filter_by([x for x in data_str['dms_id'].unique() if x in ['gb1']], 'dms_id')\n",
    "gb1_str = gb1_str.filter_by([x for x in gb1_str['mut_type'].unique() if x not in ['synonymous']], 'mut_type')\n",
    "del gb1_str['id']\n",
    "gb1_str = gb1_str.add_row_number('id',0)\n",
    "\n",
    "pred = Nogb1_model.predict(gb1_str)\n",
    "\n",
    "stats.pearsonr(pred, gb1_str['scaled_effect1'])\n",
    "#stats.spearmanr(pred, gb1_str['scaled_effect1'])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 32,
   "metadata": {
    "collapsed": false,
    "deletable": true,
    "editable": true
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(0.38433112371687506, 4.0442515474933188e-38)"
      ]
     },
     "execution_count": 32,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "## Missing evolutionary features\n",
    "gb1_evol = data_evol.filter_by([x for x in data_evol['dms_id'].unique() if x in ['gb1']], 'dms_id')\n",
    "gb1_evol = gb1_evol.filter_by([x for x in gb1_evol['mut_type'].unique() if x not in ['synonymous']], 'mut_type')\n",
    "del gb1_evol['id']\n",
    "gb1_evol = gb1_evol.add_row_number('id',0)\n",
    "\n",
    "pred = Nogb1_model.predict(gb1_evol)\n",
    "stats.pearsonr(pred, gb1_evol['scaled_effect1'])\n",
    "\n",
    "stats.pearsonr(pred, gb1_evol['scaled_effect1'])\n",
    "#stats.spearmanr(pred, gb1_evol['scaled_effect1'])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true,
    "deletable": true,
    "editable": true
   },
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 2",
   "language": "python",
   "name": "python2"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 2
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython2",
   "version": "2.7.3"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 0
}