notebooks/jupyter/classification/classify_by_drawingSelection.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "c6860758-e614-4e0b-b9bb-3cc55af78e9c",
   "metadata": {},
   "source": [
    "# Intro"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "82a37590-30ae-4d92-a889-0c173759fae9",
   "metadata": {},
   "source": [
    "Welcome to the interactive 'classification by drawing' notebook!\n",
    "This notebook goes through each step and allows you to tune parameters and view how it changes the results.\n",
    "\n",
    "The notebook proceeds as follows:\n",
    "\n",
    "1. **Import** libraries\n",
    "2. Define **paths** to data\n",
    "3. Run data through the **pipeline**. (ROInet embedding + UMAP)\n",
    "4. **Draw** a circle around a region of the UMAP embedding to select as 'good ROIs to keep'\n",
    "5. **Visualize** results\n",
    "6. **Save** results\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "c4a1f6ef-7075-4cc3-a046-dd170bb957b9",
   "metadata": {},
   "source": [
    "# Import libraries"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "452c74a7-c5bb-4b8b-ad89-b52323b131b0",
   "metadata": {},
   "source": [
    "Widen the notebook"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "bab4b303",
   "metadata": {},
   "outputs": [],
   "source": [
    "# widen jupyter notebook window\n",
    "from IPython.display import display, HTML\n",
    "display(HTML(\"<style>.container {width:95% !important; }</style>\"))\n",
    "display(HTML(\"<style>:root { --jp-notebook-max-width: 100% !important; }</style>\"))"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "d0fe2758-1ac9-4a38-986e-5fd0dc26746b",
   "metadata": {},
   "source": [
    "Import basic libraries"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "aaeb2017-7988-47b1-96d6-cc9b0c68edf3",
   "metadata": {},
   "outputs": [],
   "source": [
    "from pathlib import Path\n",
    "import tempfile\n",
    "\n",
    "import numpy as np\n",
    "from umap import UMAP"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "3e77df19-80a4-4fa9-b7df-e19647d0727a",
   "metadata": {},
   "source": [
    "Import `roicat`"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "9904f613-2094-48bd-8265-cfd7e0aed001",
   "metadata": {},
   "outputs": [],
   "source": [
    "import roicat"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "670046e1-e139-4a5b-a946-57ba7e59794c",
   "metadata": {},
   "source": [
    "# Find paths to data"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "872ed40f-14d1-47d5-82cb-27eb991e49e0",
   "metadata": {},
   "source": [
    "##### Prepare list of paths to data"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "469b96a5-db27-418e-8f40-f62a6760f6d6",
   "metadata": {},
   "source": [
    "In this example we are using suite2p output files, but other data types can be used (CaImAn, etc.) \\\n",
    "See the notebook on ingesting diverse data: https://github.com/RichieHakim/ROICaT/blob/main/notebooks/jupyter/other/demo_custom_data_importing.ipynb\n",
    "\n",
    "Make a list containing the paths to all the input files.\n",
    "\n",
    "In this example we are using suite2p, so the following are defined:\n",
    "1. `paths_allStat`: a list to all the stat.npy files\n",
    "2. `paths_allOps`: a list with ops.npy files that correspond 1-to-1 with the stat.npy files"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "4753932d-9180-4ead-a4b3-58f4269c0aa7",
   "metadata": {},
   "outputs": [],
   "source": [
    "dir_allOuterFolders = r'/media/rich/bigSSD/analysis_data/face_rhythm/mouse_0916N/'\n",
    "\n",
    "pathSuffixToStat = 'stat.npy'\n",
    "pathSuffixToOps = 'ops.npy'\n",
    "\n",
    "paths_allStat = roicat.helpers.find_paths(\n",
    "    dir_outer=dir_allOuterFolders,\n",
    "    reMatch=pathSuffixToStat,\n",
    "    depth=7,\n",
    ")[:]\n",
    "paths_allOps  = np.array([Path(path).resolve().parent / pathSuffixToOps for path in paths_allStat])[:]\n",
    "\n",
    "print(f'paths to all stat files:');\n",
    "[print(path) for path in paths_allStat];\n",
    "print('');\n",
    "print(f'paths to all ops files:');\n",
    "[print(path) for path in paths_allOps];\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "02e22710-d8ad-4345-a4d3-267a9c4455aa",
   "metadata": {},
   "source": [
    "**Important parameters**:\n",
    "\n",
    "- `um_per_pixel` (float):\n",
    "    - Resolution. 'micrometers per pixel' of the imaging field of view."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "a375034f-80ca-4b68-a4c7-13dff75a71ca",
   "metadata": {},
   "outputs": [],
   "source": [
    "data = roicat.data_importing.Data_suite2p(\n",
    "    paths_statFiles=paths_allStat,\n",
    "    paths_opsFiles=paths_allOps,\n",
    "    um_per_pixel=1.6,  \n",
    "    new_or_old_suite2p='new',\n",
    "    type_meanImg='meanImgE',\n",
    "    verbose=True,\n",
    ")\n",
    "\n",
    "assert data.check_completeness(verbose=False)['classification_inference'], f\"Data object is missing attributes necessary for tracking.\""
   ]
  },
  {
   "cell_type": "markdown",
   "id": "1fc72dd3-0b6a-448e-acf1-3935dbd86aff",
   "metadata": {},
   "source": [
    "# ROInet embedding\n",
    "\n",
    "This step passes the images of each ROI through the ROInet neural network. The inputs are the images, the output is an array describing the visual properties of each ROI."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "61b81261-1fba-4b0c-a248-bb6997623812",
   "metadata": {},
   "source": [
    "##### 1. Initialize ROInet\n",
    "\n",
    "Initialize the ROInet object. The `ROInet_embedder` class will automatically download and load a pretrained ROInet model. If you have a GPU, this step will be much faster."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "df4bbc7e-2517-4fb6-a59a-d31e46a2e30c",
   "metadata": {},
   "outputs": [],
   "source": [
    "DEVICE = roicat.helpers.set_device(use_GPU=True, verbose=True)\n",
    "dir_temp = tempfile.gettempdir()\n",
    "\n",
    "roinet = roicat.ROInet.ROInet_embedder(\n",
    "    device=DEVICE,  ## Which torch device to use ('cpu', 'cuda', etc.)\n",
    "    dir_networkFiles=dir_temp,  ## Directory to download the pretrained network to\n",
    "    download_method='check_local_first',  ## Check to see if a model has already been downloaded to the location (will skip if hash matches)\n",
    "    download_url='https://osf.io/c8m3b/download',  ## URL of the model\n",
    "    download_hash='357a8d9b630ec79f3e015d0056a4c2d5',  ## Hash of the model file\n",
    "    forward_pass_version='head',  ## How the data is passed through the network\n",
    "    verbose=True,  ## Whether to print updates\n",
    ")\n",
    "\n",
    "roinet.generate_dataloader(\n",
    "    ROI_images=data.ROI_images,  ## Input images of ROIs\n",
    "    um_per_pixel=data.um_per_pixel,  ## Resolution of FOV\n",
    "    pref_plot=False,  ## Whether or not to plot the ROI sizes\n",
    ");"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "0bbc5714",
   "metadata": {},
   "source": [
    "##### 2. Check ROI_images sizes\n",
    "In general, you want to see that a neuron fills roughly 25-50% of the area of the image. \\\n",
    "**Adjust `um_per_pixel` above to rescale image size**"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "5ffc9cd1-407b-4ee8-adbf-49dbe5664298",
   "metadata": {},
   "outputs": [],
   "source": [
    "roicat.visualization.display_toggle_image_stack(roinet.ROI_images_rs[:1000], image_size=(200,200))"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "7ee2b723",
   "metadata": {},
   "source": [
    "##### 3. Pass data through network\n",
    "\n",
    "Pass the data through the network. Expect for large datasets (~40,000 ROIs) that this takes around 15 minutes on CPU or 1 minute on GPU."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "71dbef15-4dee-4143-b4ec-2dc704bc1842",
   "metadata": {},
   "outputs": [],
   "source": [
    "roinet.generate_latents();"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "591f80f8-33bc-44ce-8bc8-bbf352bb0837",
   "metadata": {},
   "source": [
    "# UMAP embedding\n",
    "\n",
    "Reduce the dimensionality of the output of ROInet (~100 dims) to 2 dimensions so that we can visualize it. Feel free to use any settings here that do a good job of clustering your data as you see fit."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "4ce910f7",
   "metadata": {},
   "outputs": [],
   "source": [
    "umap = UMAP(\n",
    "    n_neighbors=25,\n",
    "    n_components=2,\n",
    "    n_epochs=400,\n",
    "    verbose=True,\n",
    "    densmap=False,\n",
    ")\n",
    "emb = umap.fit_transform(roinet.latents)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "cd9e6c9d-c4aa-4a20-ad15-d0329bbad201",
   "metadata": {},
   "source": [
    "# Draw selection"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "4a3debd6-e29f-497e-8aab-5c73b957899b",
   "metadata": {},
   "source": [
    "In order to visualize the kinds of ROIs at each region of the plot, we need to select a subset of points to overlay ROI images onto."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "802791ab-e0fc-424a-8320-d6a1e641f61d",
   "metadata": {},
   "outputs": [],
   "source": [
    "idx_images_overlay = roicat.visualization.get_spread_out_points(\n",
    "    emb,\n",
    "    n_ims=min(emb.shape[0], 1500),  ## Select number of overlayed images here\n",
    "    dist_im_to_point=0.8,\n",
    "#     border_frac=0.05,\n",
    ")\n",
    "\n",
    "images_overlay = roinet.ROI_images_rs[idx_images_overlay]"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "05cb7b22-c206-41e7-b1d6-55e58593f915",
   "metadata": {},
   "source": [
    "Now we can use an interactive plot (using the holoviews library) to select our region of the scatterplot to circle.\\\n",
    "This plot works as follows:\n",
    "- Use the **LASSO TOOL** to circle a region on the plot containing the images of ROIs that you'd like to keep/extract/mark.\n",
    "    - You can circle multiple times, but only the last one will be saved\n",
    "- The saved indices are saved in a temporary file that can be recovered using the `fn_get_indices` function output below. Just call `fn_get_indices()` and it will return a list of the integer indices.\n",
    "- If it is difficult to see the images, do the following:\n",
    "    - adjust the number of images in the above function (`roicat.visualization.get_spread_out_points`) using the `n_ims` argument\n",
    "    - adjust the overlap of the images in the below function (`roicat.visualization.select_region_scatterPlot`) using the `frac_overlap_allowed` argument"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "abcfedfa-6b15-42f2-a0ed-d244c61e7546",
   "metadata": {},
   "outputs": [],
   "source": [
    "fn_get_indices, layout, path_tempFile = roicat.visualization.select_region_scatterPlot(\n",
    "    data=emb,\n",
    "    idx_images_overlay=idx_images_overlay,\n",
    "    images_overlay=images_overlay[:, 6:30][:,:,6:30],\n",
    "    size_images_overlay=0.35,\n",
    "    frac_overlap_allowed=0.5,\n",
    "    figsize=(1200,1200),\n",
    "    alpha_points=1.0,\n",
    "    size_points=10,\n",
    "    color_points='b',\n",
    ");"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "ec7e6060-cd39-4666-9a39-e3b540f34562",
   "metadata": {},
   "source": [
    "Drop the results into easier to use output variables"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "9005978f-a9b4-4a56-87d1-fe6844e238b9",
   "metadata": {},
   "outputs": [],
   "source": [
    "# roicat.helpers.export_svg_hv_bokeh(layout, '/home/rich/Desktop/umap_with_labels_dotsOnly.svg')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "b9b82100",
   "metadata": {},
   "outputs": [],
   "source": [
    "n_sessions = len(data.ROI_images)\n",
    "idx_session_cat = np.concatenate([[ii]*data.ROI_images[ii].shape[0] for ii in range(n_sessions)])\n",
    "bool_good_cat = roicat.helpers.idx2bool(fn_get_indices(), length=len(idx_session_cat))\n",
    "preds_good_sessions = [np.int64((bool_good_cat * (idx_session_cat==ii))[idx_session_cat==ii]) for ii in range(data.n_sessions)]\n",
    "\n",
    "classification_output = {\n",
    "    'preds': preds_good_sessions,\n",
    "    'spatialFootprints': data.spatialFootprints,\n",
    "    'FOV_height': data.FOV_height,\n",
    "    'FOV_width': data.FOV_width,\n",
    "    'input_data': data.paths_stat,\n",
    "}"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "ec4640cd-047b-4bf3-90f8-6f4d72eda620",
   "metadata": {},
   "source": [
    "# Visualize outputs"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "be1f30c1-6013-4515-b137-04ee1020f87c",
   "metadata": {},
   "outputs": [],
   "source": [
    "print(f\"Number of 'good' and 'bad' ROIs from each session:\")\n",
    "print([f\"good: {p.sum()} / bad: {(p!=1).sum()}\" for p in preds_good_sessions])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "380e79cf-9722-4c24-9cc5-bb39d27a0f88",
   "metadata": {},
   "outputs": [],
   "source": [
    "FOVs_colored = roicat.visualization.compute_colored_FOV(\n",
    "    spatialFootprints=data.spatialFootprints,\n",
    "    FOV_height=data.FOV_height,\n",
    "    FOV_width=data.FOV_width,\n",
    "    labels=preds_good_sessions,\n",
    "    cmap=roicat.helpers.simple_cmap([[1,0,0],[0,1,0]]),\n",
    ")\n",
    "\n",
    "roicat.visualization.display_toggle_image_stack(\n",
    "    FOVs_colored,\n",
    "    image_size=(FOVs_colored[0].shape[0]*2, FOVs_colored[0].shape[1]*2)\n",
    ")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "53667b6f-0bb0-4738-8878-754998be0a55",
   "metadata": {},
   "source": [
    "# Save results\n",
    "\n",
    "The results file can be opened using any of the following methods:\n",
    "1. `roicat.helpers.pickle_load(path)`\n",
    "2. `np.load(path)`\n",
    "3. \n",
    "```\n",
    "    import pickle\n",
    "    with open(path_save, mode='rb') as f:\n",
    "        test = pickle.load(f)\n",
    "  ```"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "fe17605e-96da-4d5d-b111-0d96bb0888b4",
   "metadata": {},
   "outputs": [],
   "source": [
    "mouse = 'mouse_0916N'"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "ca2f1539",
   "metadata": {},
   "outputs": [],
   "source": [
    "dir_save = f'/media/rich/bigSSD/analysis_data/face_rhythm/{mouse}/'\n",
    "filename_save = mouse\n",
    "\n",
    "path_save = str(Path(dir_save).resolve() / (filename_save + '.ROICaT.classification_drawn.results' + '.pkl'))\n",
    "print(f'path_save: {path_save}')\n",
    "\n",
    "roicat.helpers.pickle_save(classification_output, path_save, mkdir=True)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "d3a77d9f-a95c-4533-9df7-07f105595630",
   "metadata": {},
   "source": [
    "# Thank you\n",
    "If you encountered any difficulties, please let us know at the issues page: https://github.com/RichieHakim/ROICaT/issues"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "0950587c-a4d1-4959-a167-e06265003bf6",
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "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.11.3"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}