2e-DT2Vec_develope_ML_model.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import random\n",
    "from collections import defaultdict\n",
    "from pprint import pprint\n",
    "import matplotlib.pyplot as plt\n",
    "from sklearn.model_selection import KFold\n",
    "import itertools\n",
    "from sklearn.manifold import TSNE\n",
    "from sknetwork.data import karate_club, painters, movie_actor\n",
    "from sknetwork.clustering import Louvain, BiLouvain, modularity, bimodularity\n",
    "from sknetwork.linalg import normalize\n",
    "from sknetwork.utils import bipartite2undirected, membership_matrix\n",
    "from sknetwork.visualization import svg_graph, svg_digraph, svg_bigraph\n",
    "from sknetwork.hierarchy import LouvainHierarchy, BiLouvainHierarchy\n",
    "from sknetwork.hierarchy import dasgupta_score, dasgupta_cost\n",
    "from sklearn.cluster import MiniBatchKMeans, KMeans, FeatureAgglomeration, AffinityPropagation, MeanShift, SpectralClustering, AgglomerativeClustering, AgglomerativeClustering, DBSCAN, Birch\n",
    "import graphviz\n",
    "from xgboost import plot_tree\n",
    "import csv\n",
    "import os\n",
    "import json\n",
    "import gensim\n",
    "import numpy as np\n",
    "from pprint import pprint\n",
    "import pandas as pd\n",
    "import seaborn as sns\n",
    "from sklearn.decomposition import PCA\n",
    "from ggplot import *\n",
    "from scipy.stats import chi2_contingency\n",
    "from matplotlib import rc\n",
    "import scipy\n",
    "import pickle\n",
    "import shap\n",
    "import requests\n",
    "import json\n",
    "import re\n",
    "from os.path import isfile, join\n",
    "from os import listdir\n",
    "from bs4 import BeautifulSoup\n",
    "from collections import Counter\n",
    "from ggplot import *\n",
    "from sklearn import metrics\n",
    "from sklearn.preprocessing import MinMaxScaler\n",
    "from matplotlib import pyplot\n",
    "import networkx as nx\n",
    "import networkx\n",
    "import community\n",
    "from yellowbrick.cluster import KElbowVisualizer\n",
    "from community import community_louvain\n",
    "from xgboost import XGBClassifier, plot_importance\n",
    "from sklearn.metrics import (\n",
    "    precision_score,\n",
    "    recall_score,\n",
    "    f1_score,\n",
    "    confusion_matrix,\n",
    "    matthews_corrcoef,\n",
    "    auc,\n",
    "    average_precision_score,\n",
    "    fbeta_score)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# for an imbalance data change it to True\n",
    "np.random.seed(42)\n",
    "SCALE_POS_RATIO = False\n",
    "HIGHLY_POS = 5.5 #positive threshold\n",
    "HIGHLY_NEG = 0\n",
    "MAIN_DIR = '.'"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Functions"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def community_layout(g, partition):\n",
    "    \"\"\"\n",
    "    Compute the layout for a modular graph.\n",
    "\n",
    "\n",
    "    Arguments:\n",
    "    ----------\n",
    "    g -- networkx.Graph or networkx.DiGraph instance\n",
    "        graph to plot\n",
    "\n",
    "    partition -- dict mapping int node -> int community\n",
    "        graph partitions\n",
    "\n",
    "\n",
    "    Returns:\n",
    "    --------\n",
    "    pos -- dict mapping int node -> (float x, float y)\n",
    "        node positions\n",
    "\n",
    "    \"\"\"\n",
    "\n",
    "    pos_communities = _position_communities(g, partition, scale=3.)\n",
    "\n",
    "    pos_nodes = _position_nodes(g, partition, scale=1.)\n",
    "\n",
    "    # combine positions\n",
    "    pos = dict()\n",
    "    for node in g.nodes():\n",
    "        pos[node] = pos_communities[node] + pos_nodes[node]\n",
    "\n",
    "    return pos\n",
    "\n",
    "def _position_communities(g, partition, **kwargs):\n",
    "\n",
    "    # create a weighted graph, in which each node corresponds to a community,\n",
    "    # and each edge weight to the number of edges between communities\n",
    "    between_community_edges = _find_between_community_edges(g, partition)\n",
    "\n",
    "    communities = set(partition.values())\n",
    "    hypergraph = nx.DiGraph()\n",
    "    hypergraph.add_nodes_from(communities)\n",
    "    for (ci, cj), edges in between_community_edges.items():\n",
    "        hypergraph.add_edge(ci, cj, weight=len(edges))\n",
    "\n",
    "    # find layout for communities\n",
    "    pos_communities = nx.spring_layout(hypergraph, **kwargs)\n",
    "\n",
    "    # set node positions to position of community\n",
    "    pos = dict()\n",
    "    for node, community in partition.items():\n",
    "        pos[node] = pos_communities[community]\n",
    "\n",
    "    return pos\n",
    "\n",
    "def _find_between_community_edges(g, partition):\n",
    "\n",
    "    edges = dict()\n",
    "\n",
    "    for (ni, nj) in g.edges():\n",
    "        ci = partition[ni]\n",
    "        cj = partition[nj]\n",
    "\n",
    "        if ci != cj:\n",
    "            try:\n",
    "                edges[(ci, cj)] += [(ni, nj)]\n",
    "            except KeyError:\n",
    "                edges[(ci, cj)] = [(ni, nj)]\n",
    "\n",
    "    return edges\n",
    "\n",
    "def _position_nodes(g, partition, **kwargs):\n",
    "    \"\"\"\n",
    "    Positions nodes within communities.\n",
    "    \"\"\"\n",
    "\n",
    "    communities = dict()\n",
    "    for node, community in partition.items():\n",
    "        try:\n",
    "            communities[community] += [node]\n",
    "        except KeyError:\n",
    "            communities[community] = [node]\n",
    "\n",
    "    pos = dict()\n",
    "    for ci, nodes in communities.items():\n",
    "        subgraph = g.subgraph(nodes)\n",
    "        pos_subgraph = nx.spring_layout(subgraph, **kwargs)\n",
    "        pos.update(pos_subgraph)\n",
    "\n",
    "    return pos"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def embedding(DDS, PPS_seq, mapping=True):\n",
    "    \"\"\"\n",
    "    This function maps drug-drug similarities and protein-protein similarities to vectors\n",
    "    \n",
    "    Args: \n",
    "    PPS_seq: A DataFrame of protein-protein smilarities \n",
    "    DDS: A DataFrame of drug-drug smilarities \n",
    "    mapping: Binary (False/True), read from the saved files (mapped before)\n",
    "    \n",
    "    Returns: A dataframe vectors of drugs and proteins \n",
    "    \"\"\"\n",
    "    \n",
    "    PPS_seq = PPS_seq[PPS_seq['weight']!=0]\n",
    "    DDS = DDS[DDS['weight']!=0]\n",
    "    \n",
    "    DDS.to_csv('edglist_drugs.edgelist', sep=' ', index=False, header=False)\n",
    "    PPS_seq.to_csv('edglist_proteins.edgelist', sep=' ', index=False, header=False)\n",
    "    \n",
    "    if mapping:\n",
    "        # nod2vec (for drug)\n",
    "        os.system(f'PYTHONHASHSEED=10 python2 node2vec/src/main.py --workers 8 --input edglist_drugs.edgelist --output dim100_drugs.emb --weighted --dimensions 100')\n",
    "         # nod2vec (for proteins)\n",
    "        os.system(f'PYTHONHASHSEED=10 python2 node2vec/src/main.py --workers 8 --input edglist_proteins.edgelist --output dim100_proteins.emb --weighted --dimensions 100')\n",
    "\n",
    "    embeddings_seq_drug = pd.read_csv('dim100_drugs.emb', sep=' ', skiprows=[0], header=None, index_col=0)  \n",
    "    embeddings_seq_protein = pd.read_csv('dim100_proteins.emb', sep=' ', skiprows=[0], header=None, index_col=0)  \n",
    "\n",
    "    embeddings_seq = embeddings_seq_drug.append(embeddings_seq_protein)\n",
    "    embeddings_seq.index.name = 'ID'\n",
    "    \n",
    "    return embeddings_seq, embeddings_seq_drug, embeddings_seq_protein"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def calculate_metrics(y, y_pred):\n",
    "    \"\"\"\n",
    "    This function calculates machine learning metrics\n",
    "    \n",
    "    Args: Real lables and predicted labels\n",
    "\n",
    "    Returns: A dictionary containing all the results\n",
    "    \"\"\"\n",
    "    \n",
    "    result = {\n",
    "        'accuracy': metrics.accuracy_score(y, y_pred),\n",
    "        'precision': precision_score(y, y_pred, average='binary'),\n",
    "        'recall': recall_score(y, y_pred, average='binary'),\n",
    "        'f1_score': f1_score(y, y_pred, average='binary'),\n",
    "        'average_precision_score': average_precision_score(y, y_pred),\n",
    "        'f2_score':fbeta_score(y, y_pred, beta=2)\n",
    "    }\n",
    "    \n",
    "    try:\n",
    "        result['ROC'] = metrics.roc_auc_score(y, y_pred)\n",
    "    except:\n",
    "        pass\n",
    "\n",
    "    return result"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def calculate_metrics_multiclass(y, y_pred):\n",
    "    \"\"\"\n",
    "    This function calculates machine learning metrics\n",
    "    \n",
    "    Args: Real lables and predicted labels\n",
    "\n",
    "    Returns: A dictionary containing all the results\n",
    "    \"\"\"\n",
    "    \n",
    "    result = {\n",
    "        'accuracy': metrics.accuracy_score(y, y_pred),\n",
    "        'precision': precision_score(y, y_pred, average='samples'),\n",
    "        'recall': recall_score(y, y_pred, average='samples'),\n",
    "        'f1_score': f1_score(y, y_pred, average='samples'),\n",
    "    }\n",
    "    \n",
    "\n",
    "    return result"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def concating(df, embedding):\n",
    "    \"\"\"\n",
    "    This function concats drugs and targets vectors\n",
    "    \n",
    "    Args: \n",
    "        df: A DataFrame of drug-target interactions\n",
    "        embedding: A DataFrame of embedded vectors for each drug and target\n",
    "        \n",
    "    Returns: A DataFrame of drug-target interaction vectors\n",
    "    \"\"\"\n",
    "    dataset = []\n",
    "    for idx, row in df.iterrows():\n",
    "        from_vector = embedding.loc[row['from']]\n",
    "        to_vector = embedding.loc[row['to']]\n",
    "        features = from_vector.append(to_vector).reset_index(drop=True)\n",
    "        features = features.append(row)\n",
    "        dataset.append(features)\n",
    "\n",
    "    df_final = pd.DataFrame(dataset)\n",
    "    df_final.drop(columns=['from', 'to'], inplace=True)\n",
    "    return df_final"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def get_data(df):\n",
    "    \"\"\"\n",
    "    This function separates features (as x) and labels (as y)\n",
    "    \n",
    "    Args: A Datafram containing features and labels\n",
    "        \n",
    "    Returns: \n",
    "        x: Features\n",
    "        y: labels\n",
    "    \"\"\"\n",
    "    \n",
    "    df['label'] = df['weight']\n",
    "    df.drop(columns = ['weight'], inplace=True)\n",
    "    print(df['label'].value_counts())\n",
    "    \n",
    "    X = df.drop(columns=['label'])\n",
    "    y = pd.DataFrame(df['label'])\n",
    "    \n",
    "    return X, y.values.ravel()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def get_sample_weight(y_train):\n",
    "    \"\"\"\n",
    "    This function calculates the weigh_ratio for training imbalance dataset\n",
    "    \n",
    "    Args: Labels of train-set\n",
    "    \n",
    "    Returns: The weight ratio\n",
    "    \"\"\"\n",
    "    weight_ratio = float(len(y_train[y_train == 0]))/float(len(y_train[y_train == 1]))\n",
    "    \n",
    "    return weight_ratio"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def normalized_data(df):\n",
    "    df_n = df.copy()\n",
    "    mean_weight = df_n['weight'].mean()\n",
    "    std_weight = df_n['weight'].std()\n",
    "    df_n['weight_new'] = (df_n['weight']-mean_weight)/std_weight\n",
    "    df_n = df_n.drop(columns=[\"weight\"])\n",
    "    df_n = df_n.rename(columns={'weight_new':'weight'})\n",
    "    return df_n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def scale_data(df):\n",
    "    scaler = MinMaxScaler(feature_range=(0, 1))\n",
    "    data = df[['weight']]\n",
    "    scaler.fit(data)\n",
    "    data = scaler.transform(data)\n",
    "    df['weight'] = data\n",
    "    return df"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "with open('all_nod2vec_new.pkl', 'rb') as f:\n",
    "    all_nod2vec_new = pickle.load(f)\n",
    "all_nod2vec_new"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 1- Reading input data "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# DTI\n",
    "DTI= pd.read_csv(f'{MAIN_DIR}edgelis/DTI.csv', index_col=0)\n",
    "num_drug = len(set(list(DTI['from'])))\n",
    "num_protein = len(set(list(DTI['to'])))\n",
    "\n",
    "print(f'Number of drug: {num_drug} and number of protein: {num_protein}')\n",
    "DTI.head()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# PPS\n",
    "PPS_seq = pd.read_csv(f'{MAIN_DIR}edgelis/PPS_seq.csv', index_col=0)\n",
    "#PPS_seq = normalized_data(PPS_seq)\n",
    "#PPS_seq = scale_data(PPS_seq)\n",
    "PPS_seq.head()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "with open('protein2num.pkl', 'rb') as f:\n",
    "    protein2num = pickle.load(f)\n",
    "num2protein = dict((str(int(y)),x) for x,y in protein2num.items())"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "PPS_seq['weight'].mean()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "PPS_seq_new = PPS_seq[PPS_seq['weight']>=0.001]\n",
    "PPS_seq_new = PPS_seq.copy()\n",
    "PPS_seq_new['weight'].mean()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# clusters PPS network\n",
    "edgeList_pps = PPS_seq_new.values.tolist()\n",
    "G = networkx.Graph()\n",
    "weights = []\n",
    "\n",
    "for i in range(len(edgeList_pps)):\n",
    "    G.add_edge(edgeList_pps[i][0], edgeList_pps[i][1], weight=edgeList_pps[i][2])\n",
    "    weights.append(edgeList_pps[i][2])\n",
    "    \n",
    "A = networkx.adjacency_matrix(G).A\n",
    "PPS_adj = A.copy()\n",
    "\n",
    "louvain = Louvain()\n",
    "labels = louvain.fit_transform(PPS_adj)\n",
    "\n",
    "labels_unique, counts = np.unique(labels, return_counts=True)\n",
    "\n",
    "PPS_cluster_label= pd.DataFrame({'target':list(G.nodes()), 'label':labels})\n",
    "PPS_cluster_label['label']= PPS_cluster_label['label'].astype(str)\n",
    "print(labels_unique, counts)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": false
   },
   "outputs": [],
   "source": [
    "target2cluster = dict(zip(PPS_cluster_label.target, PPS_cluster_label.label))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "modularity2 = community.modularity(target2cluster, G, weight='weight')\n",
    "print(\"The modularity Q based on networkx is {}\".format(modularity2))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "color_dict = {k: v for k, v in enumerate(['#58ACFA','#FF1493', 'yellow','orange', '#00CED1','#5F9EA0','#006400','#96bf65','#fcc808','#7b2b48',\n",
    " '#e96957','#e06000','#173679','#d2dd49','#684a6b','#096eb2','#ce482a', 'red', 'lime', 'lightslategray',\n",
    "                                      'olive', 'rosybrown', 'sienna', 'darkmagenta','midnightblue','maroon',\n",
    "                                      'lightcoral','gold','sandybrown','tomato','lawngreen','lightgreen','darkorchid',\n",
    "                                      'lightskyblue','darkgreen'])}\n",
    "color_dict= {str(key): value for key, value in color_dict.items()}\n",
    "\n",
    "with open('color_dict_target.pkl', 'wb') as fp:\n",
    "    pickle.dump(color_dict, fp)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "partition = community_louvain.best_partition(G)\n",
    "pos = community_layout(G, target2cluster)\n",
    "plt.figure(figsize=(9,5))\n",
    "nx.draw(G, pos, node_color=[color_dict[v] for v in target2cluster.values()], edge_color=weights, node_size=[15]*len(G.nodes()))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "weight = [element * 1000 for element in weights]\n",
    "weight = [40 if i>=40 else i for i in weight]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#partition = community_louvain.best_partition(G)\n",
    "pos = nx.spring_layout(G, scale=2)\n",
    "\n",
    "plt.figure(figsize=(10,6))\n",
    "nx.draw(G, pos, node_color=[color_dict[v] for v in target2cluster.values()], edge_color=weight, node_size=[20]*len(G.nodes()))\n",
    "plt.savefig('target_graph.png')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# DDS\n",
    "DDS = pd.read_csv(f'{MAIN_DIR}DDS_known_ChEMBLid_T{num_drug}.csv')\n",
    "DDS.head()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "DDS_new = DDS.copy()\n",
    "DDS_new"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "edgeList_dds = DDS_new.values.tolist()\n",
    "G = networkx.Graph()\n",
    "weights = []\n",
    "for i in range(len(edgeList_dds)):\n",
    "    G.add_edge(edgeList_dds[i][0], edgeList_dds[i][1], weight=edgeList_dds[i][2])\n",
    "    weights.append(edgeList_dds[i][2])\n",
    "    \n",
    "A = networkx.adjacency_matrix(G)\n",
    "DDS_adj = A.copy()\n",
    "\n",
    "louvain = Louvain()\n",
    "labels = louvain.fit_transform(DDS_adj)\n",
    "\n",
    "DDS_cluster_label= pd.DataFrame({'drug':list(G.nodes()), 'label':labels})\n",
    "DDS_cluster_label['label']= DDS_cluster_label['label'].astype(str)\n",
    "labels_unique, counts = np.unique(labels, return_counts=True)\n",
    "print(labels_unique, counts)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "drug2cluster = dict(zip(DDS_cluster_label.drug, DDS_cluster_label.label))\n",
    "color_dict_drug = {k: v for k, v in enumerate( ['mediumvioletred','mediumblue','gold','green','violet','mediumturquoise','mediumvioletred','darkgoldenrod','pink','dimgray'])}\n",
    "color_dict_drug= {str(key): value for key, value in color_dict_drug.items()}\n",
    "\n",
    "with open('color_dict_drug.pkl', 'wb') as fp:\n",
    "    pickle.dump(color_dict_drug, fp)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "modularity2 = community.modularity(drug2cluster, G, weight='weight')\n",
    "print(\"The modularity Q based on networkx is {}\".format(modularity2))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "partition = community_louvain.best_partition(G)\n",
    "pos = community_layout(G, drug2cluster)\n",
    "\n",
    "plt.figure(figsize=(9,5))\n",
    "nx.draw(G, pos, node_color=[color_dict_drug[v] for v in drug2cluster.values()], edge_color=weights, node_size=[15]*len(G.nodes()))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "pos = nx.spring_layout(G, scale=2)\n",
    "\n",
    "plt.figure(figsize=(8,6))\n",
    "nx.draw(G, pos, node_color=[color_dict_drug[v] for v in drug2cluster.values()], edge_color=weights, node_size=[30]*len(G.nodes()))\n",
    "plt.savefig('drug_graph.png')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "DDS_boxplot = DDS[['weight']]\n",
    "DDS_boxplot['type'] = 'DDS'\n",
    "\n",
    "PPS_boxplot = PPS_seq[['weight']]\n",
    "PPS_boxplot['type'] = 'PPS'\n",
    "\n",
    "df_boxplot = DDS_boxplot.append(PPS_boxplot)\n",
    "\n",
    "boxplot = df_boxplot.boxplot(by='type',fontsize=15, figsize=(6,7))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": true
   },
   "outputs": [],
   "source": [
    "vals, names, xs = [],[],[]\n",
    "for i, col in enumerate(DTI[[\"weight\"]]):\n",
    "    vals.append(DTI[col].values)\n",
    "    names.append(col)\n",
    "    xs.append(np.random.normal(i + 1, 0.04, DTI[col].values.shape[0]))\n",
    "\n",
    "plt.boxplot(vals, labels=[\"Drug-target interaction\"])\n",
    "palette = ['y']\n",
    "for x, val, c in zip(xs, vals, palette):\n",
    "    plt.scatter(x, val, alpha=0.2, color='#72BEB7')\n",
    "    \n",
    "    plt.ylabel(\"pChEMBL Value\", fontweight='normal', fontsize=14)\n",
    "    sns.despine(bottom=True) # removes right and top axis lines\n",
    "    plt.axhline(y=5.50, color='b', linestyle='--', linewidth=1, label='Activity thresholds')\n",
    "    plt.axhline(y=0, color='r', linestyle='--', linewidth=1, label='Negative interaction')\n",
    "\n",
    "    plt.legend(bbox_to_anchor=(0.31, 1.06), loc=2, borderaxespad=0, framealpha=1, facecolor ='white', frameon=True)\n",
    "    plt.legend(loc='upper right')\n",
    "    plt.savefig('scatter.png')\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 3- Creating test/train intearctions in ML"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# define postive/negative interactions (positive--> higher than threshold)\n",
    "df_high_pos = DTI[DTI['weight']>= HIGHLY_POS]\n",
    "df_high_pos['L'] = 1\n",
    "\n",
    "df_high_neg = DTI[DTI['weight']< HIGHLY_NEG]\n",
    "df_high_neg['L'] = 0\n",
    "\n",
    "# neutral DTI (positive but less than threshold)\n",
    "df_neutral= DTI[(DTI['weight']<HIGHLY_POS) & (DTI['weight']>HIGHLY_NEG)]\n",
    "df_neutral['L'] = -1\n",
    "\n",
    "print(f'Total number of DTI for creating graph: {len(DTI)}')\n",
    "print(f'Num of gene/protein: {num_protein}')\n",
    "print(f'Num of drug: {num_drug}')\n",
    "print(f'Num of +/- DTI for developing ML-model: {len(df_high_pos)+len(df_high_neg)}')\n",
    "print(f'Num of highly + interactions >= {HIGHLY_POS}: {len(df_high_pos)}')\n",
    "print(f'Num of - interactions = {HIGHLY_NEG}: {len(df_high_neg)}')\n",
    "print(f'Num of positive interactions: {len(df_neutral)}')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# dataframe with real value of interactions as \"weight\" and binary value of interactions as \"label\"\n",
    "DTI_for_ML = df_high_pos.append(df_high_neg)\n",
    "print(DTI_for_ML['L'].value_counts())\n",
    "DTI_for_ML.head()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# without normalization and scaling\n",
    "# Normalizing based on DDS, PPS, DTI_neutral\n",
    "# positive and negative DTI for machine learning are binaray and do not need scaling and normalization\n",
    "\n",
    "DTI_boxplot = DTI[(DTI[\"weight\"]<HIGHLY_POS) & (DTI[\"weight\"]>HIGHLY_NEG)][['weight']]\n",
    "DTI_boxplot['type'] = 'DTI'\n",
    "\n",
    "DDS_boxplot = DDS[['weight']]\n",
    "DDS_boxplot['type'] = 'DDS'\n",
    "\n",
    "PPS_boxplot = PPS_seq[['weight']]\n",
    "PPS_boxplot['type'] = 'PPS'"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#  RUN ALL\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def plot_pca_2d_with_Louvain_clusters(embeddings_seq, drug2cluster, target2cluster, type_data= 'drug'):\n",
    "    \n",
    "    \"\"\"\n",
    "    This function plot PCA of drugs and proteins vectors \n",
    "    \n",
    "    Args: \n",
    "    PPS_seq: A DataFrame of protein-protein smilarities \n",
    "    DDS: A DataFrame of drug-drug smilarities \n",
    "    df_total: A dataframe vectors of drugs and proteins \n",
    "    \n",
    "    \"\"\"\n",
    "    df = embeddings_seq.copy()\n",
    "    df['cluster'] = df.index\n",
    "    \n",
    "    if type_data== 'drug':\n",
    "        df = df[df['cluster'].isin(list(drug2cluster.keys()))]\n",
    "        df['cluster'] = df['cluster'].replace(drug2cluster)\n",
    "    else:\n",
    "        df = df[df['cluster'].isin(list(target2cluster.keys()))]\n",
    "        df['cluster'] = df['cluster'].replace(target2cluster)\n",
    "\n",
    "    \n",
    "    pca = PCA(n_components=2)\n",
    "    #pca = TSNE(n_components=2, random_state=42, perplexity=50, n_iter= 400 )\n",
    "    pca_result = pca.fit_transform(df.drop(columns=['cluster']).values)\n",
    "  \n",
    "    \n",
    "    df['PCA-1'] = pca_result[:,0]\n",
    "    df['PCA-2'] = pca_result[:,1]\n",
    "    \n",
    "    if type_data== 'target':\n",
    "        color_ids = ['#58ACFA','#FF1493', 'yellow','orange', '#00CED1','#5F9EA0','#006400','#96bf65','#fcc808','#7b2b48',\n",
    " '#e96957','#e06000','#173679','#d2dd49','#684a6b','#096eb2','#ce482a', 'red', 'lime', 'lightslategray',\n",
    "                                      'olive', 'rosybrown', 'sienna', 'darkmagenta','midnightblue','maroon',\n",
    "                                      'lightcoral','gold','sandybrown','tomato','lawngreen','lightgreen','darkorchid',\n",
    "                                      'lightskyblue','darkgreen']\n",
    "    if type_data== 'drug':\n",
    "        color_ids = ['mediumvioletred','mediumblue','gold','green','violet','mediumturquoise','mediumvioletred','darkgoldenrod','pink','dimgray']\n",
    "\n",
    "    chart = ggplot(df, aes(x='PCA-1', y='PCA-2',  color='factor(cluster)') ) \\\n",
    "        + geom_point(size=120, alpha=0.8) \\\n",
    "        + scale_color_manual(values= color_ids)    \n",
    "    \n",
    "    chart.save(f'{type_data}_pca_cluster.png', width=12, height=8)  \n",
    "        \n",
    "   \n",
    "    return chart, df"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def plot_pca_2d_with_clusters(df, DDS, PPS_seq, pca_t='both'):\n",
    "    \n",
    "    \"\"\"\n",
    "    This function plot PCA of drugs and proteins vectors \n",
    "    \n",
    "    Args: \n",
    "    PPS_seq: A DataFrame of protein-protein smilarities \n",
    "    DDS: A DataFrame of drug-drug smilarities \n",
    "    df_total: A dataframe vectors of drugs and proteins \n",
    "    \n",
    "    \"\"\"\n",
    "    df_ChEMBL2DrugBank = pd.read_csv('ChEMBL2DrugBank.csv')\n",
    "    ChEMBL2DrugBank = pd.Series(df_ChEMBL2DrugBank['DrugBank'].values, index=df_ChEMBL2DrugBank['ChEMBL']).to_dict()\n",
    "    DrugBank2ChEMBL = dict((y,x) for x,y in ChEMBL2DrugBank.items())\n",
    "\n",
    "    list_protein = ['unclassified_protein', 'Trascription_factor', 'Transporter', 'Secreted_protein',\n",
    "                   'other_categories', 'membrance_receptor', 'Ion_Channel', 'Epigenetic_regulation', \n",
    "                    'Enzyme_Cytochrome P450', 'Cytosolic_protein', 'Enzyme_Hydrolase', 'Enzyme_Kinase',\n",
    "                    'Enzyme_Ligase','Enzyme_Lyase','Enzyme_Oxidoreductase','Enzyme_Phosphatase',\n",
    "                    'Enzyme_Protease', 'Enzyme_rest', 'Enzyme_Transferase']\n",
    "    \n",
    "\n",
    "    protein_type = {}\n",
    "    for p in list_protein:\n",
    "        df_p = pd.read_csv(f'Protein_class/{p}.csv')\n",
    "        df_p['type_protein'] = p\n",
    "\n",
    "        dict_protein = pd.Series(df_p['type_protein'].values, index = df_p['ChEMBL ID']).to_dict()\n",
    "        dict_protein['CHEMBL2055'] = 'Trascription_factor'\n",
    "        protein_type.update(dict_protein)\n",
    "\n",
    "    df_total = df.copy()\n",
    "    df_total['name'] = df_total.index\n",
    "    \n",
    "    drug_name = list(set(list(DDS['from'])+list(DDS['to'])))\n",
    "    protein_name = list(set(list(PPS_seq['from'])+list(PPS_seq['to'])))\n",
    "\n",
    "    drug_vec = df_total.loc[drug_name]\n",
    "    drug_vec['type'] = 'Drug'\n",
    "    drug_vec['name'] = 'CHEMBL' + drug_vec['name'].astype(str)\n",
    "    \n",
    "    with open('protein2num.pkl', 'rb') as f:\n",
    "            protein2num = pickle.load(f)\n",
    "    num2protein = dict((str(int(y)),x) for x,y in protein2num.items())\n",
    "   \n",
    "            \n",
    "    protein_vec = df_total.loc[protein_name]\n",
    "    protein_vec['name'] = protein_vec['name'].astype(str).replace(num2protein)\n",
    "    protein_vec['type']= protein_vec['name'].replace(protein_type)\n",
    "    \n",
    "    \n",
    "    if pca_t=='drug':\n",
    "        df = drug_vec.copy()\n",
    "    elif pca_t== 'target':\n",
    "        df = protein_vec.copy()\n",
    "    \n",
    "    pca = PCA(n_components=2)\n",
    "   # pca = TSNE(n_components=2, random_state=42)\n",
    "    pca_result = pca.fit_transform(df.drop(columns=['type', 'name']).values)\n",
    "  \n",
    "    \n",
    "    df['PCA-1'] = pca_result[:,0]\n",
    "    df['PCA-2'] = pca_result[:,1]\n",
    "        \n",
    "    chart = ggplot(df, aes(x='PCA-1', y='PCA-2',  color='factor(type)', label = 'factor(name)') ) \\\n",
    "        + geom_point(size=100, alpha=0.8) \\\n",
    "        + scale_color_manual(values = ['#58ACFA','#FF1493', '#00BFFF','#00CED1','#5F9EA0','#006400','#96bf65','#fcc808','#7b2b48',\n",
    " '#e96957','#e06000','#173679','#e8a1a2','#d2dd49','#684a6b','#096eb2','#bde1e9','#d2dd49','#ce482a'])\\\n",
    "      #  + geom_text(aes(label='factor(name)'), size=10, color='#464546')\n",
    "    \n",
    "    \n",
    "    chart.save(f'{pca_t}_pca.png', width=12, height=8)  \n",
    "        \n",
    "   \n",
    "    return chart, df"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def plot_pca(df):\n",
    "    \n",
    "    \"\"\"\n",
    "    This function plot PCA of drugs and proteins vectors \n",
    "    \n",
    "    Args: \n",
    "    PPS_seq: A DataFrame of protein-protein smilarities \n",
    "    DDS: A DataFrame of drug-drug smilarities \n",
    "    df_total: A dataframe vectors of drugs and proteins \n",
    "    \n",
    "    \"\"\"\n",
    "    \n",
    "    pca = PCA(n_components=2)\n",
    "    pca_result = pca.fit_transform(df.drop(columns=['label']).values)\n",
    "  \n",
    "    \n",
    "    df['PCA-1'] = pca_result[:,0]\n",
    "    df['PCA-2'] = pca_result[:,1]\n",
    "        \n",
    "    chart = ggplot(df, aes(x='PCA-1', y='PCA-2',  color='factor(label)') ) \\\n",
    "        + geom_point(size=120, alpha=0.8) \\\n",
    "        \n",
    "    \n",
    "   chart.save(f'{pca_t}_pca_GSD.png', width=12, height=8)  \n",
    "   \n",
    "    return chart"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def DT2Vec(DTI, embeddings_seq, HIGHLY_POS, HIGHLY_NEG, n=10, multiclass=False):\n",
    "    \"\"\"\n",
    "    This function ceates the machine learning model (train and test the model) \n",
    "    \n",
    "    Args: \n",
    "        DTI: A DataFrame of drug-target interactions\n",
    "        embeddings_seq: A dataframe vectors of drugs and proteins \n",
    "        HIGHLY_POS: threshold for positive interactions \n",
    "        HIGHLY_NEG: threshold for negative interactions\n",
    "        n: Number of folds\n",
    "        \n",
    "    Returns: \n",
    "        all_train_results: Traning results for n folds\n",
    "        all_test_results: Testing results for n folds\n",
    "    \"\"\"\n",
    "    \n",
    "    df_high_pos = DTI[DTI[\"weight\"]>= HIGHLY_POS]\n",
    "    df_high_pos[\"L\"] = 1\n",
    "    \n",
    "    df_high_neg = DTI[DTI[\"weight\"]< HIGHLY_NEG]\n",
    "    df_high_neg[\"L\"] = 0\n",
    "\n",
    "    df_neutral = DTI[(DTI[\"weight\"]<HIGHLY_POS) & (DTI[\"weight\"]>HIGHLY_NEG)]\n",
    "    df_neutral[\"L\"] = 2\n",
    "    \n",
    "    if multiclass:\n",
    "         DTI_for_ML = (df_high_pos.append(df_high_neg)).append(df_neutral)\n",
    "    else:\n",
    "        DTI_for_ML = df_high_pos.append(df_high_neg)\n",
    "    \n",
    "    \n",
    "    num_crossVal = 0 \n",
    "    all_train_results, all_test_results = [], []\n",
    "    kf = KFold(n_splits=n, random_state=42, shuffle=True)\n",
    "    kf.get_n_splits(DTI_for_ML)\n",
    "    \n",
    "    for train_index, test_index in kf.split(DTI_for_ML):\n",
    "        DTI_for_ML_train, DTI_for_ML_test = DTI_for_ML.iloc[train_index], DTI_for_ML.iloc[test_index]\n",
    "        num_crossVal = num_crossVal + 1 \n",
    "        print(f'\\n\\n\\n  KFold: {num_crossVal}')\n",
    "        print(f'# all DTI: {DTI_for_ML.shape[0]}')\n",
    "        print(f'# DTI train-set: {DTI_for_ML_train.shape[0]}')\n",
    "        print(f'# DTI test-set: {DTI_for_ML_test.shape[0]}')\n",
    "\n",
    "        df_train_seq = concating(DTI_for_ML_train, embeddings_seq).drop(columns=[\"weight\"])\n",
    "        df_train_seq[\"weight\"] = df_train_seq[\"L\"]\n",
    "        final_train_seq = df_train_seq.drop(columns=[\"L\"])\n",
    "\n",
    "        df_test_seq = concating(DTI_for_ML_test , embeddings_seq).drop(columns=[\"weight\"])\n",
    "        df_test_seq [\"weight\"] = df_test_seq [\"L\"]\n",
    "        final_test_seq = df_test_seq.drop(columns=[\"L\"])\n",
    "        final_test_seq.head()\n",
    "\n",
    "        X, y = get_data(final_train_seq)\n",
    "\n",
    "        if multiclass:\n",
    "            model = XGBClassifier(objective='multi:softmax', learning_rate=0.9, max_depth=2, min_child_weight=2)\n",
    "        else:\n",
    "            model = XGBClassifier(learning_rate= 0.4, max_depth= 4, min_child_weight=2) \n",
    "\n",
    "        model.fit(X, y)\n",
    "\n",
    "        pred_train = model.predict(X)\n",
    "        \n",
    "        if multiclass:\n",
    "            train_results = calculate_metrics_multiclass(y, pred_train)\n",
    "        else:\n",
    "            train_results = calculate_metrics(y, pred_train)\n",
    "            \n",
    "        train_results['model'] = model\n",
    "        train_results['embbeding_seq'] = embeddings_seq\n",
    "        train_results['Kfold'] = num_crossVal\n",
    "        all_train_results.append(train_results)\n",
    "        del X\n",
    "        del y\n",
    "\n",
    "        print('\\n')\n",
    "\n",
    "        X, y = get_data(final_test_seq)\n",
    "        pred_test = model.predict(X)\n",
    "        test_results = calculate_metrics(y, pred_test)\n",
    "        test_results['model'] = model\n",
    "        test_results['Kfold'] = num_crossVal\n",
    "        all_test_results.append(test_results)\n",
    "        del X\n",
    "        del y\n",
    "            \n",
    "    return pd.DataFrame(all_train_results), pd.DataFrame(all_test_results), final_train_seq"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def DT2Vec_external(DTI_new, DTI_external_test, embeddings_seq, HIGHLY_POS, HIGHLY_NEG, multiclass=False):\n",
    "    \"\"\"\n",
    "    This function ceates the machine learning model  \n",
    "    \n",
    "    Args: \n",
    "        DTI: A DataFrame of drug-target interactions\n",
    "        embeddings_seq: A dataframe vectors of drugs and proteins \n",
    "        HIGHLY_POS: threshold for positive interactions \n",
    "        HIGHLY_NEG: threshold for negative interactions\n",
    "        \n",
    "    Returns: \n",
    "        all_train_results: Traning results for n folds\n",
    "        all_test_results: Testing results for n folds\n",
    "        df_DTI_proba: A dataframe of DTI with proba\n",
    "\n",
    "    \"\"\"\n",
    "    \n",
    "    df_high_pos = DTI_new[DTI_new[\"weight\"]>= HIGHLY_POS]\n",
    "    df_high_pos[\"L\"] = 1\n",
    "    \n",
    "    df_high_neg = DTI_new[DTI_new[\"weight\"]<= HIGHLY_NEG]\n",
    "    df_high_neg[\"L\"] = 0\n",
    "\n",
    "    df_neutral = DTI_new[(DTI_new[\"weight\"]<HIGHLY_POS) & (DTI_new[\"weight\"]>HIGHLY_NEG)]\n",
    "    df_neutral[\"L\"] = 2\n",
    "    \n",
    "    if multiclass:\n",
    "        DTI_for_ML = (df_high_pos.append(df_high_neg)).append(df_neutral)\n",
    "    else:\n",
    "        DTI_for_ML = df_high_pos.append(df_high_neg)\n",
    "\n",
    "    df_seq = concating(DTI_for_ML, embeddings_seq).drop(columns=[\"weight\"])\n",
    "    df_seq[\"weight\"] = df_seq[\"L\"]\n",
    "    final_seq = df_seq.drop(columns=[\"L\"])\n",
    "        \n",
    "    X, y = get_data(final_seq)\n",
    "\n",
    "\n",
    "    if multiclass:\n",
    "        model = XGBClassifier(objective='multi:softmax', learning_rate=0.9, max_depth=2, min_child_weight=2)\n",
    "    else:\n",
    "        model = XGBClassifier(learning_rate=0.9, max_depth=2, min_child_weight=2)\n",
    "\n",
    "    model.fit(X, y)\n",
    "    \n",
    "    ## for the test-set ##\n",
    "    df_high_pos_test = DTI_external_test[DTI_external_test[\"weight\"]>= HIGHLY_POS]\n",
    "    df_high_pos_test[\"L\"] = 1\n",
    "    \n",
    "    df_high_neg_test = DTI_external_test[DTI_external_test[\"weight\"]<=HIGHLY_NEG]\n",
    "    df_high_neg_test[\"L\"] = 0\n",
    "\n",
    "    df_neutral = DTI_external_test[(DTI_external_test[\"weight\"]<HIGHLY_POS) & (DTI_external_test[\"weight\"]>HIGHLY_NEG)]\n",
    "    df_neutral[\"L\"] = 2\n",
    "    \n",
    "    if multiclass:\n",
    "         DTI_for_ML_test = (df_high_pos_test.append(df_high_neg_test)).append(df_neutral_test)\n",
    "    else:\n",
    "        DTI_for_ML_test = df_high_pos_test.append(df_high_neg_test)\n",
    "\n",
    "    df_seq_test = concating(DTI_for_ML_test, embeddings_seq).drop(columns=[\"weight\"])\n",
    "    df_seq_test[\"weight\"] = df_seq_test[\"L\"]\n",
    "    final_seq_test = df_seq_test.drop(columns=[\"L\"])\n",
    "    \n",
    "    X_test, y_test = get_data(final_seq_test)\n",
    "    pred_test = model.predict(X_test)\n",
    "    test_results = calculate_metrics(y_test, pred_test)\n",
    "    del X\n",
    "    del y\n",
    "    \n",
    "    \n",
    "    return model, test_results, final_seq_test"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def DT2Vec_external_best_model(model, DTI_external_test, embeddings_seq, HIGHLY_POS, HIGHLY_NEG, multiclass=False):\n",
    "    \"\"\"\n",
    "    This function ceates the machine learning model  \n",
    "    \n",
    "    Args: \n",
    "        DTI: A DataFrame of drug-target interactions\n",
    "        embeddings_seq: A dataframe vectors of drugs and proteins \n",
    "        HIGHLY_POS: threshold for positive interactions \n",
    "        HIGHLY_NEG: threshold for negative interactions\n",
    "        \n",
    "    Returns: \n",
    "        all_train_results: Traning results for n folds\n",
    "        all_test_results: Testing results for n folds\n",
    "        df_DTI_proba: A dataframe of DTI with proba\n",
    "\n",
    "    \"\"\"\n",
    "    \n",
    "    ## for the test-set ##\n",
    "    df_high_pos_test = DTI_external_test[DTI_external_test[\"weight\"]>= HIGHLY_POS]\n",
    "    df_high_pos_test[\"L\"] = 1\n",
    "    \n",
    "    df_high_neg_test = DTI_external_test[DTI_external_test[\"weight\"]<= HIGHLY_NEG]\n",
    "    df_high_neg_test[\"L\"] = 0\n",
    "\n",
    "    df_neutral = DTI_external_test[(DTI_external_test[\"weight\"]<HIGHLY_POS) & (DTI_external_test[\"weight\"]>HIGHLY_NEG)]\n",
    "    df_neutral[\"L\"] = 2\n",
    "    \n",
    "    if multiclass:\n",
    "         DTI_for_ML_test = (df_high_pos_test.append(df_high_neg_test)).append(df_neutral_test)\n",
    "    else:\n",
    "        DTI_for_ML_test = df_high_pos_test.append(df_high_neg_test)\n",
    "\n",
    "    df_seq_test = concating(DTI_for_ML_test, embeddings_seq).drop(columns=[\"weight\"])\n",
    "    df_seq_test[\"weight\"] = df_seq_test[\"L\"]\n",
    "    final_seq_test = df_seq_test.drop(columns=[\"L\"])\n",
    "    \n",
    "    X_test, y_test = get_data(final_seq_test)\n",
    "    pred_test = model.predict(X_test)\n",
    "    test_results = calculate_metrics(y_test, pred_test)\n",
    "    \n",
    "    return test_results, final_seq_test"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##  DT2Vec "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Embed drugs and proteins to vector \n",
    "embeddings_seq, embeddings_seq_drug, embeddings_seq_protein = embedding(DDS, PPS_seq, mapping=False)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": false
   },
   "outputs": [],
   "source": [
    "chart, df_target_seq = plot_pca_2d_with_clusters(embeddings_seq, DDS, PPS_seq, pca_t='target')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": false
   },
   "outputs": [],
   "source": [
    "chart, df_Lou_target = plot_pca_2d_with_Louvain_clusters(embeddings_seq, drug2cluster, target2cluster, type_data= 'target')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "chart, df_Lou_drug = plot_pca_2d_with_Louvain_clusters(embeddings_seq, drug2cluster, target2cluster, type_data= 'drug')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## K-means Clustering"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "cluster_value = {}\n",
    "from sklearn.metrics import precision_score, recall_score, accuracy_score, f1_score, confusion_matrix, silhouette_score, silhouette_samples, calinski_harabasz_score, davies_bouldin_score, roc_curve, auc\n",
    "\n",
    "\n",
    "for n_clusters in range(4,5):\n",
    "    clusterer = MiniBatchKMeans(n_clusters, random_state=42)\n",
    "    cluster_labels = clusterer.fit_predict(embeddings_seq_drug)\n",
    "    \n",
    "\n",
    "    # The silhouette_score gives the average value for all the samples.\n",
    "    # This gives a perspective into the density and separation of the formed\n",
    "    # clusters\n",
    "    silhouette_avg = silhouette_score(embeddings_seq_drug, cluster_labels, random_state=42)\n",
    "    print(\"For n_clusters =\", n_clusters,\n",
    "          \"The average silhouette_score is :\", silhouette_avg)\n",
    "\n",
    "    calinski_harabaz_avg = calinski_harabasz_score(embeddings_seq_drug, cluster_labels)\n",
    "    print(\"The average calinski_harabaz is :\", calinski_harabaz_avg)\n",
    "\n",
    "    davies_bouldin_avg= davies_bouldin_score(embeddings_seq_drug, cluster_labels)\n",
    "    print(\"The average davies_bouldin_score is :\", davies_bouldin_avg)\n",
    "    \n",
    "    cluster_value[n_clusters] = {'S':silhouette_avg, 'CH':calinski_harabaz_avg, 'DB':davies_bouldin_avg, 'labels':cluster_labels}\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "len(cluster_value[4]['labels'])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "df_drug_kmeans = embeddings_seq_drug.copy()\n",
    "df_drug_kmeans['label'] = cluster_value[4]['labels']\n",
    "df_drug_kmeans"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "chart, df_drug_seq = plot_pca_2d_with_clusters(embeddings_seq, DDS, PPS_seq, pca_t='drug')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Adding clusters to DTI"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "protein_dict = pd.Series(df_target_seq['type'].values, index = df_target_seq.index).to_dict()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "dict_Lou_target = df_Lou_target['cluster'].to_dict()\n",
    "dict_Lou_drug = df_Lou_drug['cluster'].to_dict()\n",
    "\n",
    "DTI_cluster = DTI.copy()\n",
    "DTI_cluster['target_C']= DTI_cluster['to'].map(dict_Lou_target) \n",
    "DTI_cluster['drug_C']= DTI_cluster['from'].map(dict_Lou_drug)\n",
    "DTI_cluster['target_type']= DTI_cluster['to'].map(protein_dict) \n",
    "DTI_cluster['label']= DTI_cluster['weight']\n",
    "DTI_cluster.loc[DTI_cluster['label'] <= HIGHLY_NEG, 'label'] = -100\n",
    "DTI_cluster.loc[DTI_cluster['label'] >= HIGHLY_POS, 'label'] = 100\n",
    "DTI_cluster['label'] = np.where((DTI_cluster['label']<HIGHLY_POS )&(DTI_cluster['label']>HIGHLY_NEG), -1 , DTI_cluster['label'])\n",
    "DTI_cluster['label'] = DTI_cluster['label'].map({-100:'Negative',100:'Positive', -1:'Weak'})\n",
    "\n",
    "DTI_cluster['from_ChEMBL'] = 'CHEMBL' + DTI_cluster['from'].astype(str)\n",
    "DTI_cluster['to_ChEMBL'] = DTI_cluster['to'].astype(str).map(num2protein)\n",
    "\n",
    "\n",
    "DTI_cluster.to_csv('DTI_cluster.csv')\n",
    "\n",
    "with open('Lou_target.pkl', 'wb') as fp:\n",
    "    pickle.dump(dict_Lou_target, fp)\n",
    "\n",
    "with open('Lou_drug.pkl', 'wb') as fp:\n",
    "    pickle.dump(dict_Lou_drug, fp)\n",
    "    \n",
    "with open('protein_dict.pkl', 'wb') as fp:\n",
    "    pickle.dump(protein_dict, fp)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# split data to train-validation and external test \n",
    "DTI_external_test = DTI.sample(frac = 0.1)\n",
    "DTI_new = DTI.drop(DTI_external_test.index)\n",
    "print (f'DTI: {len(DTI)}')\n",
    "print (f'DTI_external_test: {len(DTI_external_test)}')\n",
    "print (f'DTI_train_internal_test: {len(DTI_new)}')\n",
    "DTI_external_test.head()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": true
   },
   "outputs": [],
   "source": [
    "# XGBoost\n",
    "resuts_all_runs_list = []\n",
    "resuts_external_list = []\n",
    "\n",
    "for i in range(0,1):\n",
    "    DTI_external_test = DTI.sample(frac = 0.1)\n",
    "    DTI_new = DTI.drop(DTI_external_test.index)\n",
    "    \n",
    "    all_train_results, all_test_results, final_train_seq = DT2Vec(DTI_new, embeddings_seq, HIGHLY_POS, HIGHLY_NEG, n=10, multiclass=False)\n",
    "    resuts_all_runs_list.append(all_test_results) \n",
    "    # select the best model based on \"average_precision_score\"\n",
    "    model = all_test_results.sort_values(by=['average_precision_score'], ascending=False)['model'][0]\n",
    "    test_results, final_seq_test = DT2Vec_external_best_model(model, DTI_external_test, embeddings_seq, HIGHLY_POS, HIGHLY_NEG, multiclass=False)\n",
    "    resuts_external_list.append(test_results) \n",
    "\n",
    "resuts_all_runs = pd.concat(resuts_all_runs_list, ignore_index=True)\n",
    "resuts_external = pd.DataFrame(resuts_external_list)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# results on external-testset (5 runs)--> best model (of 10 fold CV) was selected and was applied on external testsets in each run\n",
    "resuts_external"
   ]
  },
  {
   "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.8.10"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}