PoincareMSA_colab_examples.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {
    "colab_type": "text",
    "id": "view-in-github"
   },
   "source": [
    "<a href=\"https://colab.research.google.com/github/DSIMB/PoincareMSA/blob/master/PoincareMSA_colab_examples.ipynb\" target=\"_parent\"><img src=\"https://colab.research.google.com/assets/colab-badge.svg\" alt=\"Open In Colab\"/></a>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "W5FTaqy0xqw2"
   },
   "source": [
    "<img src=\"https://github.com/DSIMB/PoincareMSA/blob/master/.github/PoincareMSA_small_logo.png?raw=true\" height=\"100\" style=\"height:100px;margin-left: 0px;\">\n",
    "\n",
    "# Poincaré maps for visualization of large protein famillies\n",
    "\n",
    "**Authors**: Anna Klimovskaia Susmelj <a href=\"https://orcid.org/0000-0002-1057-6690\" target=\"_blank\"><img alt=\"ORCID logo\" src=\"https://info.orcid.org/wp-content/uploads/2019/11/orcid_16x16.png\"></a>, Yani Ren, Yann Vander Meersche <a href=\"https://orcid.org/0000-0002-6680-8104\" target=\"_blank\"><img alt=\"ORCID logo\" src=\"https://info.orcid.org/wp-content/uploads/2019/11/orcid_16x16.png\"></a>, Jean-Christophe Gelly <a href=\"https://orcid.org/0000-0001-5138-361X\" target=\"_blank\"><img alt=\"ORCID logo\" src=\"https://info.orcid.org/wp-content/uploads/2019/11/orcid_16x16.png\"></a> and Tatiana Galochkina <a href=\"https://orcid.org/0000-0002-3608-5208\" target=\"_blank\"><img alt=\"ORCID logo\" src=\"https://info.orcid.org/wp-content/uploads/2019/11/orcid_16x16.png\"></a>\n",
    "\n",
    "<p>\n",
    "\tPoincaréMSA is a computational tool for visualization of large protein families developed by <a href=\"https://doi.org/10.1093/bib/bbad103\" target=\"_blank\" data-bs-toggle=\"tooltip\" style=\"text-decoration:none;\" data-bs-placement=\"top\" title=\"Susmelj A.K., Ren Y., Vander Meersche Y., Gelly J.C, Galochkina T. (2023). Poincaré maps for visualization of large protein families. Briefings in bioinformatics, -. https://doi.org/10.1093/bib/bbad103\">Susmelj et al.</a> PoincaréMSA builds an interactive projection of an input protein multiple sequence alignment provided by the user or built on-the-flow from the target sequence. The underlying algorithm is based on Poincaré maps introduced by <a href=\"https://doi.org/10.1038/s41467-020-16822-4\" target=\"_blank\" data-bs-toggle=\"tooltip\" style=\"text-decoration:none;\" data-bs-placement=\"top\" title=\"Klimovskaia, A., Lopez-Paz, D., Bottou, L. et al. Poincaré maps for analyzing complex hierarchies in single-cell data. Nat Commun 11, 2966 (2020).\">Klimovskaia et al.</a>. It successfully reproduces both local proximities between protein sequences as well as global hierarchical structure of the data.\n",
    "</p>\n",
    "<p>\n",
    "\tYou can find a step-by-step explication of PoincaréMSA construction on the <a href=\"https://www.dsimb.inserm.fr/POINCARE_MSA/tutorial.html\" target=\"_blank\">Tutorial</a> page and several examples of PoincaréMSA projection for protein families on the <a href=\"https://www.dsimb.inserm.fr/POINCARE_MSA/examples.html\" target=\"_blank\">Examples</a> page. The source code is available at <a href=\"https://github.com/DSIMB/PoincareMSA\" target=\"_blank\">https://github.com/DSIMB/PoincareMSA</a>.<br>\n",
    "</p>\n",
    "<p class=\"mb-1\"><strong>Please use the following reference when citing the PoincaréMSA:</strong></p>\n",
    "<p>Susmelj A.K., Ren Y., Vander Meersche Y., Gelly J.C., &amp; Galochkina T. (2023). Poincaré maps for visualization of large protein families. <i>Briefings in bioinformatics</i>, bbad103, <a href=\"https://doi.org/10.1093/bib/bbad103\" target=\"_blank\">https://doi.org/10.1093/bib/bbad103</a></p>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "Oa7Rrq-YzfQB"
   },
   "source": [
    "# Notebook initialization"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "cellView": "form",
    "id": "X_SO37g-kQq7"
   },
   "outputs": [],
   "source": [
    "#@title ### Load PoincaréMSA Github repository & install dependencies\n",
    "print(\"1. Load PoincaréMSA Github repository\")\n",
    "import os\n",
    "if os.getcwd() == \"/content\":\n",
    "    !git clone https://github.com/DSIMB/PoincareMSA.git\n",
    "    %cd PoincareMSA\n",
    "    !git lfs pull\n",
    "\n",
    "# Check if the GPU is activated\n",
    "import torch\n",
    "device = torch.device('cuda' if torch.cuda.is_available() else 'cpu')\n",
    "\n",
    "print('\\nUsing device:', device)\n",
    "\n",
    "#Install missing module\n",
    "print(\"\\n2. Install dependencies\")\n",
    "!pip install adjustText\n",
    "!pip install -U kaleido\n",
    "\n",
    "import numpy as np\n",
    "import pandas as pd\n",
    "import warnings\n",
    "warnings.filterwarnings('ignore')\n",
    "\n",
    "#File import\n",
    "from google.colab import files\n",
    "import io\n",
    "\n",
    "#Import visualization functions\n",
    "from scripts.visualize_projection.pplots_new import read_embeddings, plot_embedding, plot_embedding_interactive, rotate, get_colors\n",
    "from scripts.build_poincare_map.poincare_maps_rotation import *\n",
    "%matplotlib inline"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "cellView": "form",
    "id": "GcbgnooHwIGh"
   },
   "outputs": [],
   "source": [
    "#@title # Data selection\n",
    "dataset = \"thioredoxins\" #@param [\"kinases\", \"globins\", \"thioredoxins\"]\n",
    "\n",
    "#Select annotations\n",
    "if dataset == \"kinases\":\n",
    "    mfasta = \"examples/kinases/kinases.mfasta\"\n",
    "    path_annotation = \"examples/kinases/kinase_group_new.csv\"\n",
    "elif dataset == \"globins\":\n",
    "    mfasta = \"examples/globins/glob.mfasta\"\n",
    "    path_annotation = \"examples/globins/globin_colors_new.csv\"\n",
    "elif dataset == \"thioredoxins\":\n",
    "    mfasta = \"examples/thioredoxins/thioredoxins.mfasta\"\n",
    "    path_annotation = \"examples/thioredoxins/thioredoxin_annotation.csv\"\n",
    "\n",
    "#Read mfasta file\n",
    "nb_seq = 0\n",
    "with open(mfasta, \"r\") as f_in:\n",
    "    for line in f_in:\n",
    "        if line[0] == \">\":\n",
    "            nb_seq += 1\n",
    "print(f\"\\nNumber of sequences found: {nb_seq}.\")\n",
    "\n",
    "df_annotation = pd.read_csv(path_annotation)\n",
    "print(\"\\nAnnotation file correctly loaded.\")\n",
    "annotation_names = [\"proteins_id\"] + list(df_annotation.loc[:, df_annotation.columns != 'proteins_id'].columns)\n",
    "print(f\"{len(annotation_names)} annotations found: {annotation_names}.\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "_H2vnIt_3kUz"
   },
   "source": [
    "# Data preparation\n",
    "Here we clean the input .mfasta alignment and translate each sequence to a vector ready for projection."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "cellView": "form",
    "id": "9hLq_MgX6hFC"
   },
   "outputs": [],
   "source": [
    "#@title ## Parameters for data preparation\n",
    "#@markdown ### Job name\n",
    "#@markdown Name for the output folder\n",
    "out_name = \"poincareMSA\" #@param {type:\"string\"}\n",
    "\n",
    "#@markdown ### Threshold for filtering gapped positions\n",
    "#@markdown Positions with proportion of gaps above the given threshold are removed from the alignment. If your alignment is very gapped, you may want to increase this value.\n",
    "gapth = 0.9 #@param {type:\"number\"}\n",
    "\n",
    "#@markdown ## Run data preparation\n",
    "#@markdown Data preparation consists in `.mfasta` cleaning according to a gap threshold and translation of each sequence to the PSSM profile.\n",
    "prep_parameters = \"scripts/prepare_data\" + \" \" + mfasta + \" \" + out_name + \" \" + out_name + \" \" + str(gapth)\n",
    "bash_projection = \"bash scripts/prepare_data/create_projection.sh \" + prep_parameters\n",
    "!{bash_projection}"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "Dv0lg4Jn-dfG"
   },
   "source": [
    "# Projection"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "cellView": "form",
    "collapsed": true,
    "id": "B2_Cm3P9DJjD"
   },
   "outputs": [],
   "source": [
    "#@title ### Projection parameters\n",
    "#@markdown Here you control different parameters of Poincaré maps. In our computational experiments the best results were achieved for the following values provided by default. The impact of different parameters is analyzed in the original paper [1].\n",
    "knn = 5 #@param {type:\"number\"}\n",
    "gamma = 2 #@param {type:\"number\"}\n",
    "sigma = 1 #@param {type:\"number\"}\n",
    "cospca = 0 #@param {type:\"number\"}\n",
    "batchs = 4 #@param {type:\"number\"}\n",
    "epochs = 1000 #@param {type:\"number\"}\n",
    "seed = 0 #@param {type:\"number\"}\n",
    "\n",
    "#@markdown ## Building projection and preparing data for visualization\n",
    "#@markdown This step creates a projection of encoded sequences to a Poincaré disk.\n",
    "bash_pm = \"python3 \"+ \"scripts/build_poincare_map/main.py --input_path \" + out_name + \"/fasta\" + str(gapth) + \" --output_path \" + out_name + \"/projections/ --gamma \"+ str(gamma) +\" --pca \"+ str(cospca) + \" --epochs \"+ str(epochs) +\" --seed \"+ str(seed) + \" --knn \" + str(knn)\n",
    "!{bash_pm}"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "gRxHHfjX7CHF"
   },
   "source": [
    "# Projection visualization"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "cellView": "form",
    "id": "hm_w-BlvTSVY"
   },
   "outputs": [],
   "source": [
    "#@title ### Prepare data for visualization\n",
    "#@markdown Check if the annotation file is provided and prepares a dataframe for visualization.\n",
    "\n",
    "#Check that an annotation file was provided. Create a dummy one instead\n",
    "if not path_annotation:\n",
    "    df_annotation = pd.DataFrame(list(range(1,nb_seq+1)), columns=[\"proteins_id\"], dtype=\"float\")\n",
    "    df_annotation.to_csv(\"dummy_annotation.csv\", index=False)\n",
    "    path_annotation = \"dummy_annotation.csv\"\n",
    "    annotation_names = [\"proteins_id\"]\n",
    "\n",
    "#Create the DataFrame for visualization\n",
    "path_embedding = f\"{out_name}/projections/PM{knn:1.0f}sigma={sigma:2.2f}gamma={gamma:2.2f}cosinepca={cospca:1.0f}_seed{seed:1.0f}.csv\"\n",
    "df_embedding = read_embeddings(path_embedding, path_annotation, withroot=False)\n",
    "\n",
    "# save original coordinate before any rotation\n",
    "df_embedding[\"pm1_com\"] = df_embedding[\"pm1\"]\n",
    "df_embedding[\"pm2_com\"] = df_embedding[\"pm2\"]\n",
    "\n",
    "#@markdown Here are different labels found in your annotation file (if one uploaded):\n",
    "print(f\"{len(annotation_names)} annotations found: {annotation_names}.\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "cellView": "form",
    "id": "zLb9h9S4eKUD"
   },
   "outputs": [],
   "source": [
    "#@title ### Create interactive plot\n",
    "#@markdown Here you can set different parameters to color & annotate the resulting projection.\n",
    "#@markdown One can also provide a protein ID used to center your projection.\n",
    "\n",
    "\n",
    "title = \"\" #@param {type:\"string\"}\n",
    "\n",
    "#Labels name\n",
    "#@markdown ---\n",
    "#@markdown #### Select the coloring from annotation .csv file:\n",
    "labels_name = \"\" #@param {type:\"string\"}\n",
    "if labels_name == \"\":\n",
    "    labels_name = None\n",
    "elif labels_name not in annotation_names:\n",
    "    raise NameError(f\"labels_name {labels_name} is not in the availables annotations.\\nAvailables annotations: {annotation_names}\")\n",
    "\n",
    "#Labels text\n",
    "#@markdown #### Select classes to label among the \"labels_name\" or \"second_labels_name\" column (comma separated list):\n",
    "second_labels_name = \"\" #@param {type:\"string\"}\n",
    "if second_labels_name == \"\":\n",
    "    second_labels_name = None\n",
    "elif second_labels_name not in annotation_names:\n",
    "    raise NameError(f'\"second_labels_name\" {second_labels_name} is not in the availables annotations.\\nAvailables annotations: {annotation_names}')\n",
    "\n",
    "labels_text = \"\" #@param {type:\"string\"}\n",
    "if labels_text:\n",
    "    try:\n",
    "        labels_text = [s.strip() for s in labels_text.split(\",\")]\n",
    "    except:\n",
    "        print('Error: \"label_text\" field is not a valid list.')\n",
    "else:\n",
    "    labels_text = [\"\"]\n",
    "\n",
    "#Convert labels_text to labels_name dtype\n",
    "if labels_name and second_labels_name is None:\n",
    "    if labels_name and labels_text != [\"\"]:\n",
    "        try:\n",
    "            labels_text_dtype = df_annotation[labels_name].dtypes\n",
    "            labels_text = list(np.array(labels_text).astype(labels_text_dtype))\n",
    "        except:\n",
    "            raise TypeError(f'\"labels_text\" is not compatible with {labels_name}\" data format ({labels_text_dtype}).')\n",
    "else:\n",
    "    if second_labels_name and labels_text != [\"\"]:\n",
    "        try:\n",
    "            labels_text_dtype = df_annotation[second_labels_name].dtypes\n",
    "            labels_text = list(np.array(labels_text).astype(labels_text_dtype))\n",
    "        except:\n",
    "            raise TypeError(f'\"labels_text\" is not compatible with {second_labels_name}\" data format ({labels_text_dtype}).')\n",
    "\n",
    "show_text = True #@param {type:\"boolean\"}\n",
    "#@markdown ---\n",
    "\n",
    "#@markdown #### Use a custom color palette:\n",
    "color_palette = None #@param {type:\"raw\"}\n",
    "use_custom_palette = False #@param {type:\"boolean\"}\n",
    "\n",
    "if not use_custom_palette:\n",
    "    color_palette = None\n",
    "\n",
    "#@markdown ---\n",
    "\n",
    "#@markdown #### Select projection center:\n",
    "\n",
    "id_for_rotation = None #@param {type:\"number\"}\n",
    "use_selected_center = False #@param {type:\"boolean\"}\n",
    "\n",
    "tmp_id_for_rotation = -1\n",
    "if use_selected_center and id_for_rotation is not None:\n",
    "    if id_for_rotation != tmp_id_for_rotation:\n",
    "        tmp_id_for_rotation = id_for_rotation\n",
    "\n",
    "        poincare_coord = df_embedding[[\"pm1\",\"pm2\"]].sort_values(by=[\"proteins_id\"])\n",
    "        poincare_coord = poincare_coord.values[:,:2]\n",
    "        model = PoincareMaps(poincare_coord)\n",
    "        model.get_distances()\n",
    "\n",
    "        model.iroot = id_for_rotation - 1  #Zero indexing\n",
    "        model.rotate()\n",
    "\n",
    "        df_embedding[\"pm1\"] = model.coordinates_rotated[:,0]\n",
    "        df_embedding[\"pm2\"] = model.coordinates_rotated[:,1]\n",
    "else:\n",
    "    df_embedding[\"pm1\"] = df_embedding[\"pm1_com\"]\n",
    "    df_embedding[\"pm2\"] = df_embedding[\"pm2_com\"]\n",
    "\n",
    "#Plot graph\n",
    "fig = plot_embedding_interactive(df_embedding, \n",
    "                                 labels_name = labels_name,\n",
    "                                 second_labels_name = second_labels_name, \n",
    "                                 show_text = show_text,\n",
    "                                 labels_text = labels_text,\n",
    "                                 color_palette = color_palette, \n",
    "                                 title = title, \n",
    "                                 fontsize = 11)\n",
    "fig.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "cellView": "form",
    "id": "CGaefRbSWA0d"
   },
   "outputs": [],
   "source": [
    "#@title Save plot to file\n",
    "output_name = \"figure\" #@param {type:\"string\"}\n",
    "output_format = \"html\" #@param [\"png\", \"html\", \"pdf\", \"svg\"]\n",
    "\n",
    "if output_format != \"html\":\n",
    "    fig.write_image(f\"{output_name}.{output_format}\", engine=\"kaleido\")\n",
    "else:\n",
    "    fig.write_html(f\"{output_name}.{output_format}\")\n",
    "files.download(f\"{output_name}.{output_format}\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "cellView": "form",
    "id": "KB6XdhXCWA0d"
   },
   "outputs": [],
   "source": [
    "#@title Download intermediate data\n",
    "bash_command = f\"zip -r -q {out_name}.zip {out_name}\"\n",
    "!{bash_command}\n",
    "\n",
    "files.download(f\"{out_name}.zip\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "nCodoZxXWA0d"
   },
   "source": [
    "# Help"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "TlwImUjbN_CG"
   },
   "source": [
    "### Enabling the GPU\n",
    "\n",
    "To enable GPU in your notebook, select the following menu options −\n",
    "```\n",
    "Runtime / Change runtime type\n",
    "```\n",
    "\n",
    "<figure>\n",
    "<center>\n",
    "<img src=\"https://github.com/DSIMB/PoincareMSA/blob/master/.github/colab_gpu.png?raw=true\" width=500>\n",
    "</center>\n",
    "</figure>\n",
    "\n"
   ]
  }
 ],
 "metadata": {
  "accelerator": "GPU",
  "colab": {
   "include_colab_link": true,
   "provenance": []
  },
  "gpuClass": "standard",
  "kernelspec": {
   "display_name": "Python 3 (ipykernel)",
   "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.10.12"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 1
}