CellTypesNotebook-CodingFocused.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Defining cell types by their electrophysiology\n",
    "\n",
    "### The brain has thousands of different types of cells. How do we even begin to tease them apart?\n",
    "\n",
    "We can define neurons by their <b>gene expression patterns</b>, <b>electrophysiology features</b>, and <b>structure</b>. Here, we'll use those three features to compare and contrast cell types in the brain.\n",
    "\n",
    "This notebook will help us investigate specific features in the electrophysiology dataset from the Allen Brain Atlas. Along the way, you'll encounter some concepts about coding in Python.\n",
    "\n",
    "### Before diving into this notebook:\n",
    "- Complete [the web-based activity](https://sites.google.com/ucsd.edu/neuroedu/cell-types/web-instructions?authuser=0), which asks you to look at this data on the Allen Institute website.\n",
    "- Work through [this Introduction to Jupyter Notebooks tutorial](https://github.com/ajuavinett/CellTypesLesson/blob/master/Introduction%20to%20Jupyter%20Notebooks.ipynb) which will introduce you to Python and the Jupyter Notebook environment.\n",
    "\n",
    "<hr>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Tutorial Steps\n",
    "1. [Set up our coding environment](#setup)\n",
    "\n",
    "**Data from a single cell**. *Steps 2-4 are an optional demonstration of the same cell from the [web-based activity](https://sites.google.com/ucsd.edu/neuroedu/cell-types?authuser=0). These steps can be skipped.*\n",
    "\n",
    "2. [Import data for a single cell](#import)\n",
    "3. [Plot a raw sweep of data](#plotsweep)\n",
    "4. [Plot the morphology of the cell](#morphology)\n",
    "\n",
    "**Data from many cells**. *Steps 5-7 will help you look at pre-computed features for all of the cells in the database.*\n",
    "\n",
    "5. [Analyze computed features](#metrics)\n",
    "6. [Compare action potential waveforms](#waveforms)\n",
    "7. [Compare cell types](#compare)\n",
    "<hr>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<a id=\"setup\"></a>\n",
    "## Step 1. Set up our coding environment\n",
    "Whenever we start an analysis in Python, we need to be sure to import the necessary code packages. If you're running this notebook on Colab or Binder, the cells below will install packages into your coding environment -- we are *not* installing anything on your computer.\n",
    "\n",
    "### Install the AllenSDK\n",
    "The Allen Institute has compiled a set of code and tools called a **Software Development Kit** (SDK). These tools will help us import and analyze the cell types data. See [Technical Notes](#technical) at the end of this notebook for more information about working with the AllenSDK.\n",
    "\n",
    "><b>Task</b>: Run the cell below, which will install the allensdk into your coding environment.\n",
    "\n",
    "**Technical notes about installing the allensdk**\n",
    "- If you're running this in Colab, you'll also be prompted to **Restart Runtime** after this is completed. Click the **Restart Runtime** button (nothing will happen), and then you're ready to proceed.\n",
    "- If you receive an error or are running this notebook on your local computer, there are additional instructions on how to install the SDK locally <a href=\"http://alleninstitute.github.io/AllenSDK/install.html\">here</a>."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# This will ensure that the AllenSDK is installed.\n",
    "try:\n",
    "    import allensdk\n",
    "    print('allensdk imported')\n",
    "except ImportError as e:\n",
    "    !pip install allensdk"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "><b>Task</b>: We also need to make sure that our coding environment has [NumPy](https://numpy.org/), [Pandas](https://pandas.pydata.org/), and [Matplotlib](https://matplotlib.org/) already installed. Run the cell below -- any packages that are missing will be installed for you."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "numpy already installed\n",
      "pandas already installed\n",
      "matplotlib already installed\n"
     ]
    }
   ],
   "source": [
    "# This will ensure that NumPy, Pandas, and Matplotlib are installed.\n",
    "try:\n",
    "    import numpy\n",
    "    print('numpy already installed')\n",
    "except ImportError as e:\n",
    "    !pip install numpy\n",
    "try:\n",
    "    import pandas\n",
    "    print('pandas already installed')\n",
    "except ImportError as e:\n",
    "    !pip install pandas  \n",
    "try:\n",
    "    import matplotlib\n",
    "    print('matplotlib already installed')\n",
    "except ImportError as e:\n",
    "    !pip install matplotlib"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Import common packages\n",
    "Below, we'll `import` a common selection of packages that will help us analyze and plot our data. We'll also configure the plotting in our notebook.\n",
    "><b>Task</b>: Just like you did in the [introductory notebook](https://github.com/ajuavinett/CellTypesLesson/blob/master/Introduction%20to%20Jupyter%20Notebooks.ipynb), import the numpy module nicknamed as <code>np</code>. Add a <code>print</code> message at the end that says \"Packages imported!\" so that you know the code ran."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Import our plotting package from matplotlib\n",
    "import matplotlib.pyplot as plt\n",
    "\n",
    "# Specify that all plots will happen inline & in high resolution\n",
    "%matplotlib inline  \n",
    "%config InlineBackend.figure_format = 'retina'\n",
    "\n",
    "# Import pandas for working with databases\n",
    "import pandas as pd\n",
    "\n",
    "# Import numpy below\n",
    "\n",
    "\n",
    "# Add your print() statement below\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Import the CellTypesModule from the allensdk\n",
    "With the allensdk installed, we can `import` the **CellTypesCache module**.\n",
    "\n",
    "The CellTypesCache that we're importing provides tools to allow us to get information from the cell types database. We're giving it a **manifest** filename as well. CellTypesCache will create this manifest file, which contains metadata about the cache. If you want, you can look in the cell_types folder in your code directory and take a look at the file.\n",
    "\n",
    "><b>Task</b>: Run the cell below. If you'd like a slightly faster way of working with the raw data files required by this notebook, you can also **uncomment** (remove the `#`) from the `!git` and `%cd` lines below to clone the data from Github."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Import the \"Cell Types Cache\" from the AllenSDK core package\n",
    "from allensdk.core.cell_types_cache import CellTypesCache\n",
    "\n",
    "#If you're short on time, clone the repository to get the raw data\n",
    "#!git clone https://github.com/ajuavinett/CellTypesLesson.git\n",
    "#%cd CellTypesLesson\n",
    "\n",
    "#Initialize the cache as 'ctc' (cell types cache)\n",
    "ctc = CellTypesCache(manifest_file='cell_types/manifest.json')\n",
    "\n",
    "print('CellTypesCache imported.')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<a id=\"import\"></a>\n",
    "\n",
    "## Step 2. Import Cell Types data\n",
    "**Note**:  *Steps 2-4 are an optional demonstration of the same cell from the [web-based activity](https://sites.google.com/ucsd.edu/neuroedu/cell-types?authuser=0). These steps can be skipped if you'd like to dive into analyzing populations of cells below.*\n",
    "\n",
    "Now that we have the module that we need, let's import a raw sweep of the data. The cell below will grab the data for the same experiment you just looked at on the website. This data is in the form of a [**Neuroscience Without Borders** (NWB)](https://www.nwb.org/) file.\n",
    "\n",
    "><b>Task</b>: Find the cell specimen ID for the first cell you looked at on the Allen Institute website (hint: it's in the URL), and assign this to <code>cell_id</code> below by replacing the <code>...</code>. Run the cell.\n",
    "\n",
    "This might take a minute or two. You should wait until the circle in the upper right is <i>not</i> filled (Jupyter Notebook) or you have a green checkmark (Colab) to continue."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Enter your cell_id below\n",
    "cell_id = ...\n",
    "\n",
    "# Get the electrophysiology (ephys) data for that cell\n",
    "data = ctc.get_ephys_data(cell_id)\n",
    "print('Data retrieved')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Thankfully, our NWB file has some built-in **methods** to enable us to pull out a recording sweep. We can access methods of objects like our `data` object by adding a period, and then the method. That's what we're doing below, with `data.get_sweep()`.\n",
    "\n",
    "\n",
    "><b>Task:</b> Choose your favorite sweep below. (<u>Hint</u>: go back to the website to see what the sweep numbers are.)\n",
    "    \n",
    "<i>Note</i>: You may get an `H5pyDeprecationWarning`, but don't worry about it. This is out of our control. :)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Assign your favorite sweep number to a variable \"sweep_number\" below.\n",
    "sweep_number = ...\n",
    "\n",
    "sweep_data = data.get_sweep(sweep_number) \n",
    "print('Sweep obtained')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<a id=\"plotsweep\"></a>\n",
    "## Step 3. Plot a raw sweep of data\n",
    "Now that you've pulled down some data, chosen a cell, and chosen a sweep number, let's plot that data.\n",
    "\n",
    "><b>Task:</b> Run the cell below to get the stimulus and recorded response information from the dataset."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": false
   },
   "outputs": [],
   "source": [
    "# Get the stimulus trace (in amps) and convert to pA\n",
    "stim_current = sweep_data['stimulus'] * 1e12\n",
    "\n",
    "# Get the voltage trace (in volts) and convert to mV\n",
    "response_voltage = sweep_data['response'] * 1e3\n",
    "\n",
    "# Get the sampling rate and can create a time axis for our data\n",
    "sampling_rate = sweep_data['sampling_rate'] # in Hz\n",
    "timestamps = (np.arange(0, len(response_voltage)) * (1.0 / sampling_rate))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "><b>Task</b>: In the cell below, use the <code>plt.plot(x,y)</code> to plot our voltage trace.\n",
    ">- You will need to give it two arguments, which are variables we created above: <code>timestamps</code> (x axis) and <code>response_voltage</code>(y).\n",
    "- Without changing the limits on the x-axis, you won't be able to see individual action potentials.\n",
    "- Modify the x-axis using <code>plt.xlim([min,max])</code> to specify the limits (replace <code>min</code> and <code>max</code> with numbers that make sense for this x-axis.\n",
    "- If you'd like to plot the current that was injected into the cell, you can plot <code>stim_current</code> instead of <code>response_voltage.</code>"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Plot the raw recording here\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<a id=\"morphology\"></a>\n",
    "## Step 4. Plot the morphology of the cell\n",
    "The Cell Types Database also contains **3D reconstructions** of neuronal morphologies. Here, we'll plot the reconstruction of our cell's morphology.\n",
    "\n",
    "*Note*: It will take several minutes to run the cell below, possibly longer over a slow internet connection. Running this cell is optional and can be skipped if necessary."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Import necessary toolbox\n",
    "from allensdk.core.swc import Marker\n",
    "\n",
    "# Download and open morphology and marker files\n",
    "morphology = ctc.get_reconstruction(cell_id) \n",
    "markers = ctc.get_reconstruction_markers(cell_id) \n",
    "\n",
    "# Set up our plot\n",
    "fig, axes = plt.subplots(1, 2, sharey=True, sharex=True, figsize=(10,10))\n",
    "axes[0].set_aspect('equal')\n",
    "axes[1].set_aspect('equal')\n",
    "\n",
    "# Make a line drawing of x-y and y-z views\n",
    "for n in morphology.compartment_list:\n",
    "    for c in morphology.children_of(n):\n",
    "        axes[0].plot([n['x'], c['x']], [n['y'], c['y']], color='black')\n",
    "        axes[1].plot([n['z'], c['z']], [n['y'], c['y']], color='black')\n",
    "\n",
    "# cut dendrite markers\n",
    "dm = [ m for m in markers if m['name'] == Marker.CUT_DENDRITE ]\n",
    "axes[0].scatter([m['x'] for m in dm], [m['y'] for m in dm], color='#3333ff')\n",
    "axes[1].scatter([m['z'] for m in dm], [m['y'] for m in dm], color='#3333ff')\n",
    "\n",
    "# no reconstruction markers\n",
    "nm = [ m for m in markers if m['name'] == Marker.NO_RECONSTRUCTION ]\n",
    "axes[0].scatter([m['x'] for m in nm], [m['y'] for m in nm], color='#333333')\n",
    "axes[1].scatter([m['z'] for m in nm], [m['y'] for m in nm], color='#333333')\n",
    "axes[0].set_ylabel('y')\n",
    "axes[0].set_xlabel('x')\n",
    "axes[1].set_xlabel('z')\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<a id=\"metrics\"></a>\n",
    "## Step 5. Analyze pre-computed features\n",
    "\n",
    "The Cell Types Database contains a set of features that have already been computed, which could serve as good starting points for analysis. We can query the database to get these features. Below, we'll use the Pandas package that we imported above to create a **[dataframe](https://pandas.pydata.org/pandas-docs/stable/user_guide/dsintro.html#dataframe)** of our data.\n",
    "\n",
    "><b>Task</b>: Run the cell below. It'll take ~10 seconds. It will print a list of all of the feature available, as well as produce a dataframe, which looks something like an Excel spreadsheet. You can scroll to the right to see many of the different features available in this dataset."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Download all electrophysiology features for all cells\n",
    "ephys_features = ctc.get_ephys_features()\n",
    "dataframe = pd.DataFrame(ephys_features).set_index('specimen_id')\n",
    "\n",
    "print('Ephys features available for %d cells:' % len(dataframe))\n",
    "dataframe.head() # Just show the first 5 rows (the head) of our dataframe "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "As you can see in the dataframe above, there are many **pre-computed features** available in this dataset. [Here's a glossary](https://docs.google.com/document/d/1YGLwkMTebwrXd_1E817LFbztMjSTCWh83Mlp3_3ZMEo/edit?usp=sharing), in case you're curious.\n",
    "\n",
    "![](https://github.com/ajuavinett/CellTypesLesson/blob/master/docs/ap_features.png?raw=true)\n",
    "Image from the <a href=\"http://help.brain-map.org/download/attachments/8323525/CellTypes_Ephys_Overview.pdf\">Allen Institute Cell Types Database Technical Whitepaper.</a>\n",
    "<br><br>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let's first look at the speed of the trough, and the ratio between the upstroke and downstroke of the action potential:\n",
    "- **Action potential fast trough** (<code>fast_trough_v_long_square</code>): Minimum value of the membrane potential in the interval lasting 5 ms after the peak.\n",
    "- **Upstroke/downstroke ratio** (<code>upstroke_downstroke_ratio_long_square</code>)</b>: The ratio between the absolute values of the action potential peak upstroke and the action potential peak downstroke.</div>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Above, we created a pandas dataframe of all of these features. Here, we'll assign the columns we're interested to two different **variables**, so that they will contain all of the datapoints we're interested in. We can access different columns of the dataframe by using the syntax `dataframe['column_of_interest']`. The columns of interest here are `fast_trough_v_long_square` and `upstroke_downstroke_ratio_long_square`].\n",
    "\n",
    "><b>Task:</b> Edit and run the cell below to store these columns into our two new variables."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "fast_trough = ...\n",
    "upstroke_downstroke = ..."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "><b>Task:</b> Create a scatterplot that plots the fast trough (x axis) versus the upstroke-downstroke ratio (y axis). Label your axes accordingly using `plt.xlabel()` and `plt.ylabel()`.\n",
    "    \n",
    "If you need help, see the [`plt.scatter()` documentation](https://matplotlib.org/api/_as_gen/matplotlib.pyplot.scatter.html). "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Your scatterplot here\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "It looks like there may be roughly two clusters in the data above. Maybe they relate to whether the cells are presumably excitatory (spiny) cells or inhibitory (aspiny) cells. Let's query the API and split up the two sets to see.\n",
    "\n",
    "><b>Task:</b> The cell below will dig up the dendrite type of these cells and add that to our dataframe. Then, it'll create our same scatterplot, where each dot is colored by dendrite type. All you need to do is run the cell!"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Get information about our cells' dendrites\n",
    "cells = ctc.get_cells()\n",
    "full_dataframe = dataframe.join(pd.DataFrame(cells).set_index('id'))\n",
    "\n",
    "# Create a dataframe for spiny cells, and a dataframe for aspiny cells\n",
    "spiny_df = full_dataframe[full_dataframe['dendrite_type'] == 'spiny']\n",
    "aspiny_df = full_dataframe[full_dataframe['dendrite_type'] == 'aspiny']\n",
    "\n",
    "# Create our plot! Calling scatter twice like this will draw both of these on the same plot.\n",
    "plt.scatter(spiny_df['fast_trough_v_long_square'],spiny_df['upstroke_downstroke_ratio_long_square'])\n",
    "plt.scatter(aspiny_df['fast_trough_v_long_square'],aspiny_df['upstroke_downstroke_ratio_long_square'])\n",
    "\n",
    "plt.ylabel('upstroke-downstroke ratio')\n",
    "plt.xlabel('fast trough depth (mV)')\n",
    "plt.legend(['Spiny','Aspiny'])\n",
    "    \n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Looks like these two clusters do partially relate to the dendritic type. Cells with spiny dendrites (which are typically excitatory cells) have a big ratio of upstroke:downstroke, and a more shallow trough (less negative). Cells with aspiny dendrites (typically inhibitory cells) are a little bit more varied. But </i>only</i> aspiny cells have a low upstroke:downstroke ratio and a deeper trough (more negative)."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<a id=\"waveforms\"></a>\n",
    "\n",
    "## Step 6. Compare waveforms\n",
    "Let's take a closer look at the action potentials of these cells to see what these features actually mean for the action potential waveform by choosing one of the cells with the highest upstroke:downstroke ratio. Our first line of code, where it says `dataframe.sort_values()` is the code that will arrange our dataframe by the **upstroke_downstroke_ratio_long_square** column.\n",
    "\n",
    "This first time around, we'll organize it so that the highest ratio is at the top (`ascending=False`). This is an example of a **boolean** in Python. You can change this to say `ascending=True` if you want to sort with *lowest* ratio at the top.\n",
    "\n",
    "><b>Task</b>: Run the cell below."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Sort the dataframe and reassign\n",
    "sorted_dataframe = dataframe.sort_values('upstroke_downstroke_ratio_long_square',ascending=False)\n",
    "\n",
    "# Assign one of the top cells in our dataframe (default = 2) and the ratio to different variables\n",
    "specimen_id = sorted_dataframe.index[2]\n",
    "ratio = sorted_dataframe.iloc[2]['upstroke_downstroke_ratio_long_square']\n",
    "\n",
    "# Print our results so that we can see them\n",
    "print('Specimen ID: ' + str(specimen_id) + ' with upstroke-downstroke ratio: ' + str(ratio))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now we can take a closer look at the action potential for that cell by grabbing a raw sweep of recording from it, just like we did above.\n",
    "\n",
    "><b>Task:</b> Run the cell below. This may take a minute or so. *Note*: You may receive a \"H5pyDeprecationWarning,\" but you can ignore this!"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": true
   },
   "outputs": [],
   "source": [
    "# Get the data for our specimen\n",
    "upstroke_data = ctc.get_ephys_data(specimen_id)\n",
    "\n",
    "# Get one sweep for our specimen (I've already handselected a gorgeous one for you, 45)\n",
    "upstroke_sweep = upstroke_data.get_sweep(45) \n",
    "\n",
    "# Get the current & voltage traces\n",
    "current = upstroke_sweep['stimulus'] * 1e12 # in A, converted to pA\n",
    "voltage = upstroke_sweep['response'] * 1e3 # converted to mV!\n",
    "\n",
    "# Get the time stamps for our voltage trace\n",
    "timestamps = (np.arange(0, len(voltage)) * (1.0 / upstroke_sweep['sampling_rate']))\n",
    "\n",
    "print('Sweep obtained')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "><b>Task:</b> Plot the sweep we obtained above. <i>Hint</i>: You'll want to use `plt.plot(x,y)` where `x` is the `timestamps` and `y` is the `voltage`. Be sure to give your plot accurate labels as well."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Plot the new sweep here\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "><b>Task</b> Generate a similar plot for a cell with a <b>low</b> upstroke ratio. Similiar to above, zoom in on the x axis so that you can actually see the shape of the action potential waveform.\n",
    "\n",
    "><b>Hint</b>: You only need to change <i>one</i> value in all of the code in this step in order to make this change. How did we arrange our dataframe at first?"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "As you'll hopefully see, even that one feature, upstroke:downstroke ratio, means the shape of the action potential is dramatically different. The other feature we looked at above, size of the trough, is highly correlated with upstroke:downstroke. You can see that by comparing the two cells here. Cells with high upstroke:downstroke tend to have less negative troughs (undershoots) after the action potential."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<a id=\"compare\"></a>\n",
    "\n",
    "## Step 7. Compare cell types\n",
    "Let's get out of the action potential weeds a bit. What if we want to know a big picture thing, such as are *human cells different than mouse cells?* Or *how are excitatory cells different from inhibitory cells?*  To ask these questions, we can pull out the data for two different cell types, defined by their species, dendrite type, or transgenic line.\n",
    "\n",
    "**About Transgenic Cre Lines.** The Allen Institute for Brain Science uses transgenic mouse lines that have Cre-expressing cells to mark specific types of cells in the brain. This technology is called the **Cre-Lox system**, and is a common way in neuroscience (and some other fields) to target cells based on their expression of specific genetic promotors. For more information about Cre/Lox technology, see [this website](https://old.abmgood.com/marketing/knowledge_base/Cre-Lox_Recombination.php). Information about the different Cre lines that are available can be found in [this glossary](https://docs.google.com/document/d/1ZMMZgc7cS5BHhoWNqzjw95BdxOuj5wrYl9I7PV2HeUI/edit?usp=sharing) or on the [Allen Institute's website](http://connectivity.brain-map.org/transgenic).\n",
    "\n",
    "**For this final step, it's up to you to choose which cell types to compare.** You'll also decide which pre-computed feature to compare between these cell types.\n",
    "\n",
    "- If you'd like to compare cells from different **species**, the column name is `species`.\n",
    "- If you'd like to compare **spiny vs. aspiny cells**, the column name is `dendrite_type`.\n",
    "- If you'd like to compare two **transgenic lines** (mouse cells only), the column name is `transgenic_line`. What if we want to know whether different genetically-identified cells have different intrinsic physiology?\n",
    "\n",
    "\n",
    ">**Task**: Assign `column_name` below to the name of your column to see the unique values in that column. Make sure your column is a **string**, in other words, it should be in single  quotes."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Define your column name below\n",
    "column_name = ...\n",
    "\n",
    "print(full_dataframe[column_name].unique())"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "> Using the possible values in your column, create two separate dataframes by **subsetting** the dataframe below. Assign `celltype_1` and `celltype_2` to the names of your cell types, for example 'spiny' and 'aspiny'. Make sure your cell type names are in quotes (they should be strings) and *exactly* match what is found in the dataframe."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Define your cell type variables below\n",
    "celltype_1 = ...\n",
    "celltype_2 = ...\n",
    "\n",
    "celltype_1_df = full_dataframe[full_dataframe[column_name] == celltype_1]\n",
    "celltype_2_df = full_dataframe[full_dataframe[column_name] == celltype_2]\n",
    "\n",
    "print(\"Type 1 # Cells: %d\" % len(celltype_1_df))\n",
    "print(\"Type 2 # Cells: %d\" % len(celltype_2_df))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let's start by plotting a distribution of the recorded resting membrane potential (`vrest`) for one cell type versuss the other cell type.\n",
    "\n",
    ">**Task**: Run the cell below to plot a histogram to compare one pre-computedd feature of your choice between your two cell types.\n",
    "\n",
    "- Note that the distribution is normalized by the total count (`density=True`), since there may be very different numbers of cells for your two cell types. You can set `density` to false to plot the raw numbers of cells.\n",
    "- You can also specify the number of bins with `bins= < #bins > `.\n",
    "- Look through the [`plt.hist()` documentation](https://matplotlib.org/3.3.1/api/_as_gen/matplotlib.pyplot.hist.html) for more information."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "plt.figure()\n",
    "\n",
    "# Change your feature below.\n",
    "feature = 'vrest'\n",
    "\n",
    "# Plot the histogram, with density = True \n",
    "plt.hist([celltype_1_df[feature],celltype_2_df[feature]],density = True)\n",
    "\n",
    "# Change the x label below:\n",
    "plt.xlabel('< feature name >')\n",
    "plt.ylabel('Normalized Number of Cells')\n",
    "plt.legend(['Cell Type 1','Cell Type 2'])\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "> <b>Task</b>: Choose a different feature to compare between your cell types, and rerun the plot above. Use the documentation below to get the exact name of the feature (in parentheses), and change the x label axis so that we know what you're plotting.  Right click to save your image when you're done.\n",
    "\n",
    "Here are a few additional pre-computed features you might consider comparing (you can find a complete glossary [here](https://docs.google.com/document/d/1YGLwkMTebwrXd_1E817LFbztMjSTCWh83Mlp3_3ZMEo/edit?usp=sharing)):\n",
    "\n",
    "- <b>Tau (<code>tau</code>)</b>: time constant of the membrane in milliseconds\n",
    "- <b>Adapation ratio (<code>adaptation</code>)</b>: The rate at which firing speeds up or slows down during a stimulus<br>\n",
    "- <b>Average ISI (<code>avg_isi</code>)</b>: The mean value of all interspike intervals in a sweep<br>\n",
    "- **Slope of f/I curve** (<code>f_i_curve_slope</code>)</b>: slope of the curve between firing rate (f) and current injected<br>\n",
    "- **Input Resistance** (<code>input_resistance_mohm</code>)</b>: The input resistance of the cell, in megaohms.<br>\n",
    "- **Voltage of after-hyperpolarization** (<code>trough_v_short_square</code>)</b>: minimum value of the membrane potential during the after-hyperpolarization"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "><b>Task</b>: It's more common to plot summary statistics like a mean or median, so let's compare our two cell types with a boxplot. To do so, we can use `plt.boxplot()` ([Documentation here](https://matplotlib.org/3.3.2/api/_as_gen/matplotlib.pyplot.boxplot.html)). The code below is already set up for you -- just run it and edit your labels as necessary. Right click to save the plot when you're done."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Boxplot creation lines below\n",
    "plt.boxplot([celltype_1_df[feature],celltype_2_df[feature]])\n",
    "plt.ylabel('< your label here > ') # y-axis label\n",
    "plt.xticks([1, 2], ['Cell Type 1','Cell Type 2'])\n",
    "\n",
    "# Plot title -- be sure to update!\n",
    "plt.title('< your title here > ')\n",
    "\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "That's it for today -- great work!"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from IPython.display import HTML\n",
    "print('Nice work!')\n",
    "HTML('<img src=\"https://media.giphy.com/media/xUOwGhOrYP0jP6iAy4/giphy.gif\">')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "-----------\n",
    "\n",
    "## Technical notes & credits\n",
    "\n",
    "This notebook demonstrates most of the features of the AllenSDK that help manipulate data in the Cell Types Database.  The main entry point will be through the `CellTypesCache` class. `CellTypesCache` is responsible for downloading Cell Types Database data to a standard directory structure on your hard drive.  If you use this class, you will not have to keep track of where your data lives, other than a root directory.\n",
    "\n",
    "Much more information can be found in the <a href=\"http://help.brain-map.org/download/attachments/8323525/CellTypes_Ephys_Overview.pdf\">Allen Brain Atlas whitepaper</a> as well as in their <a href=\"http://alleninstitute.github.io/AllenSDK/cell_types.html\">GitHub documentation</a>.\n",
    "\n",
    "This file modified from <a href='https://alleninstitute.github.io/AllenSDK/_static/examples/nb/cell_types.html'>this</a> notebook.\n",
    "\n",
    "In case you're curious, <a href=\"https://pandas.pydata.org/pandas-docs/stable/reference/api/pandas.Series.plot.html \">here's documentation</a> for plotting pandas series (which we do quite a bit above). You can always Google questions you have!)"
   ]
  }
 ],
 "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.7.3"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 1
}