S_DISCO_2024_Machine_Learning.ipynb

{
  "nbformat": 4,
  "nbformat_minor": 0,
  "metadata": {
    "colab": {
      "provenance": [],
      "authorship_tag": "ABX9TyN4yt9PRXS6JCGTZZTkLeJV",
      "include_colab_link": true
    },
    "kernelspec": {
      "name": "python3",
      "display_name": "Python 3"
    },
    "language_info": {
      "name": "python"
    }
  },
  "cells": [
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "view-in-github",
        "colab_type": "text"
      },
      "source": [
        "<a href=\"https://colab.research.google.com/github/cbedart/S-DISCO/blob/2024/S_DISCO_2024_Machine_Learning.ipynb\" target=\"_parent\"><img src=\"https://colab.research.google.com/assets/colab-badge.svg\" alt=\"Open In Colab\"/></a>"
      ]
    },
    {
      "cell_type": "markdown",
      "source": [
        "**<center><h1>S-DISCO 2024 - Machine Learning</h1></center>**\n",
        "\n",
        "---\n",
        "\n"
      ],
      "metadata": {
        "id": "0IZ4QKzDKt2W"
      }
    },
    {
      "cell_type": "markdown",
      "source": [
        "## <u>**Installation of mandatory packages**</u>"
      ],
      "metadata": {
        "id": "2htILZ6gK1pl"
      }
    },
    {
      "cell_type": "code",
      "source": [
        "# @title Installation of prerequisites\n",
        "!pip install rdkit\n",
        "\n",
        "import pandas as pd\n",
        "import numpy as np\n",
        "\n",
        "from rdkit import Chem\n",
        "from rdkit.Chem import AllChem, Descriptors, rdFingerprintGenerator, Draw\n",
        "from rdkit import RDLogger\n",
        "RDLogger.DisableLog('rdApp.*')\n",
        "\n",
        "from sklearn.decomposition import PCA, KernelPCA\n",
        "from sklearn import preprocessing\n",
        "from sklearn import decomposition\n",
        "from sklearn import datasets\n",
        "from sklearn.cluster import KMeans\n",
        "\n",
        "from sklearn import tree\n",
        "from sklearn.tree import export_graphviz, DecisionTreeClassifier\n",
        "from sklearn.metrics import confusion_matrix, ConfusionMatrixDisplay\n",
        "from sklearn.ensemble import RandomForestClassifier\n",
        "from sklearn.neural_network import MLPClassifier\n",
        "\n",
        "from sklearn.model_selection import train_test_split\n",
        "from sklearn.model_selection import cross_val_score\n",
        "from sklearn.model_selection import RepeatedStratifiedKFold\n",
        "from sklearn.model_selection import LeaveOneOut\n",
        "\n",
        "from sklearn.preprocessing import StandardScaler\n",
        "\n",
        "import ipywidgets as widgets\n",
        "from IPython.display import Image, display, clear_output, Javascript\n",
        "\n",
        "import matplotlib.pyplot as plt\n",
        "import seaborn as sns; sns.set()\n",
        "import plotly.express as px\n",
        "import pydotplus\n",
        "\n",
        "################################################################################\n",
        "\n",
        "launch_rf_hyperparameters = 0\n",
        "\n",
        "################################################################################\n",
        "\n",
        "def tree_graph_to_png(tree, feature_names, class_names, png_file_to_save):\n",
        "    tree_str = export_graphviz(tree, feature_names=feature_names, class_names=class_names,\n",
        "                                     filled=True, out_file=None)\n",
        "    graph = pydotplus.graph_from_dot_data(tree_str)\n",
        "    graph.write_png(png_file_to_save)\n",
        "\n",
        "\n",
        "def getMolDescriptors(mol, missingVal=None):\n",
        "    res = {}\n",
        "    for nm,fn in Descriptors._descList:\n",
        "        try:\n",
        "            val = fn(mol)\n",
        "        except:\n",
        "            val = missingVal\n",
        "        res[nm] = val\n",
        "    return res\n",
        "\n",
        "################################################################################\n",
        "\n",
        "print(\"\\n\\n\\033[1mPrerequisites succesfully installed !\\033[0m\")"
      ],
      "metadata": {
        "cellView": "form",
        "id": "PVzJ3M8hKnum"
      },
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "code",
      "source": [
        "# @title Loading data from the previous session\n",
        "\n",
        "!wget https://github.com/cbedart/S-DISCO/raw/refs/heads/2023/qsar_dict.pkl > /dev/null 2>&1\n",
        "qsar_dict = pd.read_pickle(\"qsar_dict.pkl\")\n",
        "\n",
        "################################################################################\n",
        "\n",
        "desc_GR_DTs = qsar_dict[\"desc_clustered_prep\"]\n",
        "\n",
        "pic50_limit = 8\n",
        "\n",
        "ygroup2 = []\n",
        "ycat2 = []\n",
        "for i in desc_GR_DTs.iloc[:,2]:\n",
        "  if i > pic50_limit:\n",
        "    ygroup2.append(\"Active\")\n",
        "    ycat2.append(1)\n",
        "  else:\n",
        "    ygroup2.append(\"Less Active\")\n",
        "    ycat2.append(0)\n",
        "\n",
        "desc_GR_activity2 = pd.concat([desc_GR_DTs.iloc[:,:4], pd.DataFrame(ygroup2, columns=[\"Activity_Label\"]), pd.DataFrame(ycat2, columns=[\"Activity\"]), desc_GR_DTs.iloc[:,4:]], axis = \"columns\")\n",
        "\n",
        "X2 = desc_GR_activity2.iloc[:,6:]\n",
        "y2 = desc_GR_activity2.loc[:,\"Activity\"]\n",
        "\n",
        "X_train2, X_test2, y_train2, y_test2 = train_test_split(X2, y2, test_size=0.2, random_state=42)\n",
        "\n",
        "################################################################################\n",
        "\n",
        "print(\"Consistency in splitting your data for training and validation is crucial\\nfor reliable model evaluation and comparison. Using the same split ensures\\nthat your evaluation metrics are consistent across different model iterations,\\nmaking it easier to compare their performance.\")\n",
        "print(\"\\nHere, we will do a single data split, with :\")\n",
        "print(\"- A train set containing 80% of the data ({} compounds)\".format(len(X_train2)))\n",
        "print(\"- A test set containing 20% of the data ({} compounds)\".format(len(X_test2)))\n",
        "\n",
        "print(\"\\nWe will also setup the pIC50 limit to {}, with :\".format(pic50_limit))\n",
        "print(\"- A total of {} active compounds (pIC50 < {}) - Coded as 1\".format(desc_GR_activity2[\"Activity\"].value_counts()[1], pic50_limit))\n",
        "print(\"- A total of {} less active compounds (pIC50 > {}) - Coded as 0\".format(desc_GR_activity2[\"Activity\"].value_counts()[0], pic50_limit))\n",
        "\n",
        "print(\"\\n\\033[1mDone !\\033[0m \\n\\n\")\n",
        "\n",
        "\n",
        "display(desc_GR_DTs)\n"
      ],
      "metadata": {
        "id": "Oqf3cx_1WX9z",
        "cellView": "form"
      },
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "source": [
        "## <u>**Model \\#3 - Decision tree**</u>"
      ],
      "metadata": {
        "id": "GNIfnvRAK3Jc"
      }
    },
    {
      "cell_type": "code",
      "source": [
        "# @title First decision tree\n",
        "\n",
        "dt1 = DecisionTreeClassifier()\n",
        "dt1.fit(X_train2, y_train2)\n",
        "\n",
        "y_pred_dt1 = dt1.predict(X_test2)\n",
        "# dt1_confusion = pd.DataFrame(confusion_matrix(y_test2, y_pred_dt1), columns=[\"Less Actives\", \"Actives\"], index=[\"Less Actives\", \"Actives\"])\n",
        "\n",
        "print(\"Accuracy on the train set = {}\\n\".format(round(dt1.score(X_train2, y_train2),3)))\n",
        "print(\"Accuracy on the test set = {}\".format(round(dt1.score(X_test2, y_test2),3)))\n",
        "\n",
        "cvk_dt1 = RepeatedStratifiedKFold(n_splits=5, n_repeats=3, random_state=1)\n",
        "scores_cvk_dt1 = cross_val_score(dt1, X2, y2, scoring='accuracy', cv=cvk_dt1, n_jobs=-1)\n",
        "\n",
        "print(\"Accuracy using 10-fold Cross Validation = {}\\n\".format(round(np.mean(scores_cvk_dt1),3)))\n",
        "\n",
        "\n",
        "tree_graph_to_png(dt1, feature_names=X2.columns, class_names=[\"Less Actives\", \"Actives\"], png_file_to_save='dt1.png')\n",
        "\n",
        "print(\"Confusion matrix:\")\n",
        "fig, axn = plt.subplots(1,2, figsize=(12,3))\n",
        "axn[0].grid(False) ; axn[1].grid(False)\n",
        "confusion1_dt1 = ConfusionMatrixDisplay.from_estimator(dt1, X_test2, y_test2, display_labels=[\"Less Active\", \"Active\"], ax=axn[0], cmap=plt.cm.Blues, normalize=None)\n",
        "confusion1_dt1.ax_.set_title(\"Confusion matrix, without normalization\")\n",
        "# confusion1_dt1\n",
        "confusion2_dt1 = ConfusionMatrixDisplay.from_estimator(dt1, X_test2, y_test2, display_labels=[\"Less Active\", \"Active\"], ax=axn[1], cmap=plt.cm.Blues, normalize=\"true\")\n",
        "confusion2_dt1.ax_.set_title(\"Normalized confusion matrix\")\n",
        "# confusion2_dt1\n",
        "plt.show()\n",
        "\n",
        "dt1_acc_pr_rec = [round(dt1.score(X_train2, y_train2),3),\n",
        "                  round(dt1.score(X_test2, y_test2),3),\n",
        "                  round(confusion2_dt1.confusion_matrix[1,1] / (confusion2_dt1.confusion_matrix[1,1] + confusion2_dt1.confusion_matrix[0,1]), 3),\n",
        "                  round(confusion2_dt1.confusion_matrix[1,1] / (confusion2_dt1.confusion_matrix[1,1] + confusion2_dt1.confusion_matrix[1,0]), 3)]"
      ],
      "metadata": {
        "cellView": "form",
        "id": "DMDL61OPZbXI"
      },
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "code",
      "source": [
        "# @title Display the first decision tree\n",
        "\n",
        "# @markdown Do not hesitate to right clic >>> \"Open in a new tab\" to better see the tree\n",
        "\n",
        "display(Image('dt1.png'))"
      ],
      "metadata": {
        "cellView": "form",
        "id": "a1IrbJMbaXxk"
      },
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "source": [
        "## <u>**Model \\#4 - Random forest**</u>"
      ],
      "metadata": {
        "id": "KKzpj86Vdamq"
      }
    },
    {
      "cell_type": "code",
      "source": [
        "# @title Build the first random forest model\n",
        "\n",
        "# @markdown Choose the number of decision trees you will generate:\n",
        "number_of_trees = 100 # @param {type:\"integer\"}\n",
        "\n",
        "# @markdown Choose if you set an early termination of the trees after X splits, or just build the bigger trees you can (= None)\n",
        "maximum_depth = \"None\" # @param [\"None\", 1, 2, 3, 4, 5]\n",
        "\n",
        "if maximum_depth == \"None\":\n",
        "  maximum_depth = None\n",
        "else:\n",
        "  maximum_depth = int(maximum_depth)\n",
        "\n",
        "print(\"\\n\\033[1mUsing {} decision trees with an early termination setup at {}\\033[0m\\n\".format(number_of_trees, maximum_depth))\n",
        "\n",
        "rf1 = RandomForestClassifier(n_estimators = number_of_trees, max_depth = maximum_depth, random_state=42)\n",
        "rf1.fit(X_train2, y_train2)\n",
        "\n",
        "y_pred_rf1 = rf1.predict(X_test2)\n",
        "rf1_confusion = pd.DataFrame(confusion_matrix(y_test2, y_pred_rf1), columns=[\"Less Actives\", \"Actives\"], index=[\"Less Actives\", \"Actives\"])\n",
        "\n",
        "print(\"Accuracy on the train set = {}\\n\".format(round(rf1.score(X_train2, y_train2),3)))\n",
        "print(\"Accuracy on the test set = {}\".format(round(rf1.score(X_test2, y_test2),3)))\n",
        "\n",
        "cvk_rf1 = RepeatedStratifiedKFold(n_splits=5, n_repeats = 3, random_state=1)\n",
        "scores_cvk_rf1 = cross_val_score(rf1, X2, y2, scoring='accuracy', cv=cvk_rf1, n_jobs=-1)\n",
        "\n",
        "print(\"Accuracy using 5-fold Cross Validation = {}\\n\".format(round(np.mean(scores_cvk_rf1),3)))\n",
        "\n",
        "print(\"Confusion matrix:\")\n",
        "fig, axn = plt.subplots(1,2, figsize=(12,3))\n",
        "axn[0].grid(False) ; axn[1].grid(False)\n",
        "confusion1_rf1 = ConfusionMatrixDisplay.from_estimator(rf1, X_test2, y_test2, display_labels=[\"Less Active\", \"Active\"], ax=axn[0], cmap=plt.cm.Blues, normalize=None)\n",
        "confusion1_rf1.ax_.set_title(\"Confusion matrix, without normalization\")\n",
        "# confusion1_rf1\n",
        "confusion2_rf1 = ConfusionMatrixDisplay.from_estimator(rf1, X_test2, y_test2, display_labels=[\"Less Active\", \"Active\"], ax=axn[1], cmap=plt.cm.Blues, normalize=\"true\")\n",
        "confusion2_rf1.ax_.set_title(\"Normalized confusion matrix\")\n",
        "# confusion2_rf1\n",
        "plt.show()\n",
        "\n",
        "rf1_acc_pr_rec = [round(rf1.score(X_train2, y_train2),3),\n",
        "                  round(rf1.score(X_test2, y_test2),3),\n",
        "                  round(confusion2_rf1.confusion_matrix[1,1] / (confusion2_rf1.confusion_matrix[1,1] + confusion2_rf1.confusion_matrix[0,1]), 3),\n",
        "                  round(confusion2_rf1.confusion_matrix[1,1] / (confusion2_rf1.confusion_matrix[1,1] + confusion2_rf1.confusion_matrix[1,0]), 3)]"
      ],
      "metadata": {
        "cellView": "form",
        "id": "iGOqesBqd6-D"
      },
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "code",
      "source": [
        "# @title Hyperparameters optimization of the random forest model\n",
        "\n",
        "if launch_rf_hyperparameters == 1:\n",
        "  raise FileExistsError('Already done')\n",
        "\n",
        "rf_accuracies_train = []\n",
        "rf_accuracies_test = []\n",
        "rf_precision_test = []\n",
        "rf_recall_test = []\n",
        "\n",
        "\n",
        "for maximum_depth_loop in [1, 2, 3, 4, 5, None]:\n",
        "  for number_of_trees_loop in [10,100,1000]:\n",
        "  # for number_of_trees_loop in [10,100,1000,10000]:\n",
        "    print(\"Currently running depth {} with {} trees\".format(maximum_depth_loop, number_of_trees_loop))\n",
        "    rfloop = RandomForestClassifier(n_estimators = number_of_trees_loop, max_depth = maximum_depth_loop, random_state=42)\n",
        "    rfloop.fit(X_train2, y_train2)\n",
        "    y_pred_rfloop = rfloop.predict(X_test2)\n",
        "    rfloop_confusion = confusion_matrix(y_test2, y_pred_rfloop, normalize=None)\n",
        "\n",
        "    rf_accuracies_train.append(round(rfloop.score(X_train2, y_train2),3))\n",
        "    rf_accuracies_test.append(round(rfloop.score(X_test2, y_test2),3))\n",
        "    rf_precision_test.append(round(rfloop_confusion[1,1] / (rfloop_confusion[1,1] + rfloop_confusion[0,1]), 3))\n",
        "    rf_recall_test.append(round(rfloop_confusion[1,1] / (rfloop_confusion[1,1] + rfloop_confusion[1,0]), 3))\n",
        "\n",
        "\n",
        "print(\"\\n\\033[1mDone !\\033[0m\")"
      ],
      "metadata": {
        "cellView": "form",
        "id": "vaSW7oEtfhVi"
      },
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "code",
      "source": [
        "# @title Plot for the hyperparameters optimization\n",
        "\n",
        "plt.figure(figsize=(10,5))\n",
        "plt.plot(rf_accuracies_train, marker=\".\")\n",
        "plt.plot(rf_accuracies_test, marker=\".\")\n",
        "plt.plot(rf_precision_test, marker=\".\")\n",
        "plt.plot(rf_recall_test, marker=\".\")\n",
        "plt.title(\"Quality assessment metrics for RF hyperparameters optimization\")\n",
        "plt.legend([\"Accuracy train\", \"Accuracy test\", \"Precision\", \"Recall\"])\n",
        "plt.ylim([-0.05,1.05])\n",
        "\n",
        "labels = []\n",
        "for maximum_depth_loop in [1, 2, 3, 4, 5, None]:\n",
        "  # for number_of_trees_loop in [10,100,1000,10000]:\n",
        "  for number_of_trees_loop in [10,100,1000]:\n",
        "    labels.append(\"{} / {}\".format(maximum_depth_loop, number_of_trees_loop))\n",
        "\n",
        "plt.xticks(range(18), labels, rotation='vertical')\n",
        "plt.show()"
      ],
      "metadata": {
        "cellView": "form",
        "id": "yv3WsRmBqLtc"
      },
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "source": [
        "## <u>**Model \\#5 - Artificial neural networks**</u>"
      ],
      "metadata": {
        "id": "ody6VUZYnlPT"
      }
    },
    {
      "cell_type": "code",
      "source": [
        "# @title Build your own artificial neural network\n",
        "display(Javascript('''google.colab.output.setIframeHeight(0, true, {maxHeight: 700})'''))\n",
        "\n",
        "# @markdown Number of hidden layers:\n",
        "nb_layers = 1 # @param [1,2,3,4,5] {type:\"raw\"}\n",
        "\n",
        "# @markdown Number of nodes per hidden layers:\n",
        "layer1 = 100 # @param {type:\"integer\"}\n",
        "layer2 = 100 # @param {type:\"integer\"}\n",
        "layer3 = 100 # @param {type:\"integer\"}\n",
        "layer4 = 100 # @param {type:\"integer\"}\n",
        "layer5 = 100 # @param {type:\"integer\"}\n",
        "\n",
        "# @markdown Activation function:\n",
        "activation_function = \"logistic\" # @param [\"identity\", \"logistic\", \"tanh\", \"relu\"] {type:\"string\"}\n",
        "\n",
        "# @markdown Solver for weight optimization :\n",
        "solver_opt = \"adam\" # @param [\"adam\", \"lbfgs\"] {type:\"string\"}\n",
        "\n",
        "################################################################################\n",
        "\n",
        "if nb_layers == 1:\n",
        "  if layer1 <= 0:\n",
        "    raise ValueError('Wrong setup - You have a wrong number nodes for some of hidden layers')\n",
        "  custom_layer_sizes=(layer1)\n",
        "elif nb_layers == 2:\n",
        "  if layer1 <= 0 or layer2 <= 0:\n",
        "    raise ValueError('Wrong setup - You have a wrong number nodes for some of hidden layers')\n",
        "  custom_layer_sizes=(layer1, layer2)\n",
        "elif nb_layers == 3:\n",
        "  if layer1 <= 0 or layer2 <= 0 or layer3 <= 0:\n",
        "    raise ValueError('Wrong setup - You have a wrong number nodes for some of hidden layers')\n",
        "  custom_layer_sizes=(layer1, layer2, layer3)\n",
        "elif nb_layers == 4:\n",
        "  if layer1 <= 0 or layer2 <= 0 or layer3 <= 0 or layer4 <= 0:\n",
        "    raise ValueError('Wrong setup - You have a wrong number nodes for some of hidden layers')\n",
        "  custom_layer_sizes=(layer1, layer2, layer3, layer4)\n",
        "else :\n",
        "  if layer1 <= 0 or layer2 <= 0 or layer3 <= 0 or layer4 <= 0 or layer5 <= 0:\n",
        "    raise ValueError('Wrong setup - You have a wrong number nodes for some of hidden layers')\n",
        "  custom_layer_sizes=(layer1, layer2, layer3, layer4, layer5)\n",
        "\n",
        "nn1 = MLPClassifier(hidden_layer_sizes=custom_layer_sizes, activation = activation_function, verbose = True, solver = solver_opt, alpha=1e-4, max_iter= 2000)\n",
        "# nn1 = MLPClassifier(hidden_layer_sizes=custom_layer_sizes, activation = 'identity', verbose = True, solver = 'lbfgs', alpha=1e-5, max_iter= 2000)\n",
        "nn1 = nn1.fit(X_train2, y_train2)\n",
        "y_pred_nn1 = nn1.predict(X_test2)\n",
        "\n",
        "################################################################################\n",
        "\n",
        "print(\"\\n\\n\\033[1mQuality assessment metrics - {} hidden layers / {} architecture:\\033[0m\\n\".format(nb_layers, custom_layer_sizes))\n",
        "print(\"Accuracy on the train set = {}\\n\".format(round(nn1.score(X_train2, y_train2),3)))\n",
        "print(\"Accuracy on the test set = {}\\n\".format(round(nn1.score(X_test2, y_test2),3)))\n",
        "\n",
        "print(\"Confusion matrix:\")\n",
        "fig, axn = plt.subplots(1,2, figsize=(12,3))\n",
        "axn[0].grid(False) ; axn[1].grid(False)\n",
        "confusion1_nn1 = ConfusionMatrixDisplay.from_estimator(nn1, X_test2, y_test2, display_labels=[\"Less Active\", \"Active\"], ax=axn[0], cmap=plt.cm.Blues, normalize=None)\n",
        "confusion1_nn1.ax_.set_title(\"Confusion matrix, without normalization\")\n",
        "# confusion1_nn1\n",
        "confusion2_nn1 = ConfusionMatrixDisplay.from_estimator(nn1, X_test2, y_test2, display_labels=[\"Less Active\", \"Active\"], ax=axn[1], cmap=plt.cm.Blues, normalize=\"true\")\n",
        "confusion2_nn1.ax_.set_title(\"Normalized confusion matrix\")\n",
        "# confusion2_nn1\n",
        "plt.show()\n",
        "\n",
        "print(\"\\nPrecision = {}\".format(round(confusion1_nn1.confusion_matrix[1,1] / (confusion1_nn1.confusion_matrix[1,1] + confusion1_nn1.confusion_matrix[0,1]), 3)))\n",
        "print(\"\\nRecall = {}\".format(round(confusion1_nn1.confusion_matrix[1,1] / (confusion1_nn1.confusion_matrix[1,1] + confusion1_nn1.confusion_matrix[1,0]), 3)))\n",
        "\n",
        "nn1_acc_pr_rec = [round(nn1.score(X_train2, y_train2),3),\n",
        "                  round(nn1.score(X_test2, y_test2),3),\n",
        "                  round(confusion1_nn1.confusion_matrix[1,1] / (confusion1_nn1.confusion_matrix[1,1] + confusion1_nn1.confusion_matrix[0,1]), 3),\n",
        "                  round(confusion1_nn1.confusion_matrix[1,1] / (confusion1_nn1.confusion_matrix[1,1] + confusion1_nn1.confusion_matrix[1,0]), 3)]"
      ],
      "metadata": {
        "id": "0myYVuR-n2BV",
        "cellView": "form"
      },
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "source": [
        "## <u>**New compounds prediction !**</u>"
      ],
      "metadata": {
        "id": "3Q2qnrU4vIR4"
      }
    },
    {
      "cell_type": "code",
      "source": [
        "# @title 4 new compounds\n",
        "\n",
        "# @markdown If you want to add a custom compound, add it here, or leave the field blank (default)\n",
        "custom_smiles = \"\" # @param {type:\"string\"}\n",
        "\n",
        "# CHEMBL4081121 - Ki:  0.0800nM  / Ki:  3.30nM\n",
        "# \"Cc1cccc(S(=O)(=O)N2CCC3=Cc4c(cnn4-c4ccc(F)cc4)C[C@]3(C(=O)c3cc(C)ccn3)C2)c1\"\n",
        "\n",
        "# CHEMBL247048 - Ki:  3.70E+3nM\n",
        "# \"O=C1C=C2CCN(Cc3ccccc3)CC2(Cc2ccccc2)CC1\"\n",
        "\n",
        "# CHEMBL4088286 - Ki: 0.120nM / Ki: 9.20nM\n",
        "# \"O=C(c1ccccn1)[C@]12Cc3cnn(-c4ccc(F)cc4)c3C=C1CCN(S(=O)(=O)c1cc(F)cc(F)c1)C2\"\n",
        "\n",
        "# CHEMBL3734774 - Ki:  0.160nM / Ki:  100nM\n",
        "# \"O=C(c1ccccn1)[C@]12Cc3cnn(-c4ccc(F)cc4)c3C=C1CCN(S(=O)(=O)c1ccc(C(F)(F)F)cc1)C2\"\n",
        "\n",
        "desc_new = pd.DataFrame(columns=[\"Compound ID\"] + [i for i,j in Descriptors._descList])\n",
        "\n",
        "if custom_smiles == \"\":\n",
        "    cpds_new = [\"Compound_1\", \"Compound_2\",\"Compound_3\",\"Compound_4\"]\n",
        "    smiles_new = [\"Cc1cccc(S(=O)(=O)N2CCC3=Cc4c(cnn4-c4ccc(F)cc4)C[C@]3(C(=O)c3cc(C)ccn3)C2)c1\",\n",
        "                \"O=C1C=C2CCN(Cc3ccccc3)CC2(Cc2ccccc2)CC1\",\n",
        "                \"O=C(c1ccccn1)[C@]12Cc3cnn(-c4ccc(F)cc4)c3C=C1CCN(S(=O)(=O)c1cc(F)cc(F)c1)C2\",\n",
        "                \"O=C(c1ccccn1)[C@]12Cc3cnn(-c4ccc(F)cc4)c3C=C1CCN(S(=O)(=O)c1ccc(C(F)(F)F)cc1)C2\"]\n",
        "else:\n",
        "    cpds_new = [\"Compound_1\", \"Compound_2\",\"Compound_3\",\"Compound_4\", \"Custom_compound\"]\n",
        "    smiles_new = [\"Cc1cccc(S(=O)(=O)N2CCC3=Cc4c(cnn4-c4ccc(F)cc4)C[C@]3(C(=O)c3cc(C)ccn3)C2)c1\",\n",
        "                \"O=C1C=C2CCN(Cc3ccccc3)CC2(Cc2ccccc2)CC1\",\n",
        "                \"O=C(c1ccccn1)[C@]12Cc3cnn(-c4ccc(F)cc4)c3C=C1CCN(S(=O)(=O)c1cc(F)cc(F)c1)C2\",\n",
        "                \"O=C(c1ccccn1)[C@]12Cc3cnn(-c4ccc(F)cc4)c3C=C1CCN(S(=O)(=O)c1ccc(C(F)(F)F)cc1)C2\",\n",
        "                custom_smiles]\n",
        "\n",
        "for i in range(len(smiles_new)):\n",
        "  mol_new = Chem.AddHs(Chem.MolFromSmiles(smiles_new[i]))\n",
        "  desc_new.loc[i] = [cpds_new[i]] + list(getMolDescriptors(mol_new, missingVal = np.NAN).values())\n",
        "\n",
        "if custom_smiles == \"\":\n",
        "    display(Draw.MolsToGridImage([Chem.MolFromSmiles(smiles_new[i]) for i in range(len(smiles_new))],molsPerRow=4,subImgSize=(200,200),legends=cpds_new))\n",
        "else:\n",
        "    display(Draw.MolsToGridImage([Chem.MolFromSmiles(smiles_new[i]) for i in range(len(smiles_new))],molsPerRow=5,subImgSize=(200,200),legends=cpds_new))\n",
        "\n",
        "display(desc_new)"
      ],
      "metadata": {
        "cellView": "form",
        "id": "1WDheb4vvOK0"
      },
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "code",
      "source": [
        "# @title Scaled using the same methods applied previously for the train/test dataset\n",
        "desc_train_normal = qsar_dict['desc_clustered']\n",
        "desc_train_normal_filtered = desc_train_normal.loc[desc_train_normal[\"BindingDB_ID\"].isin(desc_GR_DTs[\"BindingDB_ID\"])]\n",
        "desc_train_normal_filtered = desc_train_normal_filtered[desc_GR_DTs.columns]\n",
        "desc_train_normal_filtered\n",
        "\n",
        "scaler = qsar_dict[\"scaler\"]\n",
        "desc_new = desc_new[list(scaler.feature_names_in_)]\n",
        "# scaler.fit(desc_train_normal_filtered.iloc[:,4:])\n",
        "\n",
        "desc_new_scaled = pd.DataFrame(scaler.transform(desc_new))\n",
        "desc_new_scaled.columns = list(scaler.feature_names_in_)\n",
        "desc_new_scaled = desc_new_scaled[X_train2.columns.to_list()]\n",
        "\n",
        "desc_new_scaled = pd.concat([pd.DataFrame(cpds_new), desc_new_scaled], axis=\"columns\")\n",
        "desc_new_scaled.columns = [\"Compound ID\"] + X_train2.columns.to_list()\n",
        "desc_new_scaled"
      ],
      "metadata": {
        "id": "WpYWMc2Rz_eh",
        "cellView": "form"
      },
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "code",
      "source": [
        "# @title Predictions time\n",
        "\n",
        "from sklearn.discriminant_analysis import LinearDiscriminantAnalysis\n",
        "lda2 = LinearDiscriminantAnalysis()\n",
        "lda2 = lda2.fit(X_train2, y_train2)\n",
        "y_pred_lda2 = lda2.predict(X_test2)\n",
        "confusion2_lda2 = confusion_matrix(y_test2, y_pred_lda2, normalize=\"true\")\n",
        "lda2_acc_pr_rec = [round(lda2.score(X_train2, y_train2),3),\n",
        "                  round(lda2.score(X_test2, y_test2),3),\n",
        "                  round(confusion2_lda2[1,1] / (confusion2_lda2[1,1] + confusion2_lda2[0,1]), 3),\n",
        "                  round(confusion2_lda2[1,1] / (confusion2_lda2[1,1] + confusion2_lda2[1,0]), 3)]\n",
        "\n",
        "a1 = lda2.predict(desc_new_scaled.iloc[:,1:])\n",
        "a2 = dt1.predict(desc_new_scaled.iloc[:,1:])\n",
        "a3 = rf1.predict(desc_new_scaled.iloc[:,1:])\n",
        "a4 = nn1.predict(desc_new_scaled.iloc[:,1:])\n",
        "\n",
        "new_results = pd.DataFrame([a1, a2, a3, a4], index=[\"LDA\", \"DT\", \"RF\", \"NN\"], columns=cpds_new).replace(0, \"Less active\").replace(1, \"Active\")\n",
        "rf1_acc_pr_rec\n",
        "new_results[[\"Accuracy train\", \"Accuracy test\", \"Precision\", \"Recall\"]] = lda2_acc_pr_rec\n",
        "new_results.loc[\"LDA\",[\"Accuracy train\", \"Accuracy test\", \"Precision\", \"Recall\"]] = dt1_acc_pr_rec\n",
        "new_results.loc[\"DT\",[\"Accuracy train\", \"Accuracy test\", \"Precision\", \"Recall\"]] = dt1_acc_pr_rec\n",
        "new_results.loc[\"RF\",[\"Accuracy train\", \"Accuracy test\", \"Precision\", \"Recall\"]] = rf1_acc_pr_rec\n",
        "new_results.loc[\"NN\",[\"Accuracy train\", \"Accuracy test\", \"Precision\", \"Recall\"]] = nn1_acc_pr_rec\n",
        "new_results"
      ],
      "metadata": {
        "id": "4GYyfyxH-84G",
        "cellView": "form"
      },
      "execution_count": null,
      "outputs": []
    }
  ]
}