accelerating_chemical_mappings.ipynb

{
  "cells": [
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "view-in-github",
        "colab_type": "text"
      },
      "source": [
        "<a href=\"https://colab.research.google.com/github/dcolinmorgan/grph/blob/main/accelerating_chemical_mappings.ipynb\" target=\"_parent\"><img src=\"https://colab.research.google.com/assets/colab-badge.svg\" alt=\"Open In Colab\"/></a>"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "dzmJAwfiAi6k"
      },
      "source": [
        "# Accelerated Chemical Mapping with [Graphistry](graphistry.com)\n",
        "\n",
        "This notebook visualizes a chemical dataset describing Blood Brain Barrier Permeability (BBBP) from [MoleculeNet](http://moleculenet.ai/datasets-1) and [ECFPs](https://pubs.acs.org/doi/10.1021/ci100050t), then reduce those high-dimensional vectors down to 2 dimensions using conventional statistical tools. Ultimately we demonstrate how such an OPEN-SOURCE analysis can be sped-up and scaled-up massively with the [graphistry](graphistry.com) environment and toolkit\n",
        "\n",
        "### [baseline HERE](https://github.com/dcolinmorgan/grph/blob/main/generic_chemical_mappings.ipynb)"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "w89wE473URRH"
      },
      "source": [
        "# Import accelerator libraries"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "ABA00KeDT6Gx"
      },
      "outputs": [],
      "source": [
        "!pip install --extra-index-url=https://pypi.nvidia.com cuml-cu11 cudf-cu11 cugraph-cu11 pylibraft_cu11 raft_dask_cu11 dask_cudf_cu11 pylibcugraph_cu11 pylibraft_cu11\n",
        "import cuml,cudf\n",
        "print(cuml.__version__)\n",
        "\n",
        "!pip install -U --force git+https://github.com/graphistry/pygraphistry.git@cudf\n",
        "!pip install -U git+https://github.com/graphistry/cu-cat.git@DT3\n",
        "# !pip install dirty_cat"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 2,
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/"
        },
        "id": "7_t1UL7YUAx3",
        "outputId": "e9df2c1a-51f8-45f5-ea14-ae338e9a710c"
      },
      "outputs": [
        {
          "output_type": "stream",
          "name": "stdout",
          "text": [
            "0.28.7+463.gfb96400\n",
            "/usr/local/lib/python3.10/dist-packages/cu_cat/__init__.py\n"
          ]
        }
      ],
      "source": [
        "import graphistry\n",
        "graphistry.register(api=3,protocol=\"https\", server=\"hub.graphistry.com\", username='dcolinmorgan', password='fXjJnkE3Gik6BWy') ## key id, secret key\n",
        "\n",
        "# graphistry.register(api=3,protocol=\"https\", server=\"hub.graphistry.com\", username='dcolinmorgan', password='***') ## key id, secret key\n",
        "print(graphistry.__version__)\n",
        "\n",
        "import cu_cat\n",
        "print(cu_cat.__file__)\n",
        "\n",
        "import os\n",
        "from collections import Counter\n",
        "import cProfile\n",
        "import pandas as pd\n",
        "import numpy as np\n",
        "import matplotlib.pyplot as plt\n",
        "from pstats import Stats\n",
        "import cuml,cudf\n",
        "from time import time\n",
        "import warnings\n",
        "warnings.filterwarnings('ignore')\n",
        "from typing import List\n",
        "import seaborn as sns\n",
        "pd.set_option('display.max_colwidth', 200)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 3,
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/"
        },
        "id": "EtU413LOUG_S",
        "outputId": "7f88b243-22fc-48d9-93c9-cdda300161ea"
      },
      "outputs": [
        {
          "output_type": "stream",
          "name": "stdout",
          "text": [
            "Tesla V100-SXM2-16GB\n"
          ]
        }
      ],
      "source": [
        "!nvidia-smi --query-gpu=gpu_name --format=csv,noheader"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "ZfM8WRfW4gOY"
      },
      "source": [
        "# Import Basics"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "prDjzDTU384B"
      },
      "outputs": [],
      "source": [
        "!pip install rdkit\n",
        "!pip install --pre deepchem\n",
        "\n",
        "from rdkit import Chem, DataStructs\n",
        "from rdkit.Chem.rdchem import Mol\n",
        "from rdkit.Chem.MolStandardize.rdMolStandardize import LargestFragmentChooser\n",
        "\n",
        "from rdkit import RDLogger\n",
        "lg = RDLogger.logger()\n",
        "lg.setLevel(RDLogger.CRITICAL)"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "L8NcwLVK5s_i"
      },
      "source": [
        "# Embed BBBP in Global Chemical Space Approximation (Dataset-Agnostic Embedding)"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "gdXebUL45usJ"
      },
      "source": [
        "### Read in and process ChEMBL data\n",
        "_Note: We want this walk-through to be fully reproducible with publicly available data. Since the dataset we used at Reverie create our internal Dataset-Agnostic Embedding can't be shared here, instead we select a random sample of ChEMBL compounds as our approximation of global chemical space._"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 10,
      "metadata": {
        "id": "HUbDr8Dp5uOB"
      },
      "outputs": [],
      "source": [
        "# Read in data from MoleculeNet\n",
        "chembl = pd.read_csv(\"https://deepchemdata.s3-us-west-1.amazonaws.com/datasets/chembl_sparse.csv.gz\", compression='gzip')\n",
        "\n",
        "# Sample a random 10k\n",
        "chembl = chembl.sample(n=20000)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 11,
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/"
        },
        "id": "DnXz47yx9MdU",
        "outputId": "f486c872-b250-4d80-9054-40023c6ce82d"
      },
      "outputs": [
        {
          "output_type": "stream",
          "name": "stderr",
          "text": [
            "WARNING:root:! Failed featurize speedup attempt. Continuing without memoization speedups.WARNING:root:Target is not of type(DataFrame) and has no columns"
          ]
        },
        {
          "output_type": "stream",
          "name": "stdout",
          "text": [
            "Using GPU: cu_cat\n"
          ]
        },
        {
          "output_type": "stream",
          "name": "stderr",
          "text": [
            "WARNING:root:! Failed umap speedup attempt. Continuing without memoization speedups.WARNING:graphistry.umap_utils:* Ignoring target column of shape (19951, 0) in UMAP fit, as it is not one dimensional"
          ]
        },
        {
          "output_type": "stream",
          "name": "stdout",
          "text": [
            "umap: \n",
            "46.321837425231934\n"
          ]
        }
      ],
      "source": [
        "\n",
        "chemblA = pd.read_csv(\"https://deepchemdata.s3-us-west-1.amazonaws.com/datasets/chembl_sparse.csv.gz\", compression='gzip')\n",
        "# p=''\n",
        "# for i in [150000,230000]:\n",
        "  # chembl = chemblA.sample(n=i)\n",
        "\n",
        "g2 = graphistry.nodes(cudf.from_pandas(chembl[\"smiles\"][chembl.smiles.str.len()<500]))\n",
        "\n",
        "  # t=time()\n",
        "  # g2=g2.featurize(feature_engine='cu_cat',memoize=True)\n",
        "  # print(\"featurize: \\n\"+str(time()-t))\n",
        "\n",
        "t=time()\n",
        "g4=g2.umap(feature_engine='cu_cat',engine='cuml',metric = \"jaccard\",\n",
        "                      n_neighbors = 25,\n",
        "                      n_components = 2,\n",
        "                      dbscan=True,\n",
        "                      # low_memory = True,\n",
        "                      min_dist = 0.001)\n",
        "j=time()-t\n",
        "print(\"umap: \\n\"+str(j))\n",
        "  # p=p+'_'+str(i)+':'+str(j)\n",
        "# g4.plot()"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "xb1xyWqn9pfM"
      },
      "outputs": [],
      "source": [
        "emb2=g4._node_embedding\n",
        "g22=graphistry.nodes(emb2.reset_index(),'ID').edges(g4._edges,'_src_implicit','_dst_implicit').bind(point_x=\"x\",point_y=\"y\").settings(url_params={\"play\":0})\n",
        "# g22.plot()"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "jovz3xVc4M8X"
      },
      "source": [
        "# Embed BBBP with UMAP"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "ZlQqzKMB4Qq7"
      },
      "source": [
        "### Read in and process small data"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "HR2S6BbP4JlB"
      },
      "outputs": [],
      "source": [
        "# Read in data from MoleculeNet\n",
        "bbbp = pd.read_csv(\"https://deepchemdata.s3-us-west-1.amazonaws.com/datasets/BBBP.csv\")\n",
        "\n",
        "# Clean up columnn names so they are easier to interpret\n",
        "bbbp = bbbp[[\"smiles\", \"p_np\", \"name\"]].reset_index(drop=True).rename({\"p_np\": \"permeable\"}, axis=1)\n",
        "\n",
        "# Remove extra fragments in SMILES (typically salts, which are irrelevant to BBB permeability)\n",
        "bbbp[\"smiles\"] = bbbp[\"smiles\"].apply(get_largest_fragment_from_smiles).dropna()\n",
        "t=time()\n",
        "# Compute descriptors and keep track of which failed to featurize\n",
        "ecfp_descriptors, keep_idx = compute_ecfp_descriptors(bbbp[\"smiles\"])\n",
        "\n",
        "# Only keep those that sucessfully featurized\n",
        "bbbp = bbbp.iloc[keep_idx]\n",
        "print(time()-t)\n",
        "print(bbbp.shape)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "i9iRe44glrPq"
      },
      "outputs": [],
      "source": [
        "BBBP=bbbp[~bbbp.name.duplicated(keep='first')]\n",
        "BBBP[['name','permeable']][BBBP.smiles.str.len()>3]#.reset_index(drop=True)"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "6L8TZ5qx530H"
      },
      "source": [
        "### ... and with graphistry"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "fdvAXUCi9oGZ"
      },
      "outputs": [],
      "source": [
        "ecfp_descriptors"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "vA8LG2dFgKrc"
      },
      "outputs": [],
      "source": [
        "BBBP=bbbp[~bbbp.name.duplicated(keep='first')]\n",
        "\n",
        "g = graphistry.nodes(cudf.from_pandas(BBBP[['smiles','permeable']][BBBP.smiles.str.len()>3]))\n",
        "t=time()\n",
        "# g2=g.featurize(feature_engine='cu_cat',memoize=True)\n",
        "g3=g.umap(feature_engine='cu_cat',engine='cuml',metric = \"jaccard\",\n",
        "                      n_neighbors = 25,\n",
        "                      n_components = 2,\n",
        "                      low_memory = False,\n",
        "                      min_dist = 0.001)\n",
        "print(\"\\n\"+str(time()-t))"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "WMzz0EfaqR03"
      },
      "outputs": [],
      "source": [
        "# g2.bind(node=BBBP.permeable) #[BBBP.smiles.str.len()>3])\n",
        "g3.encode_point_color('permeable',palette=[\"hotpink\", \"dodgerblue\"],as_continuous=True).plot()\n"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "gI0YsM12V3Nt"
      },
      "source": [
        "## using a [universal fingerprint](https://jcheminf.biomedcentral.com/articles/10.1186/s13321-020-00445-4) rather than [morgan fingerprint from 1965](https://pubs.acs.org/doi/abs/10.1021/c160017a018)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "m26BEG2HVBdv"
      },
      "outputs": [],
      "source": [
        "import locale\n",
        "def getpreferredencoding(do_setlocale = True):\n",
        "    return \"UTF-8\"\n",
        "locale.getpreferredencoding = getpreferredencoding\n",
        "\n",
        "!pip install tmap\n",
        "!pip install git+https://github.com/reymond-group/map4@v1.0\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "TIn7kPdTCT5P"
      },
      "outputs": [],
      "source": [
        "# !pip install -q condacolab\n",
        "# import condacolab\n",
        "condacolab.install()\n",
        "\n",
        "!conda install -c tmap tmap"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "Ep_UAw7CQ0S0"
      },
      "outputs": [],
      "source": [
        "from rdkit import Chem\n",
        "import tmap as tm\n",
        "from map4 import MAP4Calculator\n",
        "\n",
        "dim = 1024\n",
        "\n",
        "MAP4 = MAP4Calculator(dimensions=dim)\n",
        "ENC = tm.Minhash(dim)\n",
        "\n",
        "smiles_a = 'c1ccccc1'\n",
        "mol_a = Chem.MolFromSmiles(smiles_a)\n",
        "map4_a = MAP4.calculate(mol_a)\n",
        "\n"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "v4Ov_xoI51tG"
      },
      "source": [
        "### ... and with graphistry"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "t7qqMDBbQ0VH"
      },
      "outputs": [],
      "source": [
        "# PPB = PPB[PPB.columns.drop(list(PPB.filter(regex='UMAP')))]\n",
        "g = graphistry.nodes(cudf.from_pandas(PPB[['smiles']][PPB.smiles.str.len()>3]))\n",
        "\n",
        "t=time()\n",
        "# g2=g.featurize(feature_engine='cu_cat',memoize=True)\n",
        "g3=g.umap(feature_engine='cu_cat',engine='cuml',metric = \"jaccard\",\n",
        "                      n_neighbors = 25,\n",
        "                      n_components = 2,\n",
        "                      low_memory = False,\n",
        "                      min_dist = 0.001)\n",
        "print(\"\\n\"+str(time()-t))\n",
        "emb2=g3._node_embedding\n",
        "g222=graphistry.nodes(emb2.reset_index(),'index').edges(g3._edges,'_src_implicit','_dst_implicit').bind(point_x=\"x\",point_y=\"y\").settings(url_params={\"play\":0})\n",
        "g222.plot()"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "f4X6aRSQUkBn"
      },
      "source": [
        "## and larger data"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "S4rnM7NmM17T"
      },
      "outputs": [],
      "source": [
        "PPB = pd.read_csv(\"https://deepchemdata.s3-us-west-1.amazonaws.com/datasets/qm9.csv\")\n",
        "t=time()\n",
        "#Compute ECFP Descriptors\n",
        "ecfp_descriptors, keep_idx = compute_ecfp_descriptors(PPB[\"smiles\"])\n",
        "PPB = PPB.iloc[keep_idx]\n",
        "print(PPB.shape)\n",
        "print(time()-t)"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "_M-ZH2jJRD6y"
      },
      "source": [
        "### ... and with graphistry"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "5EHc3y0jVBba"
      },
      "outputs": [],
      "source": [
        "# PPB = PPB[PPB.columns.drop(list(PPB.filter(regex='UMAP')))]\n",
        "g = graphistry.nodes(cudf.from_pandas(PPB[['smiles']][PPB.smiles.str.len()>3]))\n",
        "t=time()\n",
        "# g2=g.featurize(feature_engine='cu_cat',memoize=True)\n",
        "g3=g.umap(feature_engine='cu_cat',engine='cuml',metric = \"jaccard\",\n",
        "                      n_neighbors = 25,\n",
        "                      n_components = 2,\n",
        "                      low_memory = False,\n",
        "                      min_dist = 0.001)\n",
        "\n",
        "print(\"\\n\"+str(time()-t))\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "FdbQoMjHQ8D-"
      },
      "outputs": [],
      "source": [
        "emb2=g3._node_embedding\n",
        "g222=graphistry.nodes(emb2.reset_index(),'index').edges(g3._edges,'_src_implicit','_dst_implicit').bind(point_x=\"x\",point_y=\"y\").settings(url_params={\"play\":0})\n",
        "g222.plot()"
      ]
    }
  ],
  "metadata": {
    "accelerator": "GPU",
    "colab": {
      "collapsed_sections": [
        "jovz3xVc4M8X"
      ],
      "machine_shape": "hm",
      "provenance": [],
      "gpuType": "V100",
      "include_colab_link": true
    },
    "kernelspec": {
      "display_name": "Python 3",
      "name": "python3"
    }
  },
  "nbformat": 4,
  "nbformat_minor": 0
}