diff --git a/exercises/navis/README.md b/exercises/navis/README.md
new file mode 100644
index 0000000..deed601
--- /dev/null
+++ b/exercises/navis/README.md
@@ -0,0 +1,12 @@
+## NAVis exercises
+
+This folder contains the exercises for the NAVis part of the workshop:
+
+- `navis_on_colab_example.ipynb` simply demonstrates how to install NAVis on [Google Colab](https://colab.research.google.com/). You can open any notebook directly from Github by clicking "File" -> "Open notebook" -> "Github" and pasting in the link.
+
+- `navis_exercise_1.ipynb` contains the first set of exercises that will introduce some of the basic concepts in NAVis
+
+- `navis_exercise_2.ipynb` demonstrates how to re-create some of the figures and analyses from the [Gouwens, Sorensen, et al. (2020)](https://pubmed.ncbi.nlm.nih.gov/33186530/) paper on mouse visual cortex neurons. Please note that you will have to download data if you want to follow along - see the notebook or the workshop's [website](https://navis-org.github.io/neuropython2024/) for details.
+
+> [!IMPORTANT]
+> All exercises require a full install of NAVis: `pip install "navis[all]"`.
diff --git a/exercises/navis/_static/Fig4_A.png b/exercises/navis/_static/Fig4_A.png
new file mode 100644
index 0000000..4c758d9
Binary files /dev/null and b/exercises/navis/_static/Fig4_A.png differ
diff --git a/exercises/navis/_static/Fig5_A.jpg b/exercises/navis/_static/Fig5_A.jpg
new file mode 100644
index 0000000..4df6f07
Binary files /dev/null and b/exercises/navis/_static/Fig5_A.jpg differ
diff --git a/exercises/navis/_static/tortuosity.jpg b/exercises/navis/_static/tortuosity.jpg
new file mode 100644
index 0000000..6d54a9d
Binary files /dev/null and b/exercises/navis/_static/tortuosity.jpg differ
diff --git a/exercises/navis/navis_exercise_1.ipynb b/exercises/navis/navis_exercise_1.ipynb
new file mode 100644
index 0000000..eaac875
--- /dev/null
+++ b/exercises/navis/navis_exercise_1.ipynb
@@ -0,0 +1,844 @@
+{
+ "cells": [
+ {
+ "cell_type": "markdown",
+ "id": "b340428a-66d8-4d3b-8f26-dea2fb6a4c31",
+ "metadata": {
+ "editable": true,
+ "raw_mimetype": "",
+ "tags": []
+ },
+ "source": [
+ "![navis logo](https://github.com/navis-org/navis/raw/master/docs/_static/logo_new_banner.png)\n",
+ "\n",
+ "# NeuroPython Workshop 2024 - NAVis Exercise 1\n",
+ "\n",
+ "_Author: Philipp Schlegel (University of Cambridge)_\n",
+ "\n",
+ "_Date: September 23, 2024_\n",
+ "\n",
+ "_Source: https://github.com/navis-org/neuropython2024/exercises/navis/navis_exercise_1.ipynb_\n",
+ "\n",
+ "In this first exercise we will introduce the basic concepts of NAVis. You will learn about:\n",
+ "\n",
+ "1. Different neuron types \n",
+ "2. Converting between neuron meshes and skeletons\n",
+ "3. Neuron properties\n",
+ "4. Lists of neurons \n",
+ "5. Ways to plot neurons\n",
+ "\n",
+ "### Requirements \n",
+ "\n",
+ "For this exercise, you will need a full, up-to-date install of NAVis:\n",
+ "\n",
+ "```shell \n",
+ "pip install \"navis[all]\" -U\n",
+ "```\n",
+ "### Part I: Neuron Types & Properties\n",
+ "\n",
+ "NAVis supports 4 different types of neuron representations:\n",
+ "- `navis.MeshNeuron` = meshes\n",
+ "- `navis.TreeNeuron` = skeletons \n",
+ "- `navis.Dotprops` = dotprops \n",
+ "- `navis.VoxelNeuron` = images\n",
+ "\n",
+ "For this tutorial, we will focus in skeletons and meshes as these are the types you will most likely work with.\n",
+ "\n",
+ "#### Skeletons"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "id": "69d9112d",
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "import navis\n",
+ "\n",
+ "# This should be 1.8.0 or higher\n",
+ "print(navis.__version__)"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "id": "435b9ceb-fab2-4b24-8f9f-cc33f448aba9",
+ "metadata": {
+ "editable": true,
+ "tags": []
+ },
+ "outputs": [],
+ "source": [
+ "# Load one of the example skeletons\n",
+ "n = navis.example_neurons(1, kind=\"skeleton\")\n",
+ "\n",
+ "# Show a short summary of the skeleton\n",
+ "n"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "id": "822bfa92",
+ "metadata": {},
+ "source": [
+ "For skeletons, the underlying data is the SWC node table which you can access via the `.nodes` property:"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "id": "6f0e56f6",
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "# The nodes table is a simple pandas DataFrame\n",
+ "n.nodes.head()"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "id": "4b0c716d",
+ "metadata": {},
+ "source": [
+ "Properties such as the ones seen in the summary can be accessed similarly:"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "id": "e98727ab",
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "n.cable_length"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "id": "c27596d5",
+ "metadata": {},
+ "source": [
+ "Try some other parameter - e.g. `n.n_leafs`!\n",
+ "\n",
+ "You may have noticed that our example neuron has a `units` property of `8 nanometer`. "
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "id": "55c73547",
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "n.units"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "id": "56d02307",
+ "metadata": {},
+ "source": [
+ "If you know your neurons units, it's useful to set that property because it allows us to do stuff like this:"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "id": "f0e46637",
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "n.cable_length * n.units"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "id": "28c34a29",
+ "metadata": {},
+ "source": [
+ "Units can be set manually using a string:"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 7,
+ "id": "221ab0b6",
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "n.units = \"8 nm\""
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "id": "5caee7a7",
+ "metadata": {},
+ "source": [
+ "#### Somata\n",
+ "\n",
+ "You've probably already seen the `.soma` property on our example neuron:"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "id": "8eb757e4",
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "n.soma"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "id": "c34dd15d",
+ "metadata": {},
+ "source": [
+ "For skeletons, `.soma` is the ID of the node that corresponds to the neuron's soma. \n",
+ "\n",
+ "NAVis tries its best to infer the soma from whatever data you throw at it. For example when reading SWC files, it will look for a `label` column and failing that will look for nodes with a large radius. You can also set the soma manually:"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 9,
+ "id": "3755efc6",
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "n.soma = 4177"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "id": "80706176",
+ "metadata": {},
+ "source": [
+ "The `.soma` property is used for example for plotting and a few other functions, so it's always good to check this is set correctly."
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "id": "a143a8aa",
+ "metadata": {},
+ "source": [
+ "##### Roots\n",
+ "\n",
+ "Skeletons are directed acyclic graphs (= trees) where each node has exactly one parent. The only exception is the neuron's \"root\", the one node that does not have a parent:"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "id": "d1f0fdcb",
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "# Like `.soma`, these are node IDs\n",
+ "n.root"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "id": "f1717c2b",
+ "metadata": {},
+ "source": [
+ "Our example neuron is fully connected and so it only has a single root. If there are any breaks in the structure, you would find multiple root nodes.\n",
+ "\n",
+ "Why should we care about the skeleton's root? Roots are important because they are used as points of reference.\n",
+ "For example, when calculating Strahler Indices (SI) we walk from the leafs to the root and if that's not correctly set your SI won't make much sense. \n",
+ "\n",
+ "Currently, our neuron isn't actually rooted to its soma - let's fix that!"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "id": "a524fe8e",
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "# Reroot the skeleton at the soma\n",
+ "navis.reroot_skeleton(n, n.soma, inplace=True)\n",
+ "\n",
+ "# Check that the new root is correct\n",
+ "n.root"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "id": "8c7e7996",
+ "metadata": {},
+ "source": [
+ "Above function call also introduced an important concept: the `inplace` parameter. If you are used to `pandas` this should look familiar!\n",
+ "Many functions in NAVis accept this parameter:\n",
+ "\n",
+ "- `inplace=False` (default): the neuron(s) are copied and those copies are modifid and returned \n",
+ "- `inplace=True`: the original neurons are modified \n",
+ "\n",
+ "The latter is useful if you don't need to keep the originals - it saves time and memory by avoiding having to copy data. \n",
+ "\n",
+ "#### Meshes\n"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "id": "edc72e83",
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "# Load one of the example neuron meshes\n",
+ "m = navis.example_neurons(1, kind=\"mesh\")\n",
+ "\n",
+ "# Show a short summary of the mesh\n",
+ "m"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "id": "f70da58e",
+ "metadata": {},
+ "source": [
+ "As you can see the summary for meshes (i.e. `MeshNeurons`) looks quite a bit different from those for skeletons (i.e. `TreeNeurons`).\n",
+ "That's because in contrast to skeletons, neuron meshes don't immediately give us a sense of topology to calculate leafs, branch points, cable length \n",
+ "and so on. For that, we have to skeletonize them - which we'll demo below.\n",
+ "\n",
+ "The basic data for meshes are vertices (as x/y/z positions in space) and faces (as vertices indices):"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "id": "fd42a8c9",
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "m.vertices"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "id": "8946e2e8",
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "m.faces"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "id": "0864c926",
+ "metadata": {},
+ "source": [
+ "NAVis tries to make it easy to move between representations. In particular the conversion _to_ skeletons is very useful\n",
+ "as it allows us to analyse the neuron's topology.\n",
+ "\n",
+ "For `MeshNeurons` that's fairly straight forward:"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "id": "fc0596dd",
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "# Skeletonize the mesh\n",
+ "sk = navis.skeletonize(m)\n",
+ "sk"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "id": "ed29c36e",
+ "metadata": {},
+ "source": [
+ "**Important**: while the conversion function try to use sensible defaults, you may still want to inspect the results carefully and make adjustements if necessary!\n",
+ "\n",
+ "You may also notice that some functions such as e.g. `navis.strahler_index` accept `MeshNeurons` even though they probably operate on a skeleton:"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "id": "85d908b3",
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "# Calculate the Strahler index for the mesh\n",
+ "navis.strahler_index(m)\n",
+ "\n",
+ "# The per-vertex Strahler index is attached as a new property\n",
+ "m.strahler_index"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "id": "37c6c735",
+ "metadata": {},
+ "source": [
+ "The above works because NAVis (lazily) generates a skeleton, calculates the Strahler index on it and then maps the indices\n",
+ "back onto the mesh. The skeleton is then cached for future use:"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "id": "7d2b60d3",
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "# The .skeleton property is generated lazily on first access\n",
+ "m.skeleton"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "id": "4c69a576",
+ "metadata": {},
+ "source": [
+ "We won't cover this here but you can adjust the skeleton-representation manually - e.g. by using `m.skeletonize()`!\n",
+ "\n",
+ "### Lists of Neurons\n",
+ "\n",
+ "So far, we've only looked at single neurons but in reality you will likely work with a whole bunch at a time.\n",
+ "For that NAVis has `NeuronLists`:"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "id": "fcdd6623",
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "# Ask for multiple skeletons\n",
+ "nl = navis.example_neurons(3, kind=\"skeleton\")\n",
+ "\n",
+ "# You can also manually construct NeuronLists from individual neurons like so:\n",
+ "# nl = navis.NeuronList([n1, n2, n3])\n",
+ "\n",
+ "# Show a short summary of the list of skeletons\n",
+ "nl"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "id": "e8bf2825",
+ "metadata": {},
+ "source": [
+ "`NeuronLists` behave similar to numpy arrays in that we can index them like so:"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "id": "b0914afc",
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "# Get the first neuron in the list\n",
+ "nl[0]"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "id": "0f92e363",
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "# Index by a slice\n",
+ "nl[1:3]"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "id": "8bfd71ab",
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "# Index by a boolean array\n",
+ "nl[nl.cable_length > 300_000]"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "id": "fd2bb544",
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "# Iterate over the neurons in the list\n",
+ "for neuron in nl:\n",
+ " print(neuron.id)"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "id": "ab8fe3db",
+ "metadata": {},
+ "source": [
+ "In addition to the list/array-like ways to index into the `NeuronList`, you can also index by ID similar to `.loc` in `pandas`:"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "id": "5628615e",
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "# Grab a single neuron by ID\n",
+ "nl.idx[1734350908]"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "id": "946635d6",
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "# Grab multiple neurons by ID\n",
+ "nl.idx[[1734350908, 722817260]]"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "id": "7b20ad94",
+ "metadata": {},
+ "source": [
+ "So `NeuronLists` make it fairly easy to organize your neurons but what else are they good for?\n",
+ "\n",
+ "They also make it easy to access your neurons properties:"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "id": "0c4e3b88",
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "# Get the number of nodes in each neuron\n",
+ "nl.n_nodes"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "id": "d19c706d",
+ "metadata": {},
+ "source": [
+ "This access is generic and should work for any neuron property:"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "id": "bc1c9bc3",
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "nl.id"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "id": "1cd32228",
+ "metadata": {},
+ "source": [
+ "We can also use the `NeuronList` to add/modify neuron properties:"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "id": "5c687e1f",
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "# Set a \"my_property\" attribute on the neurons\n",
+ "nl.set_neuron_attributes([\"a\", \"b\", \"c\"], name=\"my_property\")\n",
+ "\n",
+ "nl[0].my_property"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "id": "b7825e93",
+ "metadata": {},
+ "source": [
+ "See also `NeuronList.add_meta_data()` for adding a whole bunch of properties at once.\n",
+ "\n",
+ "### Part II: Plotting\n",
+ "\n",
+ "NAVis has various functions for plotting neurons:\n",
+ "\n",
+ "- `navis.plot3d`\n",
+ "- `navis.plot2d` \n",
+ "- `navis.plot1d` \n",
+ "- `navis.plot_flat`\n",
+ "\n",
+ "Here, we will focus on `plot3d` and `plot2d` but check out the online tutorials to learn about the rest!\n",
+ "\n",
+ "\n",
+ "#### Plot 2D\n",
+ "\n",
+ "`navis.plot2d` uses `matplotlib` to visualize neurons. While `matplotlib` has 3d axes which allow you to plot\n",
+ "objects at oblique angles, its capabilities to render multiple objects correctly layered are limited. It's \n",
+ "also fairly slow for large scenes which means it's hard to use interactively. \n",
+ "\n",
+ "On the plus side, we get very nice, clean-looking figures and there are tons of ways of customizing them."
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "id": "3ecd6959",
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "# Plot a single neuron\n",
+ "fig, ax = navis.plot2d(\n",
+ " nl[0],\n",
+ " color=\"k\", # set the color\n",
+ " radius=True, # use the radius column to set the width of the lines\n",
+ " view=(\"x\", \"-y\"), # set the view\n",
+ " method=\"2d\", # this determines whether we use 2d or 3d matplotlib axes\n",
+ ")"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "id": "0de96b3b-1de0-4e40-b3ae-099ba3471444",
+ "metadata": {
+ "editable": true,
+ "tags": []
+ },
+ "outputs": [],
+ "source": [
+ "# Plot a bunch of neurons\n",
+ "fig, ax = navis.plot2d(\n",
+ " nl, # pass a list of neurons\n",
+ " palette=\"tab10\",\n",
+ " radius=True,\n",
+ " view=(\"x\", \"-y\"),\n",
+ " method=\"2d\",\n",
+ ")"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "id": "90f648f1",
+ "metadata": {},
+ "source": [
+ "You may have noticed that the skeletons look fairly plastic - that's because by default\n",
+ "NAVis will use the radius (if present). That \n",
+ "typically looks better at the cost of slower rendering. To switch off radii:"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "id": "42bee4e9",
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "# Plot a bunch of neurons\n",
+ "fig, ax = navis.plot2d(\n",
+ " nl,\n",
+ " palette=\"tab10\",\n",
+ " radius=False, # do not use the radius\n",
+ " linewidth=.5, # set the width of the lines\n",
+ " view=(\"x\", \"-y\"),\n",
+ " method=\"2d\",\n",
+ ")"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "id": "4fcd9b91",
+ "metadata": {},
+ "source": [
+ "You can always change the initial camera position using the `view` parameter. If you use 3D axes via `method=\"3d\"`,\n",
+ "you can also rotate the camera:"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "id": "1293d82b",
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "# Plot a bunch of neurons\n",
+ "fig, ax = navis.plot2d(\n",
+ " nl,\n",
+ " palette=\"tab10\",\n",
+ " radius=False,\n",
+ " linewidth=1,\n",
+ " view=(\"x\", \"-y\"),\n",
+ " method=\"3d\", # set method to 3d - try out `3d_complex`!\n",
+ " non_view_axes3d=True # show the non-view axes\n",
+ ")\n",
+ "\n",
+ "# Move the camera around\n",
+ "ax.elev = -20\n",
+ "ax.azim = 45\n",
+ "ax.roll = 180"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "id": "6534df2b",
+ "metadata": {},
+ "source": [
+ "#### Plot 3D \n",
+ "\n",
+ "Interactive 3D plots can be generated with `navis.plot3d`. Depending on what's installed on your system and whether you are in a Jupyter environment or in e.g. a terminal, NAVis will use different backends:\n",
+ "\n",
+ "- `plotly` is the default backend for Jupyter \n",
+ "- `octarine3d` will be used by default in the terminal or scripts \n",
+ "- `vispy` is the fallback for `octarine` "
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "id": "7108ae37",
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "# This should produce a plotly plot\n",
+ "fig = navis.plot3d(nl, palette='tab10')"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "id": "a991a12a",
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "# This should produce an octarine plot (if that backend is installed)\n",
+ "# Importantly: this will NOT work in Google Colab\n",
+ "viewer = navis.plot3d(nl, palette='tab10', backend='octarine')"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "id": "d91d2ca0",
+ "metadata": {},
+ "source": [
+ "Unlike `plotly` figures which are static once generated, `octarine` (and `vispy`) viewers can be manipulated further.\n",
+ "\n",
+ "Run these cells one-by-one and look at the viewer above:"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 34,
+ "id": "2806bf93",
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "# Set all neurons to red\n",
+ "viewer.set_colors('r')"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 35,
+ "id": "58e98ed1",
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "# Clear the viewer\n",
+ "viewer.clear()"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 36,
+ "id": "cb3c7d12",
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "# Add a single neuron\n",
+ "viewer.add(nl[0])"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 37,
+ "id": "b57e27b2",
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "# Change camera position\n",
+ "viewer.set_view(\"XY\") # this also accepts more complex, custom views (see docstring)"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "id": "ffbf4c06",
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "# Take a screenshot and display using matplotlib\n",
+ "# Note: this may fail with a shader-related error if you simply clicked \"Run All\" in a Jupyter notebook\n",
+ "# I think that's because we need to give it a split second to initialize (compile shaders, render a first frame, etc.)\n",
+ "# If you run this cell by itself, it should work fine\n",
+ "im = viewer.screenshot(filename=None)\n",
+ "im.shape"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "id": "e9791a68",
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "import matplotlib.pyplot as plt\n",
+ "\n",
+ "plt.imshow(im)"
+ ]
+ }
+ ],
+ "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.8"
+ },
+ "rise": {
+ "footer": "
Philipp
",
+ "header": "Hello
",
+ "show_buttons_on_startup": false,
+ "slideNumber": false,
+ "transition": "zoom"
+ }
+ },
+ "nbformat": 4,
+ "nbformat_minor": 5
+}
diff --git a/exercises/navis/navis_exercise_2.ipynb b/exercises/navis/navis_exercise_2.ipynb
new file mode 100644
index 0000000..5ed58dd
--- /dev/null
+++ b/exercises/navis/navis_exercise_2.ipynb
@@ -0,0 +1,940 @@
+{
+ "cells": [
+ {
+ "cell_type": "markdown",
+ "id": "b340428a-66d8-4d3b-8f26-dea2fb6a4c31",
+ "metadata": {
+ "editable": true,
+ "raw_mimetype": "",
+ "slideshow": {
+ "slide_type": "slide"
+ },
+ "tags": []
+ },
+ "source": [
+ "![navis logo](https://github.com/navis-org/navis/raw/master/docs/_static/logo_new_banner.png)\n",
+ "\n",
+ "# NeuroPython Workshop 2024 - NAVis Exercise 2\n",
+ "\n",
+ "_Author: Philipp Schlegel (University of Cambridge)_\n",
+ "\n",
+ "_Date: September 23, 2024_\n",
+ "\n",
+ "_Source: https://github.com/navis-org/neuropython2024/exercises/navis/navis_exercise_1.ipynb_\n",
+ "\n",
+ "In this second exercise we will analyse morphological data from [\"Integrated Morphoelectric and Transcriptomic Classification of Cortical GABAergic Cells\"](https://www.cell.com/cell/pdf/S0092-8674(20)31254-X.pdf) by Gouwens, Sorensen _et al._ (2020).\n",
+ "\n",
+ "You will learn to:\n",
+ "\n",
+ "1. Load neurons from disk and attach meta data \n",
+ "2. Correct offsets\n",
+ "3. Generate publication-ready plots \n",
+ "4. Extract basic neuron properties \n",
+ "5. Run a clustering\n",
+ "\n",
+ "### Requirements \n",
+ "\n",
+ "As with the previous exercise, you will need a full, up-to-date install of NAVis:\n",
+ "\n",
+ "```shell \n",
+ "pip install \"navis[all]\" -U\n",
+ "```\n",
+ "\n",
+ "In addition, you will need the folder with the SWC files and the `20200711_patchseq_metadata_mouse.csv` with the meta data.\n",
+ "Please see the [Data](https://navis-org.github.io/neuropython2024/preparing/#data) section on the website for details!"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "id": "d56f6b5a",
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "import navis\n",
+ "\n",
+ "# This should be 1.8.0 or higher\n",
+ "print(navis.__version__)"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "id": "d2207ea8",
+ "metadata": {},
+ "source": [
+ "### Part I: Loading neurons\n",
+ "\n",
+ "Importing skeletons from disk is pretty straight forward using `navis.read_swc`:"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "id": "da5a0b83",
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "# Read the SWC files in the folder\n",
+ "nl = navis.read_swc(\n",
+ " \"~/Downloads/swc/\", # adjust the filepath as needed!\n",
+ " fmt=\"{name,id:int}_transformed.swc\", # parse the name and id from the file name\n",
+ " limit=\".*_transformed.swc\", # only read files that end with `_transformed.swc`\n",
+ ")\n",
+ "nl"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "id": "3af52f78",
+ "metadata": {},
+ "source": [
+ "Next, we want to load the meta data provided alongside the reconstructions:"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "id": "72904979",
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "import pandas as pd\n",
+ "\n",
+ "meta = pd.read_csv(\"~/Downloads/20200711_patchseq_metadata_mouse.csv\")\n",
+ "\n",
+ "# Strangely, the last column with the transcriptomic (\"t\") types has no name - let's fix that:\n",
+ "meta.columns = list(meta.columns[:-1]) + [\"t_type\"]\n",
+ "\n",
+ "meta.head()"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "id": "7137a8a2",
+ "metadata": {},
+ "source": [
+ "The `cell_specimen_id` should match our neurons' `.id` property. The columns we care about are:\n",
+ "\n",
+ "- `cell_soma_normalized_depth` to help us position the neurons along the y-axis for plotting\n",
+ "- `neuron_reconstruction_type` to let us drop incomplete reconstructions \n",
+ "- `t_type` to let us correlate our morphological analyses with the transcriptomic types"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "id": "11ebc3b0",
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "# Attach metadata to neurons in the list\n",
+ "nl.add_metadata(\n",
+ " meta,\n",
+ " id_col=\"cell_specimen_id\", # the column in the metadata that corresponds to the neuron id\n",
+ " columns=[\n",
+ " \"t_type\",\n",
+ " \"neuron_reconstruction_type\",\n",
+ " \"cell_soma_normalized_depth\",\n",
+ " ], # only get our subset of columns\n",
+ " register=True, # use this meta data for the summary table\n",
+ ")\n",
+ "\n",
+ "nl"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "id": "86d4e38e",
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "# With the meta data attached we can do stuff like this:\n",
+ "nl[nl.t_type == \"Lamp5 Ntn1 Npy2r\"]"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "id": "a8bcad66",
+ "metadata": {},
+ "source": [
+ "Next, we need to align the neurons according to their soma depth! The normalized `cell_soma_normalized_depth` should map to a physical range of `0` to `922.5861720311`.\n",
+ "\n",
+ "Let's demo with one neuron before we run this for all neurons:"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "id": "33d2924c",
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "# Grab one of the neurons\n",
+ "n = nl[0]\n",
+ "\n",
+ "# This is the normalized soma depth\n",
+ "print(f\"Normalized soma depth: {n.cell_soma_normalized_depth}\")\n",
+ "\n",
+ "# This is the physical soma depth\n",
+ "# Note that we're positioning from the bottom, i.e. 922 will be the surface of the slice\n",
+ "# That's because matplotlib uses the bottom left as origin, so we need the \"deepest\" neurons\n",
+ "# to be at y = 0\n",
+ "phys_y = (1 - n.cell_soma_normalized_depth) * 922.5861720311\n",
+ "print(f\"Physical soma depth: {phys_y}\")\n",
+ "\n",
+ "# Current soma\n",
+ "print(f\"Current soma coordinates: {n.soma_pos[0]}\")"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "id": "c150eadb",
+ "metadata": {},
+ "source": [
+ "We will now offset the neuron such that the soma is at x=0 and y=368.6:"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "id": "7a884429",
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "offset = [0, phys_y, 0] - n.soma_pos[0]\n",
+ "offset"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "id": "cdc2da98",
+ "metadata": {},
+ "source": [
+ "Moving or scaling neurons in NAVis is super straight forward: adding, subtracting, dividing or multiplying neurons by a number or an `[x, y, z]` vector will change their coordinates. "
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "id": "968d6ce4",
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "# Move the neuron to the new centered position\n",
+ "n += offset\n",
+ "\n",
+ "# Check the that the soma is in the correct position\n",
+ "n.soma_pos[0]"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "id": "01e59737",
+ "metadata": {},
+ "source": [
+ "That looks good! Let's do it for all neurons:"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 9,
+ "id": "bbaebed3",
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "for n in nl:\n",
+ " phys_y = (1 - n.cell_soma_normalized_depth) * 922.5861720311\n",
+ " offset = [0, phys_y, 0] - n.soma_pos[0]\n",
+ " n += offset"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "id": "b43b5b2a",
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "# Check that soma positions are correct\n",
+ "nl.soma_pos"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "id": "b4327922",
+ "metadata": {},
+ "source": [
+ "## Part II: Plotting\n",
+ "\n",
+ "Now that we have loaded and aligned the neurons, let's see if we can recreate a plot similar to those in Figure 4A:\n",
+ "\n",
+ "![Figure4A](_static/Fig4_A.png)\n",
+ "\n",
+ "Here is the complete code for plotting in a single function for re-use. During the workshop we will build this slowly bit by bit to make it easier to understand!"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "id": "5dfd2ea6",
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "def plot_t_type(t_type, color=\"magenta\", axon_color=\"purple\", offset=500):\n",
+ " \"\"\"Plot all neurons of a given transcriptomic type.\n",
+ "\n",
+ " Parameters\n",
+ " ----------\n",
+ " t_type : str\n",
+ " The transcriptomic type to plot.\n",
+ " color : str\n",
+ " The color of the dendrites.\n",
+ " axon_color : str\n",
+ " The color of the axon.\n",
+ " offset : int\n",
+ " The offset between neurons along the x-axis.\n",
+ "\n",
+ " Returns\n",
+ " -------\n",
+ " fig, ax\n",
+ " The matplotlib figure and axis.\n",
+ "\n",
+ " \"\"\"\n",
+ " to_plot = nl[nl.t_type == \"Lamp5 Ntn1 Npy2r\"]\n",
+ "\n",
+ " # Offset the neurons slightly along the x-axis\n",
+ " to_plot = [n + [offset * i, 0, 0] for i, n in enumerate(to_plot)]\n",
+ "\n",
+ " compartment_palette = {i: color for i in range(5)}\n",
+ " compartment_palette[2] = axon_color\n",
+ "\n",
+ " # Plot the neuron\n",
+ " fig, ax = navis.plot2d(\n",
+ " to_plot,\n",
+ " radius=False, # use thin lines\n",
+ " lw=1,\n",
+ " c=\"black\",\n",
+ " method=\"2d\",\n",
+ " soma=dict(\n",
+ " fc=\"black\",\n",
+ " ec=\"white\", # highlight the soma with a white outline\n",
+ " radius=10, # override the default soma radius\n",
+ " ),\n",
+ " color_by=\"label\",\n",
+ " figsize=(\n",
+ " len(to_plot) * 2,\n",
+ " 10,\n",
+ " ), # scale the figure size with the number of neurons\n",
+ " palette=compartment_palette,\n",
+ " )\n",
+ "\n",
+ " # Add the layer boundaries\n",
+ " layer_bounds = {\n",
+ " \"L1\": 0,\n",
+ " \"L2/3\": 115.1112491335,\n",
+ " \"L4\": 333.4658190171,\n",
+ " \"L5\": 453.6227158132,\n",
+ " \"L6\": 687.6482650269,\n",
+ " \"L6b\": 883.1308910545,\n",
+ " }\n",
+ "\n",
+ " for layer, y in layer_bounds.items():\n",
+ " y = 922.5861720311 - y\n",
+ " ax.axhline(y, color=\"gray\", ls=\"--\", lw=1)\n",
+ " ax.text(-300, y - 25, layer, color=\"gray\", va=\"center\", size=11)\n",
+ "\n",
+ " ax.axhline(0, color=\"gray\", ls=\"--\", lw=1)\n",
+ "\n",
+ " # Set the axis y limits according to the layers\n",
+ " ax.set_ylim(-10, 930)\n",
+ "\n",
+ " # Hide axis\n",
+ " ax.axis(\"off\")\n",
+ "\n",
+ " return fig, ax\n",
+ "\n",
+ "\n",
+ "fig, ax = plot_t_type(\"Lamp5 Ntn1 Npy2r\")"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "id": "f685b045",
+ "metadata": {},
+ "source": [
+ "Next, let's calculate the distribution of cable lengths:"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "id": "22b588fb",
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "# Take one of the Lamp5 Ntn1 Npy2r neurons\n",
+ "to_plot = nl[nl.t_type == \"Lamp5 Ntn1 Npy2r\"]\n",
+ "\n",
+ "# We're going to be cheap here and simply generate a histogram over the node positions\n",
+ "# To make this representative, we should make sure that the number of nodes per unit of cable\n",
+ "# is homogeneous across neurons. For that we will resample the neurons:\n",
+ "print(\n",
+ " f\"Sampling rate (nodes/cable) before resampling: {to_plot.sampling_resolution.mean():.2f}\"\n",
+ ")\n",
+ "\n",
+ "# Resample to 2 nodes per micron\n",
+ "to_analyze = navis.resample_skeleton(\n",
+ " to_plot,\n",
+ " resample_to=0.5,\n",
+ " map_columns=\"label\", # make sure label column is carried over\n",
+ ")\n",
+ "print(\n",
+ " f\"Sampling rate (nodes/cable) after resampling: {to_analyze.sampling_resolution.mean():.2f}\"\n",
+ ")"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "id": "96a8c3e7",
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "# Get the combined nodes table\n",
+ "nodes = to_analyze.nodes\n",
+ "nodes.head()"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "id": "95e0d7cb",
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "import seaborn as sns\n",
+ "from mpl_toolkits.axes_grid1 import make_axes_locatable\n",
+ "\n",
+ "# Plot the neurons using the function we defined earlier\n",
+ "fig, ax = plot_t_type(\"Lamp5 Ntn1 Npy2r\")\n",
+ "\n",
+ "# Add a new axis to the right of the main plot\n",
+ "divider = make_axes_locatable(ax)\n",
+ "ax_hist = divider.append_axes(\"right\", size=0.75, pad=0.05)\n",
+ "\n",
+ "# Add histograms\n",
+ "# For axon:\n",
+ "sns.kdeplot(\n",
+ " data=nodes[nodes.label == 2], y=\"y\", ax=ax_hist, color=\"purple\", linewidth=1.5\n",
+ ")\n",
+ "# For the rest:\n",
+ "sns.kdeplot(\n",
+ " data=nodes[nodes.label != 2], y=\"y\", ax=ax_hist, color=\"magenta\", linewidth=1.5\n",
+ ")\n",
+ "\n",
+ "# Add soma positions\n",
+ "soma_pos = to_plot.soma_pos.reshape(-1, 3)\n",
+ "ax_hist.scatter([0] * len(soma_pos), soma_pos[:, 1], color=\"black\", s=10, clip_on=False)\n",
+ "\n",
+ "# Set same axis limits as the main plot\n",
+ "ax_hist.set_ylim(-10, 930)\n",
+ "\n",
+ "# Hide axes\n",
+ "ax_hist.set_axis_off()"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "id": "c9f7c64c",
+ "metadata": {},
+ "source": [
+ "## Part III: Morphometrics\n",
+ "\n",
+ "For the last part of this exercise, we will extract a number of features from each neuron and try to run a simple clustering. The idea is to create something like Figure 5A to compare transcriptomic with morphology types:\n",
+ "\n",
+ "![Fig5a](_static/Fig5_A.jpg)"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "id": "0fcfab84",
+ "metadata": {},
+ "source": [
+ "For their morphology types, Gouwens, Sorensen _et al._ used a collection of features collectively called \"IVSCC\" (In-Vitro Single Cell Characterization). The exact set of features seems to vary between studies but here are some examples:\n",
+ "\n",
+ "- Axon-soma distance: the path distance from the axon root to the soma surface (μm)\n",
+ "- Total length: the combined length of all branches (μm)\n",
+ "- Max. path distance: the path distance from the soma to the furthest neurite tip (μm)\n",
+ "\n",
+ "You can find a short description for each feature in the Supplements Note for [\"Classification of electrophysiological and\n",
+ "morphological neuron types in the mouse visual cortex\"](https://www.nature.com/articles/s41593-019-0417-0) by Gouwens, Sorensen, Berg _et al._ (2019).\n",
+ "\n",
+ "We will start by calculating a few of these by hand before switching of to letting NAVis do the heavy lifting!\n",
+ "\n",
+ "First item on the menu: total length!"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "id": "05c89a8a",
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "# Grab one example neuron to analyze\n",
+ "n = nl[0]\n",
+ "\n",
+ "print(f\"Total length: {n.cable_length:.1f} microns\")"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "id": "6010439d",
+ "metadata": {},
+ "source": [
+ "Ok, that was easy. But what if you wanted to get the path length for axon and dendrites separately?\n",
+ "\n",
+ "For that, we need to split the neuron using `navis.subset_neuron`:"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 16,
+ "id": "3c3e97fa",
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "n_axon = navis.subset_neuron(n, subset=n.nodes.label == 2) # axon only\n",
+ "n_dend = navis.subset_neuron(n, subset=n.nodes.label != 2) # everything but the axon"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "id": "69569bd6",
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "print(f\"Axon length: {n_axon.cable_length:.1f} microns\")\n",
+ "print(f\"Dendrites length: {n_dend.cable_length:.1f} microns\")"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "id": "9f70544c",
+ "metadata": {},
+ "source": [
+ "Great! How about something a bit more complicated next? How about max path distance, i.e. the longest distance between any two points in the neuron?\n",
+ "\n",
+ "In the paper, they calculate both the max Euclidean and the max geodesic (along-the-arbor) distance:"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "id": "e5007b39",
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "from scipy.spatial.distance import pdist\n",
+ "\n",
+ "# Get the maximum pairwise distance between nodes\n",
+ "max_distance = pdist(n.nodes[[\"x\", \"y\", \"z\"]]).max()\n",
+ "\n",
+ "print(f\"Maximum Euclidean path distance: {max_distance:.1f} microns\")"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "id": "83e0ff69",
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "# Calculate the geodesic distance between leaf nodes and all other nodes\n",
+ "m = navis.geodesic_matrix(n, from_=n.leafs.node_id.values)\n",
+ "\n",
+ "# Rows are our leaf nodes; columns are all other nodes\n",
+ "m.head()"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "id": "8fb0fb10",
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "# The maximum geodesic distance\n",
+ "max_distance = m.values.max()\n",
+ "\n",
+ "print(f\"Maximum geodesic path distance: {max_distance:.1f} microns\")"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "id": "9226f63b",
+ "metadata": {},
+ "source": [
+ "One important note regarding the `geodesic_matrix` above: if you neuron consists of multiple disconnected pieces, some of the distances will be infinite (`np.inf`). If that's the case, you need to either reconnect the neuron (`navis.heal_skeleton`) or exclude infinite values (`m.values[m.values < np.inf].max()`)!\n",
+ "\n",
+ "Let's do something slightly more complicated for the last example feature:\n",
+ "\n",
+ "- Mean contraction: the average of the ratios of the summed Euclidean distance between bifurcations, and between bifurcations and tips, to the summed path distance between same\n",
+ "\n",
+ "That explanation may sound a bit complicated but it's actually fairly simple: it is asking, for each linear segment between branch/leaf nodes, what's the ratio of the direct Euclidean vs the along-the-arbor geodesic distance:\n",
+ "\n",
+ "![skeleton](_static/tortuosity.jpg)\n",
+ "\n",
+ "Let's try calculating this manually, shall we?"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "id": "2a690e7a",
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "# NAVis pre-computes two kinds of segments:\n",
+ "print('Linear segmnents that maximmize segment lengths:', len(n.segments))\n",
+ "print('Linear segmnents between branch/leaf nodes:', len(n.small_segments))\n",
+ "\n",
+ "# For the tortuosity we need to use the latter!\n",
+ "\n",
+ "# Each segment is a list of node IDs from distal to proximal\n",
+ "n.small_segments[0]"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "id": "c9378fa9",
+ "metadata": {},
+ "source": [
+ "Let's first calculate the Euclidean distance:"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "id": "0473c7ec",
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "import numpy as np\n",
+ "\n",
+ "# Calculate the Euclidean distance between the first and the last node of each segment\n",
+ "node_locs = dict(zip(n.nodes.node_id, n.nodes[[\"x\", \"y\", \"z\"]].values))\n",
+ "dist_eucl = []\n",
+ "for seg in n.small_segments:\n",
+ " dist_eucl.append(\n",
+ " np.linalg.norm(node_locs[seg[0]] - node_locs[seg[-1]])\n",
+ " )\n",
+ "dist_eucl = np.array(dist_eucl)\n",
+ "\n",
+ "# Inspect the results\n",
+ "dist_eucl[:10]"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "id": "1cb9b42e",
+ "metadata": {},
+ "source": [
+ "Next, we need the geodesic distance. Here you have multiple options! The fastest is probably using `navis.geodesic_matrix` again but since we've already demo'ed that, let's run this operation on the neuron's graph representation instead:"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "id": "534c2afe",
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "# NAVis precomputed both iGraph and networkx graphs for the neuron\n",
+ "n.graph, n.igraph"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "id": "84158a3a",
+ "metadata": {},
+ "source": [
+ "iGraph is faster but we have to contend with mapping between node indices and vertex IDs which is annoying. Let's use networkx instead:"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "id": "8509decd",
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "import networkx as nx\n",
+ "\n",
+ "dist_geo = []\n",
+ "for seg in n.small_segments:\n",
+ " dist_geo.append(\n",
+ " nx.shortest_path_length(n.graph, source=seg[0], target=seg[-1], weight='weight')\n",
+ " )\n",
+ "dist_geo = np.array(dist_geo)\n",
+ "dist_geo[:10]"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "id": "86f77d87",
+ "metadata": {},
+ "source": [
+ "Great! The only thing that's left is to calculate the ratios and take the mean:"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "id": "0b9a46be",
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "tort = (dist_geo / dist_eucl).mean()\n",
+ "print(f\"Tortuosity: {tort:.3f}\")"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "id": "83cdd873",
+ "metadata": {},
+ "source": [
+ "For the record: NAVis has a built-in function to calculate tortuosity:"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "id": "bd8fdc6e",
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "navis.tortuosity(n)"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "id": "dc609df9",
+ "metadata": {},
+ "source": [
+ "We hope that manually compiling a couple features may have given you an idea of how this is done in general. In particular the use of the graph representation (`.graph` or `.igraph`) can be very powerful. We also encourage you to search the [API documentation](https://navis-org.github.io/navis/api/) to see if there already is a function for what you have in mind.\n",
+ "\n",
+ "Moving on, let's use NAVis to collect our IVSCC features:"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "id": "a3af08a6",
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "ivscc = navis.ivscc_features(n)\n",
+ "ivscc"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "id": "e3bf1888",
+ "metadata": {},
+ "source": [
+ "The `ivscc_features` function is new in NAVis `1.8.0` and currently calculates up to 57 features depending on which neuron compartments are present.\n",
+ "\n",
+ "Let's run this for all neurons in our list:"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "id": "04d009b2",
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "ivscc = navis.ivscc_features(nl)\n",
+ "ivscc.head()"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 29,
+ "id": "f9374f1e",
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "# Manually add the normalized soma depth to the ivscc table\n",
+ "ivscc.loc[\"cell_soma_normalized_depth\"] = nl.cell_soma_normalized_depth"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "id": "e42563ac",
+ "metadata": {},
+ "source": [
+ "Now that we have the features, we can run a clustering on them. There are obviously tons of options for doing that including many ways to pre-process the data. Like in the paper, we will run a hierarchical clustering:"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "id": "4b3ed274",
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "from scipy.spatial.distance import pdist\n",
+ "from scipy.cluster.hierarchy import linkage, cut_tree\n",
+ "\n",
+ "# Drop neurons that don't have a t-type assigned\n",
+ "nl_w_ttype = nl[nl.t_type != 'nan']\n",
+ "\n",
+ "# First we will drop neurons that have not t-type assigned\n",
+ "ivscc_w_ttype = ivscc[nl_w_ttype.id]\n",
+ "\n",
+ "# Our IVSCC features contain some NaNs where neurons don't have a basal or apical dendrite\n",
+ "# We could throw them out but let's just fill them with zeros\n",
+ "ivscc_filled = ivscc_w_ttype.fillna(0)\n",
+ "\n",
+ "# Normalize the data to be between 0 and 1\n",
+ "ivscc_scaled = (ivscc_filled + ivscc_filled.min()) / (ivscc_filled.max() - ivscc_filled.min())\n",
+ "\n",
+ "# Calculate the Euclidean distance\n",
+ "dist = pdist(ivscc_scaled.T)\n",
+ "\n",
+ "# Perform hierarchical clustering\n",
+ "Z = linkage(dist, method='ward')\n",
+ "\n",
+ "# Cut the tree to get the clusters\n",
+ "labels = cut_tree(Z, n_clusters=40).flatten()\n",
+ "\n",
+ "# Combine cluster labels and t-types into a dataframe\n",
+ "clusters = pd.DataFrame()\n",
+ "\n",
+ "clusters['id'] = nl_w_ttype.id\n",
+ "clusters['t_type'] = nl_w_ttype.t_type\n",
+ "clusters['cluster'] = labels\n",
+ "clusters.head()"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "id": "7a4226a3",
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "from scipy.cluster.hierarchy import leaves_list\n",
+ "\n",
+ "# For plotting we need to create a confusion matrix between the clusters and the t-types\n",
+ "confusion = pd.crosstab(clusters.t_type, clusters.cluster)\n",
+ "\n",
+ "# Sort rows and columns by hierarchical clustering\n",
+ "row_linkage = linkage(pdist(confusion), method='ward')\n",
+ "row_order = leaves_list(row_linkage)\n",
+ "col_linkage = linkage(pdist(confusion.T), method='ward')\n",
+ "col_order = leaves_list(col_linkage)\n",
+ "\n",
+ "# Apply the sorting\n",
+ "confusion = confusion.iloc[row_order, col_order]\n",
+ "\n",
+ "confusion.head()"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 32,
+ "id": "aba2a795",
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "# For each t-type, generate a color based on how the first label\n",
+ "first_labels = np.unique(confusion.index.map(lambda x: x.split(' ')[0]))\n",
+ "\n",
+ "cmap = {}\n",
+ "for t, base_color in zip(first_labels, sns.color_palette(\"tab10\", len(first_labels))):\n",
+ " # All labels starting with `t`\n",
+ " t_labels = [ix for ix in confusion.index if ix.startswith(t)]\n",
+ " cmap.update(dict(zip(t_labels,sns.light_palette(base_color, n_colors=len(t_labels) * 2)[::-1])))"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "id": "af9874db",
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "# Plot the confusion matrix as scatter heatmap\n",
+ "confusion_pivot = confusion.unstack().reset_index().rename(columns={0: 'count'})\n",
+ "\n",
+ "\n",
+ "# Turn the cluster into a string - otherwise seaborn will re-order it\n",
+ "confusion_pivot = confusion_pivot.astype({'cluster': \"string\", 't_type': \"string\"})\n",
+ "\n",
+ "g = sns.relplot(\n",
+ " data=confusion_pivot,\n",
+ " y=\"t_type\", x=\"cluster\", size=\"count\",\n",
+ " palette=cmap, hue=\"t_type\",\n",
+ " edgecolor=\".7\",\n",
+ " col_order=confusion.index.values, row_order=confusion.columns.values,\n",
+ " height=10, sizes=(0, 400),\n",
+ " legend=False\n",
+ ")\n",
+ "\n",
+ "ax = g.axes[0, 0]\n",
+ "\n",
+ "ax.set_ylabel('Transcriptomic type')\n",
+ "ax.set_xlabel('Morphological cluster')"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "id": "f038cf50",
+ "metadata": {},
+ "source": [
+ "That looks very different from the figure in the paper. However, in the paper they:\n",
+ "\n",
+ "- used more features\n",
+ "- removed low-variance features\n",
+ "- ran multiple clusterings and extracted consensus clusters\n",
+ "- merged small clusters (n<3) into the next larger cluster\n",
+ "- etc.\n",
+ "\n",
+ "Bottom line: it's not suprising we don't get the same results off the bat. Still, we hope this has given you an idea of how to run clustering on morphological features! Also check out the CAJAL tutorial for a very different approach on calculating morphological similarity!"
+ ]
+ }
+ ],
+ "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.8"
+ },
+ "rise": {
+ "footer": "Philipp
",
+ "header": "Hello
",
+ "show_buttons_on_startup": false,
+ "slideNumber": false,
+ "transition": "zoom"
+ }
+ },
+ "nbformat": 4,
+ "nbformat_minor": 5
+}
diff --git a/exercises/navis/navis_on_colab_example.ipynb b/exercises/navis/navis_on_colab_example.ipynb
new file mode 100644
index 0000000..6c3e778
--- /dev/null
+++ b/exercises/navis/navis_on_colab_example.ipynb
@@ -0,0 +1,113 @@
+{
+ "cells": [
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "Using `Navis` on Google Collaboratory is super straight forward:"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "# First we need to install the navis package\n",
+ "!pip install navis[all] -U"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "Now try running this cell to make sure `Navis` has installed correctly:"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "# Import and run some test code\n",
+ "import navis\n",
+ "\n",
+ "nl = navis.example_neurons(3, kind='skeleton')\n",
+ "nl"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "# This should render as plotly 3d plot\n",
+ "navis.plot3d(nl)"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "from IPython.display import display\n",
+ "\n",
+ "import navis\n",
+ "import seaborn as sns\n",
+ "nl = navis.example_neurons(n=3)\n",
+ "\n",
+ "# Some mock labels\n",
+ "labels = [\"typeA\", \"typeA\", \"typeB\"]\n",
+ "\n",
+ "# Generate a color for each unique label\n",
+ "cmap = dict(zip(set(labels), sns.color_palette('tab20', len(set(labels)))))\n",
+ "\n",
+ "# Map back onto neurons\n",
+ "colors = [cmap[l] for l in labels]\n",
+ "\n",
+ "# Plot\n",
+ "v = navis.plot3d(nl, color=colors, connectors=True, backend='vispy', cn_size=None)\n",
+ "#v.show()"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "from importlib import reload\n",
+ "reload(navis.plotting.dd)\n",
+ "reload(navis.plotting)\n",
+ "reload(navis)\n",
+ "\n",
+ "\n",
+ "fig, ax = navis.plot2d(nl, color=colors, method='3d', view=('x', '-z'), radius=True)"
+ ]
+ }
+ ],
+ "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.11.8"
+ }
+ },
+ "nbformat": 4,
+ "nbformat_minor": 2
+}
diff --git a/material/navis/navis_on_colab.ipynb b/material/navis/navis_on_colab.ipynb
deleted file mode 100644
index d948ed1..0000000
--- a/material/navis/navis_on_colab.ipynb
+++ /dev/null
@@ -1,41 +0,0 @@
-{
- "cells": [
- {
- "cell_type": "markdown",
- "metadata": {},
- "source": [
- "Using `Navis` on Google Collaboratory is super straight forward:"
- ]
- },
- {
- "cell_type": "code",
- "execution_count": null,
- "metadata": {},
- "outputs": [],
- "source": [
- "# First we need to install the navis package\n",
- "!pip install navis[all] -U"
- ]
- },
- {
- "cell_type": "code",
- "execution_count": null,
- "metadata": {},
- "outputs": [],
- "source": [
- "# Import and run some test code\n",
- "import navis\n",
- "\n",
- "nl = navis.example_neurons(3, kind='skeleton')\n",
- "nl"
- ]
- }
- ],
- "metadata": {
- "language_info": {
- "name": "python"
- }
- },
- "nbformat": 4,
- "nbformat_minor": 2
-}