<jupyter_start><jupyter_text>Drug-Drug Similiarity SystemBased on the best model among the 6 graph embedding algorithms, we are going to use the embeddings, link predictions and the drug clustering to create a combined system to find the similiarities among the drug. This system will also help us find a similiar drug for any given input drug.   Setup Notebook and GPU<jupyter_code># Monunt google drive to import KEGG dataset
from google.colab import drive
drive.mount('/content/drive')
# Setting Colab's GPU
# Edit > Notebook Settings > GPU 
!nvidia-smi -L
# Detailed GPU Information Reference: https://medium.com/analytics-vidhya/the-google-colab-system-specification-check-69d159597417
# !nvidia-smi -q
!lscpu |grep 'Model name'<jupyter_output>Model name:          Intel(R) Xeon(R) CPU @ 2.30GHz<jupyter_text>Import libraries<jupyter_code>%%capture 
!pip install ampligraph;
%tensorflow_version 1.x 
import numpy as np
import pandas as pd
import ampligraph
from ampligraph.datasets import load_from_csv
from tabulate import tabulate
ampligraph.__version__<jupyter_output><empty_output><jupyter_text>Import Dataset<jupyter_code># Load the KEGG dataset
df = load_from_csv('.','drive/My Drive/datasets/KEGG/KEGG50K.txt',sep='\t')
headers = ["Head", "Relationship", "Tail"]
table = tabulate(df[:5], headers, tablefmt="fancy_grid") # tabulate data Reference: https://stackoverflow.com/questions/9712085/numpy-pretty-print-tabular-data
print(table) # output
print('Full set size: ', df.shape)
positives_filter = df<jupyter_output><empty_output><jupyter_text>Graph Embedding Algorithm The best performance was given by the ComplEx algorithms providing the best scores on link prediction tasks. The configurations are as follows:* **Embedding size**: 300* **Optimizer**: adam* **Learning rate**: 0.001* **Loss**: Multiclass NLL Loss* **Regularization type**: 3* **Regularization constant**: 0.01<jupyter_code>"""
Model configs ( Embedding size :  300 , Optimizer :  adam , with learning rate :  0.001 , Loss :  multiclass_nll ,regularization type:  3  regularization constant lambda:  0.01  )
Execution time: 136.5050220489502 seconds.
Number :  31
MRR: 0.46
Hits@10: 0.64
Hits@3: 0.49
Hits@1: 0.37
"""
from ampligraph.latent_features import ComplEx

model = ComplEx(batches_count=100, 
                seed=0, 
                epochs=200, 
                k=300, 
                eta=5,
                optimizer='adam', 
                optimizer_params={'lr':1e-3},
                loss='multiclass_nll',                     # Multiclass NLL Loss Reference: https://docs.ampligraph.org/en/latest/generated/ampligraph.latent_features.NLLMulticlass.html#ampligraph.latent_features.NLLMulticlass
                regularizer='LP', 
                regularizer_params={'p':3,'lambda':1e-2},  # The nuclear 3-norm proposed in the ComplEx-N3 from paper LUO18, proves to be effective for DistMult and ComplEx Reference: https://github.com/DeepGraphLearning/KnowledgeGraphEmbedding/issues/6 
                verbose=True)<jupyter_output><empty_output><jupyter_text>Fitting model<jupyter_code>import tensorflow as tf
import time
tf.logging.set_verbosity(tf.logging.ERROR)
start_time = time.time()
model.fit(df, early_stopping = False)
print("Execution time: %s seconds." % (time.time() - start_time))<jupyter_output>Average ComplEx Loss:   0.166244: 100%|██████████| 200/200 [02:25<00:00,  1.37epoch/s]<jupyter_text>Saving model<jupyter_code>from ampligraph.latent_features import save_model, restore_model
save_model(model, './ComplEx_Model.pkl')<jupyter_output><empty_output><jupyter_text>Generating Embeddings for Entities<jupyter_code>entities = np.unique(np.concatenate([df[:, 0], df[:, 2]]))
relations = np.unique(df[:, 1])
entity_embeddings = dict(zip(entities, model.get_embeddings(entities)))
len(entity_embeddings['D00043'])<jupyter_output><empty_output><jupyter_text>Cosine Similarity among Graphical Embeddings We will first find the the top k similiar drugs based on cosine similiarities wrt other drugs.<jupyter_code>def entityType(n):         # Utility functon to return entity type
  if n.startswith("h"):
    return "PATHWAY"
  elif n.startswith("D"):
    return "DRUG"
  elif n.startswith("N"):
    return "NETWORK"
  elif n.startswith("HS"):
    return "GENE"
  else:
    return "DISEASE"
entityTypeList = list(map(entityType,entities))
print(pd.DataFrame({"Entity": entities[:6000],"Type":entityTypeList[:6000]}))
drugList = []
geneList,diseaseList,pathwayList=[],[],[]
for i,entity in enumerate(entityTypeList):
  if entityTypeList[i]=='DRUG':
    drugList.append(entities[i])
  elif entityTypeList[i]=='GENE':
    geneList.append(entities[i])
  elif entityTypeList[i]=='DISEASE':
    diseaseList.append(entities[i])
  elif entityTypeList[i]=='PATHWAY':
    pathwayList.append(entities[i])
print("Total number of unique drugs : ",len(drugList))
#print("Name of Drug :",drugList[0]," Embedding :",entity_embeddings[drugList[0]])
from scipy import spatial
def cosineSimilarity(emd1,emd2):
  return 1 - spatial.distance.cosine(emd1,emd2)
cosineSimilarity(entity_embeddings[drugList[0]],entity_embeddings[drugList[0]])
def findSimiliarDrugs(drug,drugList,k):     # Utility function to find top k similiar drugs wrt the input drug
  similarityMap = dict()
  drugEmbedding = entity_embeddings[drug]
  for ent in drugList:
    if ent != drug:
      similarityMap[ent] = cosineSimilarity(drugEmbedding,entity_embeddings[ent])
  return sorted(similarityMap.items(), key=lambda kv: -kv[1])[:k] # return a sorted dictionary 
k = findSimiliarDrugs('D00043',drugList,10)
drugs,cosines=[],[]
for u,v in k:
  drugs.append(u)
  cosines.append(v)
# Top 10 drugs similiar based on cosine similairity between graph embeddings
pd.DataFrame({'drug':drugs,'cosine':cosines})<jupyter_output><empty_output><jupyter_text>* **D00043**: Name is IsoflurophateFluostigmine, name(TN) is Floropryl (TN). Also is part of two drug groups as per KEGG database DG01593 (Acetylcholinesterase inhibitor) Class: Neuropsychiatric agent and DG01594 (Butyrylcholinesterase inhibitor)	Class: Neuropsychiatric agent.* **D02558**: Name is Rivastigmine tartrateRivastigmine hydrogen tartrate used for Alzheimer's disease (DS:H00056). Also is a part of drug group DG00984 (Rivastigmine).  Clustering on Graphical Embeddings of Drugs<jupyter_code>%%capture
!pip install seaborn adjustText
from sklearn.decomposition import PCA
import matplotlib.pyplot as plt
import seaborn as sns
from adjustText import adjust_text
%matplotlib inline
drug_embeddings = dict(zip(drugList, model.get_embeddings(drugList)))
embeddings_2d = PCA(n_components=2).fit_transform(np.array([i for i in drug_embeddings.values()]))
embeddings_2d[:2]
def plot_clusters(hue,df,saved=False,figname='default.png'):
  np.random.seed(0)
  plt.figure(figsize=(12, 12))
  plt.title("{} embeddings".format(hue).capitalize())
  ax = sns.scatterplot(data=df, x="embedding1", y="embedding2", hue=hue)
  if saved==True:
    plt.savefig(figname,bbox_inches = 'tight')
from ampligraph.discovery import find_clusters
from sklearn.cluster import KMeans
# https://github.com/prakhargurawa/Machine-Learning-A-Z/blob/main/4_Clustering/k_means_clustering.py
# https://scikit-learn.org/stable/modules/generated/sklearn.cluster.KMeans.html
WCSS = []
for i in range(1, 11):
    km = KMeans(n_clusters=i, n_init=50, max_iter=500, random_state=0,init = 'k-means++')
    km.fit(embeddings_2d )
    WCSS.append(km.inertia_)

# plot
plt.plot(range(1, 11), WCSS, marker='o')
plt.xlabel('Number of clusters')
plt.ylabel('Distortion')
plt.show()
# 3 is perfect for now
clustering_algorithm = KMeans(n_clusters=3, n_init=50, max_iter=500, random_state=0,init = 'k-means++')
clusters = find_clusters(np.array(drugList), model, clustering_algorithm, mode='entity')
plot_df = pd.DataFrame({"Drug": drugList, 
                        "embedding1": embeddings_2d[:, 0], 
                        "embedding2": embeddings_2d[:, 1],
                        "cluster": "cluster" + pd.Series(clusters).astype(str)})
plot_df
plot_clusters("cluster",plot_df,saved=True,figname="drug-cluster.png")
def entityType(n):         # Utility functon to return entity type
  if n.startswith("h"):
    return "PATHWAY"
  elif n.startswith("D"):
    return "DRUG"
  elif n.startswith("N"):
    return "NETWORK"
  elif n.startswith("HS"):
    return "GENE"
  else:
    return "DISEASE"
entityTypeList = list(map(entityType,drugList))
entity_df = pd.DataFrame({"Drug": drugList, 
                        "embedding1": embeddings_2d[:, 0], 
                        "embedding2": embeddings_2d[:, 1],
                        "type":entityTypeList})
entity_df
entity_df['type'] = np.where((entity_df.Drug == 'D00043'),'D00043',entity_df.type)     # Reference : https://kanoki.org/2019/07/17/pandas-how-to-replace-values-based-on-conditions/
for drug in drugs:
  entity_df['type'] = np.where((entity_df.Drug == drug),'Top similair drugs',entity_df.type)
entity_df
alphas = []
for drug in entity_df['Drug']:
  if drug == 'D00043' or drug in drugs:
    alphas.append(1.5)
  else:
    alphas.append(0.01)
def plot_clusters(hue,alpha,df,saved=False,figname='default.png'):
  np.random.seed(0)
  plt.figure(figsize=(12, 12))
  plt.title("2D Embeddings")
  ax = sns.scatterplot(data=df, x="embedding1", y="embedding2", hue=hue,size=alpha,legend=False)
  if saved==True:
    plt.savefig(figname,bbox_inches = 'tight')

plot_clusters("type",alphas,entity_df,saved=True,figname="drug-plot.png")<jupyter_output><empty_output><jupyter_text>Link Prediction ErrorsDrugs only have relationships with Genes, Pathways and Diseases. We will try to predict there links with other possible entities and try to find error between the link prediction between pair of drugs wrt to there scores with other entities. Similair drug should produce nearby link prediction scores wrt to other entities.<jupyter_code>geneList,diseaseList,pathwayList = [],[],[]
for i,entity in enumerate(entityTypeList):
  if entityTypeList[i]=='GENE':
    geneList.append(entities[i])
  elif entityTypeList[i]=='DISEASE':
    diseaseList.append(entities[i])
  elif entityTypeList[i]=='PATHWAY':
    pathwayList.append(entities[i])
def generateData(drug,geneList,diseaseList,pathwayList):
  x = [[drug,'DRUG_TARGET_GENE',gene] for gene in geneList] + [[drug,'DRUG_EFFICACY_DISEASE',disease] for disease in diseaseList] + [[drug,'DRUG_TARGET_PATHWAY',pathway] for pathway in pathwayList] 
  return np.array(x)
from ampligraph.evaluation import evaluate_performance
from scipy.special import expit

def findLinkPredictionDifference(drug1,drug2,geneList,diseaseList,pathwayList,positives_filter,model):
  X_unseen_drug1 = generateData(drug1,geneList,diseaseList,pathwayList)
  X_unseen_drug2 = generateData(drug2,geneList,diseaseList,pathwayList)

  # calculate scores
  scores_drug1 = model.predict(X_unseen_drug1)
  scores_drug2 = model.predict(X_unseen_drug2)
  # calculate probablities
  probs_drug1 = expit(scores_drug1)
  probs_drug2 = expit(scores_drug2)
  # calculate errors
  mses = ((np.array(probs_drug1)-np.array(probs_drug2))**2).mean()
  return mses

error = findLinkPredictionDifference('D00043','D02558',geneList,diseaseList,pathwayList,positives_filter,model)
error
drugs,cosines,errors,ratio=[],[],[],[]
for u,v in k:
  drugs.append(u)
  cosines.append(v)
  error = findLinkPredictionDifference('D00043',u,geneList,diseaseList,pathwayList,positives_filter,model)
  errors.append(error)
  ratio.append(float(v/error))
pd.DataFrame({'Drug':drugs,'Cosine':cosines,'Error':errors,'Ratio':ratio})<jupyter_output><empty_output><jupyter_text>* **D00196**: Name Physostigmine (USP), part of drug group DG01132,	Physostigmine. Formula: C15H21N3O2.  Calculating top 10 drugs based on link prediction error<jupyter_code>drugs,cosines,errors,ratio=[],[],[],[]
drugEmbedding = entity_embeddings['D00043']
i = 0
for drug in drugList:
  i = i+1
  print(i)
  drugs.append(drug)
  if drug != 'D00043':
    v = cosineSimilarity(drugEmbedding,entity_embeddings[drug])
    cosines.append(v)
    error = findLinkPredictionDifference('D00043',drug,geneList,diseaseList,pathwayList,positives_filter,model)
    errors.append(error)
    ratio.append(float(v/error))
top10linkpredictionscore = pd.DataFrame({'Drug':drugList[1:],'Cosine':cosines,'Error':errors,'Ratio':ratio})
top10linkpredictionscore
cosine_df = top10linkpredictionscore.sort_values(by=['Cosine'], ascending=False)
cosine_df[:10]
error_df = top10linkpredictionscore.sort_values(by=['Error'], ascending=True)
error_df[:10]
ratio_df = top10linkpredictionscore.sort_values(by=['Ratio'], ascending=False)
ratio_df[:10]<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.12709851124485272, 'fraction_numerical': 0.028508077288565095, 'mean_word_length': 5.62591815320042, 'pattern_counts': {'":': 10, '<': 44, '<?xml version=': 0, '>': 45, 'https://': 7, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '11', 'dir': 'jupyter-structured-clean-dedup', 'id': '233780', 'license': 'MIT', 'n_tokens_mistral': 4675, 'n_tokens_neox': 4271, 'n_words': 907, 'path': 'Phase3_Drug_Drug_Similarity.ipynb', 'repo_name': 'prakhargurawa/Drug-Similarity-and-Link-Prediction-using-Graph-Embeddings-on-Medical-Knowledge-Graph', 'size': '12628'}