scoring_function_workflow.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Cretor: Daniel Castaneda Mogollon\n",
    "# Date: 02/09/2018 (date created)\n",
    "# Purpose: This workflow has been implemented with the purpose of converting PDB files to pdbqt formats, as well as\n",
    "# ligands. It also tests and predicts a score given different models (DL, NN, etc) and displays statistics parameters\n",
    "# leading to data validation."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "from __future__ import print_function\n",
    "\n",
    "import fileinput\n",
    "import glob\n",
    "import math\n",
    "import os\n",
    "import subprocess\n",
    "import re\n",
    "import matplotlib.pyplot as plt\n",
    "import pandas as pd\n",
    "from metk_dir.modelevaltoolkit import metk\n",
    "from metk_dir.modelevaltoolkit import metk_report\n",
    "from metk_dir.modelevaltoolkit import metk_plots\n",
    "from metk_dir.modelevaltoolkit import metk_util\n",
    "import sys\n",
    "from IPython.utils import io\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "def welcome():\n",
    "    print(\"Welcome to the first prototype of the UTEP Protein-ligand Workflow. This program contains the necessary tools\")\n",
    "    print(\"to run NN (neural network) and DL (deep learning) by using pdbqt files. In order for it to run, you require\")\n",
    "    print(\"a directory containing a folder for proteins and ligands, an experimental data file, and metk folder for the\")\n",
    "    print(\"proper metrics and statistical analysis.\")\n",
    "    print(\"\\n\")\n",
    "    print(\"Plese see the documentation for more information, or contact the developers in case a doubt arises\\n\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "def running_DL():\n",
    "    print(\"We will begin by running the NN and DL scoring functions.\")\n",
    "    path_to = input(\"Please specify the directory where your required files are located: \")\n",
    "    # Specifying the path to access to the set\n",
    "    path_prlig = (path_to+'/proteins_and_ligands')\n",
    "    os.chdir(path_to+'/proteins_and_ligands')  # Moving into that path\n",
    "    content = os.listdir(path_prlig)  # Storing all the content in that folder \n",
    "    match_list = []\n",
    "    for f in content:\n",
    "        match = re.findall('\\d\\w{3}$',f)  # This regex looks for all files that have exactly four alphanumeric\n",
    "        # characters (pdbID in my general-set folder\n",
    "        if match:\n",
    "            match_list.append(match[0])\n",
    "    os.chdir(path_to)\n",
    "    path_vina = (path_to+'/autodock_vina/bin/vina')\n",
    "           \n",
    "    for i in range(0, len(match_list), 1):\n",
    "        molecule = match_list[i]\n",
    "    \n",
    "        running_list = ['-receptor '+ path_prlig + '/' + str(molecule) + '/' + str(\n",
    "            molecule) + '_protein.pdbqt ' + '-ligand ' + path_prlig + '/' + str(molecule) + '/' + str(\n",
    "            molecule) + '_ligand.pdbqt -vina_executable '+path_vina]\n",
    "        \n",
    "        with io.capture_output() as captured:\n",
    "            %run  $command \n",
    "            command = ('NNScore2.py -receptor ' + path_prlig + '/' + str(molecule) + '/' + str(\n",
    "                molecule) + '_protein.pdbqt' + ' -ligand ' + path_prlig + '/' + str(molecule) + '/' + str(\n",
    "                molecule) + '_ligand.pdbqt -vina_executable ' + path_vina)\n",
    "        #print(captured.stdout)\n",
    "    return(path_to)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [],
   "source": [
    "def merging_csv_files(path2):\n",
    "    print(\"The DL scoring function has finalized, and the output should exist in a new folder generated \"\n",
    "          \"called output\")\n",
    "    print(\"\\n\")\n",
    "    print(\"Now we will merge the generated output files\")\n",
    "    path = path2+'/output'\n",
    "    os.chdir(path)\n",
    "    allFiles = glob.glob(path + \"/*.csv\")\n",
    "    fout = open(\"all_predicted_values.csv\", \"a\")\n",
    "    for file in allFiles:\n",
    "        f = open(file)\n",
    "        for line in f:\n",
    "            fout.write(line)\n",
    "            break\n",
    "        f.close()\n",
    "    fout.close()\n",
    "    return()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [],
   "source": [
    "def adding_header(path):\n",
    "    os.chdir(path+'/output')\n",
    "    headers = \"PDB_ID NN_Net_1 NN_Net_2 NN_Net_3 NN_Net_4 NN_Net_5 NN_Net_6 NN_Net_7 NN_Net_8 \" \\\n",
    "              \"NN_Net_9 NN_Net_10 NN_Net_11 NN_Net_12 NN_Net_13 NN_Net_14 NN_Net_15 NN_Net_16 \" \\\n",
    "                   \"NN_Net_17 NN_Net_18 NN_Net_19 NN_Net_20 DL_Net_1 DL_Net_2 DL_Net_3 DL_Net_4 \" \\\n",
    "                   \"DL_Net_5 DL_Net_6 DL_Net_7 DL_Net_8 DL_Net_9 DL_Net_10 DL_Net_11 DL_Net_12 \" \\\n",
    "                   \"DL_Net_13 DL_Net_14 DL_Net_15 DL_Net_16 DL_Net_17 DL_Net_18 DL_Net_19 DL_Net_20 \" \\\n",
    "                   \"DL_Net_21 DL_Net_22 DL_Net_23 DL_Net_24 DL_Net_25 DL_Net_26 DL_Net_27 DL_Net_28 \" \\\n",
    "                   \"DL_Net_29 DL_Net_30 DL_Net_31 DL_Net_32 DL_Net_33 DL_Net_34 DL_Net_35 DL_Net_36 \" \\\n",
    "                   \"DL_Net_37 DL_Net_38 DL_Net_39 DL_Net_40 DL_Net_41 DL_Net_42 DL_Net_43 DL_Net_44 \" \\\n",
    "                   \"DL_Net_45 DL_Net_46 DL_Net_47 DL_Net_48 DL_Net_49 DL_Net_50 DL_Net_51 DL_Net_52 \" \\\n",
    "                   \"DL_Net_53 DL_Net_54 DL_Net_55 DL_Net_56 DL_Net_57 DL_Net_58 DL_Net_59 DL_Net_60 \" \\\n",
    "                   \"DL_Net_61 DL_Net_62 DL_Net_63 DL_Net_64 DL_Net_65 DL_Net_66 DL_Net_67 DL_Net_68 \" \\\n",
    "                   \"DL_Net_69 DL_Net_70 DL_Net_71 DL_Net_72 DL_Net_73 DL_Net_74 DL_Net_75 DL_Net_76 \" \\\n",
    "                   \"DL_Net_77 DL_Net_78 DL_Net_79 DL_Net_80 DL_Net_81 DL_Net_82 DL_Net_83 DL_Net_84 \" \\\n",
    "                   \"DL_Net_85 DL_Net_86 DL_Net_87 DL_Net_88 DL_Net_89 DL_Net_90 DL_Net_91 DL_Net_92 \" \\\n",
    "                   \"DL_Net_93 DL_Net_94 DL_Net_95 DL_Net_96 DL_Net_97 DL_Net_98 DL_Net_99 DL_Net_100 \" \\\n",
    "                   \"Experimental Average_NN Average_DL StDev_NN StDev_DL MSE_NN MSE_DL RMSE_NN RMSE_DL \" \\\n",
    "                   \"Average_(fi-yi)_NN Average_(fi-yi)_DL Delta_G(kcal/mol)_NN Delta_G(kcal/mol)_DL \" \\\n",
    "              \"Delta_G(kcal/mol)_exp\".split()\n",
    "    for line in fileinput.input(files=['all_predicted_values.csv'], inplace=True):\n",
    "        if fileinput.isfirstline():\n",
    "            print(','.join(headers))\n",
    "        line = line.replace('\\n','')\n",
    "        print(line,)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [],
   "source": [
    "def adding_experimental_data(path):\n",
    "    id_list = []\n",
    "    id_exp = []\n",
    "    value_list = []\n",
    "    os.chdir(path+'/output')    \n",
    "    book = pd.read_csv(\"all_predicted_values.csv\")\n",
    "    os.chdir(path)\n",
    "    book_exp = pd.read_csv(\"experimental_pdb_binddata.csv\")\n",
    "    nrows = book.shape[0]\n",
    "    nrows_exp = book_exp.shape[0]\n",
    "    for i in range(0,nrows_exp,1):\n",
    "        id_exp.append(book_exp['pdbid'][i])             #Storing the ID's from the experimental file into a list\n",
    "        value_list.append(book_exp['act'][i])           #Storing the values matching the ID's into the list\n",
    "    for j in range(0,nrows,1):\n",
    "        id_list.append(book['PDB_ID'][j])\n",
    "        if id_list[j] in id_exp:\n",
    "            index_number = id_exp.index(id_list[j])\n",
    "            book.at[j,'Experimental'] = value_list[index_number]        #We add the experimental values to the dataframe\n",
    "    book.to_csv(\"all_predicted_values.csv\",sep=\",\")         #We write a new file with the edited dataframe"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [],
   "source": [
    "def statistics(path):\n",
    "    dataframe = pd.read_csv(path+'/all_predicted_values.csv')\n",
    "    nrows = dataframe.shape[0]\n",
    "    sum=0\n",
    "\n",
    "\n",
    "    #Average of NN Score\n",
    "    for i in range(0,nrows,1):\n",
    "        for j in range(1,21,1):                             #20 Networks given by NNScore\n",
    "            sum = sum+(dataframe['NN_Net_'+str(j)][i])\n",
    "        average = sum/20\n",
    "        dataframe.at[i,'Average_NN'] = average\n",
    "        sum=0\n",
    "\n",
    "    #Average of DL Score\n",
    "    for a in range(0,nrows,1):\n",
    "        for b in range(1,101,1):                            #100 Networks given by DL Score\n",
    "            sum = sum+(dataframe['DL_Net_'+str(b)][a])\n",
    "        average2 = sum/100\n",
    "        dataframe.at[a,'Average_DL'] = average2\n",
    "        sum=0\n",
    "\n",
    "    #Standard deviation of NN Score\n",
    "    df_NN = dataframe.iloc[:,2:22]\n",
    "    for c in range(0,nrows,1):\n",
    "        df_NN_row = df_NN.iloc[c]\n",
    "        df_NN_row_std = df_NN_row.std()\n",
    "        dataframe.at[c,'StDev_NN'] = df_NN_row_std\n",
    "\n",
    "    #Standard deviation of DL Score\n",
    "    df_DL = dataframe.iloc[:,22:122]\n",
    "    for d in range(0,nrows,1):\n",
    "        df_DL_row = df_DL.iloc[d]\n",
    "        df_DL_row_std = df_DL_row.std()\n",
    "        dataframe.at[d,'StDev_DL'] = df_DL_row_std\n",
    "\n",
    "    #Getting MSE for NN Score\n",
    "    sum_mse_NN = 0\n",
    "    for e in range(0,nrows,1):\n",
    "        for f in range(0,20,1):\n",
    "            df_NN_row = df_NN.iloc[e]\n",
    "            df_row= dataframe.iloc[e]\n",
    "            experimental_value = df_row[122]\n",
    "            if math.isnan(experimental_value):\n",
    "                dataframe.at[e,'MSE_NN'] = 'NaN'\n",
    "            else:\n",
    "                df_NN_value = df_NN_row[f]\n",
    "                df_NN_value = float(df_NN_value)\n",
    "                squared =((df_NN_value-experimental_value)*(df_NN_value-experimental_value))\n",
    "                sum_mse_NN = (sum_mse_NN + squared)\n",
    "                squared = 0\n",
    "                #print(sum_mse_NN)\n",
    "        mse_NN = sum_mse_NN*(0.05)\n",
    "        rmse_NN = math.sqrt(mse_NN)\n",
    "        dataframe.at[e,'MSE_NN'] = mse_NN\n",
    "        dataframe.at[e,'RMSE_NN'] = rmse_NN\n",
    "        sum_mse_NN=0\n",
    "\n",
    "    #Getting MSE for DL Score\n",
    "    sum_mse_DL = 0\n",
    "    for g in range(0, nrows, 1):\n",
    "        for h in range(0, 100, 1):\n",
    "            df_DL_row = df_DL.iloc[g]\n",
    "            df_row = dataframe.iloc[g]\n",
    "            experimental_value = df_row[122]\n",
    "            if math.isnan(experimental_value):\n",
    "                dataframe.at[g, 'MSE_DL'] = 'NaN'\n",
    "            else:\n",
    "                df_DL_value = df_DL_row[h]\n",
    "                df_DL_value = float(df_DL_value)\n",
    "                squared2 = ((df_DL_value - experimental_value) * (df_DL_value - experimental_value))\n",
    "                sum_mse_DL = (sum_mse_DL + squared2)\n",
    "                squared = 0\n",
    "        mse_DL = sum_mse_DL * (0.001)\n",
    "        rmse_DL = math.sqrt(mse_DL)\n",
    "        dataframe.at[g, 'MSE_DL'] = mse_DL\n",
    "        dataframe.at[g, 'RMSE_DL'] = rmse_DL\n",
    "        sum_mse_DL = 0\n",
    "\n",
    "    #Getting (fi-yi) for NN (average)\n",
    "    for k in range(0,nrows,1):\n",
    "        df_row = dataframe.iloc[k]\n",
    "        average_NN = df_row[123]\n",
    "        experimental_value = df_row[122]\n",
    "        if math.isnan(experimental_value):\n",
    "            dataframe.at[k, 'Average_(fi-yi)_NN'] = 'NaN'\n",
    "        else:\n",
    "            distance_NN = average_NN-experimental_value\n",
    "            dataframe.at[k, 'Average_(fi-yi)_NN'] = distance_NN\n",
    "\n",
    "    #Getting (fi-yi) for DL (average)\n",
    "    for m in range(0,nrows,1):\n",
    "        df_row = dataframe.iloc[m]\n",
    "        average_DL = df_row[124]\n",
    "        experimental_value = df_row[122]\n",
    "        if math.isnan(experimental_value):\n",
    "            dataframe.at[m, 'Average_(fi-yi)_DL'] = 'NaN'\n",
    "        else:\n",
    "            distance_DL = average_DL-experimental_value\n",
    "            dataframe.at[m, 'Average_(fi-yi)_DL'] = distance_DL\n",
    "\n",
    "\n",
    "    #Getting the delta_G values\n",
    "    for n in range(0,nrows,1):\n",
    "        df_row = dataframe.iloc[n]\n",
    "        experimental_value = df_row[122]\n",
    "        average_NN = df_row[123]\n",
    "        average_DL = df_row[124]\n",
    "        kd_nn = (10**(-average_NN))\n",
    "        kd_dl = (10**(-average_DL))\n",
    "        kd_experimental = (10**(-experimental_value))\n",
    "        delta_g_nn = (0.59248368)*math.log(kd_nn,math.e)\n",
    "        delta_g_dl = (0.59248368)*math.log(kd_dl,math.e)\n",
    "        delta_g_exp =(0.59248368)*math.log(kd_experimental,math.e)\n",
    "        dataframe.at[n, 'Delta_G(kcal/mol)_NN'] = delta_g_nn\n",
    "        dataframe.at[n, 'Delta_G(kcal/mol)_DL'] = delta_g_dl\n",
    "        dataframe.at[n, 'Delta_G(kcal/mol)_exp'] = delta_g_exp\n",
    "\n",
    "    #Creating two input files for Walters metrics, one for nn and one for dl scoring functions\n",
    "    file = open('input_walters_nn.csv','a')\n",
    "    file2 = open('input_walters_dl.csv','a')\n",
    "    file.write(\"Name,Exp,Pred\")\n",
    "    file2.write(\"Name,Exp,Pred\")\n",
    "    file.close()\n",
    "    file2.close()\n",
    "    walters_df_nn = pd.read_csv('input_walters_nn.csv')\n",
    "    walters_df_dl = pd.read_csv('input_walters_dl.csv')\n",
    "    nrows2 = dataframe.shape[0]\n",
    "    for m in range(0,nrows2,1):\n",
    "        df_row = dataframe.iloc[m]\n",
    "        pdb_id = df_row[1]\n",
    "        experimental_value = df_row[135]\n",
    "        pred_nn = df_row[133]\n",
    "        pred_dl = df_row[134]\n",
    "        walters_df_nn.at[m,'Name'] = pdb_id\n",
    "        walters_df_nn.at[m,'Exp'] = experimental_value\n",
    "        walters_df_nn.at[m,'Pred']=pred_nn\n",
    "        walters_df_dl.at[m,'Name'] = pdb_id\n",
    "        walters_df_dl.at[m,'Exp'] = experimental_value\n",
    "        walters_df_dl.at[m,'Pred'] = pred_dl\n",
    "    walters_df_nn.to_csv('input_walters_nn.csv',index = False,sep=',')\n",
    "    walters_df_dl.to_csv('input_walters_dl.csv',index = False,sep=',')\n",
    "\n",
    "    #Displaying plot\n",
    "    y_list = dataframe['Experimental']\n",
    "    x_list_nn = dataframe['Average_NN']\n",
    "    x_list_dl = dataframe['Average_DL']\n",
    "    axes = plt.gca()\n",
    "    y_limit = plt.ylim(2,12)\n",
    "    x_limit = plt.xlim(2,12)\n",
    "    plt.title(\"Experimental vs Predicted pKd values\")\n",
    "    plt.plot([-15,35],[-15,35],'r')\n",
    "    NN_plot = plt.scatter(x_list_nn,y_list,marker='x',color='b')\n",
    "    DL_plot = plt.scatter(x_list_dl,y_list,marker='o',color='g')\n",
    "    plt.ylabel(\"Experimental pKd values\")\n",
    "    plt.xlabel(\"Predicted pKd values\")\n",
    "    plt.legend((NN_plot,DL_plot),('NN Score','DL Score'),scatterpoints=1,loc=4,prop={'size':9})\n",
    "    plt.show()\n",
    "    dataframe.to_csv(\"all_predicted_values.csv\",sep=',')\n",
    "\n",
    "    return()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [],
   "source": [
    "def walters_metrics(path):\n",
    "    print('Metrics analyzed. Output files generated as: metrics_output_nn and metrics_output_dl')\n",
    "    metk.calling('input_walters_nn.csv','metrics_output_nn','uM')\n",
    "    metk.calling('input_walters_dl.csv','metrics_output_dl','uM')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "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.5.4"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}