PyUUL-napari.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "b244798c",
   "metadata": {},
   "source": [
    "First of all we import the libraries we need"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "717102f8",
   "metadata": {},
   "source": [
    "https://www.nature.com/articles/s41467-022-28327-3"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "id": "bb375a08",
   "metadata": {},
   "outputs": [],
   "source": [
    "from pyuul import utils, VolumeMaker\n",
    "import os, urllib, torch\n",
    "import numpy as np\n",
    "import napari\n",
    "\n",
    "\n",
    "viewer = napari.Viewer()"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "3dd917b1",
   "metadata": {},
   "source": [
    "We fetch the protein structures and we parse the PDB using the pyuul utils module"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "id": "152936d1",
   "metadata": {},
   "outputs": [],
   "source": [
    "# os.makedirs('exampleStructures', exist_ok=True)\n",
    "# urllib.request.urlretrieve('http://files.rcsb.org/download/101M.pdb', 'exampleStructures/101m.pdb')\n",
    "# urllib.request.urlretrieve('http://files.rcsb.org/download/5BMZ.pdb', 'exampleStructures/5bmz.pdb')\n",
    "# urllib.request.urlretrieve('http://files.rcsb.org/download/5BOX.pdb', 'exampleStructures/5box.pdb')\n",
    "# urllib.request.urlretrieve('http://files.rcsb.org/download/1WOU.pdb', 'exampleStructures/1wou.pdb')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "id": "d289af1f",
   "metadata": {},
   "outputs": [],
   "source": [
    "device = \"cpu\" # runs the volumes on CPU\n",
    "\n",
    "#coords, atname = utils.parsePDB(\"exampleStructures/\") # get coordinates and atom names\n",
    "coords, atname = utils.parsePDB(\"exampleStructures/1wou.pdb\") # get coordinates and atom names\n",
    "atoms_channel = utils.atomlistToChannels(atname) # calculates the corresponding channel of each atom\n",
    "radius = utils.atomlistToRadius(atname) # calculates the radius of each atom"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "id": "154fb47f",
   "metadata": {},
   "outputs": [],
   "source": [
    "coords = coords.to(device)\n",
    "radius = radius.to(device)\n",
    "atoms_channel = atoms_channel.to(device)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "id": "977060b0",
   "metadata": {},
   "outputs": [],
   "source": [
    "VoxelsObject = VolumeMaker.Voxels(device=device, sparse=True)\n",
    "PointCloudSurfaceObject = VolumeMaker.PointCloudVolume(device=device)\n",
    "PointCloudVolumeObject = VolumeMaker.PointCloudSurface(device=device)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "ea4cb6a1",
   "metadata": {},
   "source": [
    "### Volumetric representations"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "id": "0c50f396",
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "/Users/nsofroniew/opt/anaconda3/envs/mol/lib/python3.9/site-packages/pyuul/VolumeMaker.py:564: UserWarning: __rfloordiv__ is deprecated, and its behavior will change in a future version of pytorch. It currently rounds toward 0 (like the 'trunc' function NOT 'floor'). This results in incorrect rounding for negative values. To keep the current behavior, use torch.div(a, b, rounding_mode='trunc'), or for actual floor division, use torch.div(a, b, rounding_mode='floor').\n",
      "  npoints = maxpoints // padding_mask.sum(-1).min() + 1 # we ensure that the smallest protein has at least 5000 points\n",
      "/Users/nsofroniew/opt/anaconda3/envs/mol/lib/python3.9/site-packages/pyuul/VolumeMaker.py:463: UserWarning: __rfloordiv__ is deprecated, and its behavior will change in a future version of pytorch. It currently rounds toward 0 (like the 'trunc' function NOT 'floor'). This results in incorrect rounding for negative values. To keep the current behavior, use torch.div(a, b, rounding_mode='trunc'), or for actual floor division, use torch.div(a, b, rounding_mode='floor').\n",
      "  npoints = (maxpoints // padding_mask.sum(-1).min() + 1) * 2  # we ensure that the smallest protein has at least maxpoints points\n"
     ]
    }
   ],
   "source": [
    "SurfacePointCloud = PointCloudSurfaceObject(coords, radius)\n",
    "VolumePointCloud = PointCloudVolumeObject(coords, radius)\n",
    "VoxelRepresentation = VoxelsObject(coords, radius, atoms_channel)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "ae601085",
   "metadata": {},
   "source": [
    "### Visualize the surface point cloud data"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "id": "5f522f03",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<Points layer 'point_cloud' at 0x133c30f70>"
      ]
     },
     "execution_count": 8,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# seperate dim for each batch\n",
    "pts = []\n",
    "for i, pc in enumerate(SurfacePointCloud.numpy()):\n",
    "    pts_padded = np.concatenate((np.tile(i, (len(pc), 1)), pc), axis=1)\n",
    "    pts.append(pts_padded)\n",
    "point_cloud = np.concatenate(pts, axis=0)\n",
    "\n",
    "viewer.add_points(point_cloud, size=0.9, edge_width=0)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "id": "8a2edebb",
   "metadata": {},
   "outputs": [],
   "source": [
    "# seperate layer for each batch\n",
    "for spc in SurfacePointCloud.numpy():\n",
    "    viewer.add_points(spc, size=0.9, edge_width=0)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "id": "39874b87",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<Points layer 'spc' at 0x13a727280>"
      ]
     },
     "execution_count": 10,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# one 3D layer for all\n",
    "spc = SurfacePointCloud.numpy().reshape(-1, 3)\n",
    "viewer.add_points(spc, size=0.9, edge_width=0)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "ef911a30",
   "metadata": {},
   "source": [
    "### Visualize the volume point cloud data"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "id": "49dd0787",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<Points layer 'vpc' at 0x13ad39e80>"
      ]
     },
     "execution_count": 11,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "vpc = VolumePointCloud.numpy()\n",
    "viewer.add_points(vpc, size=0.9, edge_width=0, face_color='red')"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "b8813834",
   "metadata": {},
   "source": [
    "### Visualize the voxel representation"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 37,
   "id": "ab564a41",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<Image layer 'vr [3]' at 0x13c2d0820>"
      ]
     },
     "execution_count": 37,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "vr = VoxelRepresentation.to_dense().numpy()\n",
    "# translate = np.concatenate((np.zeros(vr.ndim - 3), -np.divide(vr.shape[-2:], 4)), axis=0)\n",
    "translate = [0, 0, 0, 0, 0]\n",
    "viewer.add_image(vr, translate=translate)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 32,
   "id": "acb41f1b",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<Image layer 'vr [1]' at 0x13ac92310>"
      ]
     },
     "execution_count": 32,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "viewer.add_image(vr)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 28,
   "id": "03585ff9",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([22., 24., 26.])"
      ]
     },
     "execution_count": 28,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "np.divide(vr.shape[-3:], 2)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 50,
   "id": "9b0427cc",
   "metadata": {},
   "outputs": [],
   "source": [
    "# viewer.layers[-1].translate = [0, 0, -11, 3, -13]"
   ]
  }
 ],
 "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.9.10"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}