NRSC7610_ReceptiveField.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Computing Visual Receptive Fields \n",
    "In this notebook, we'll estimate the spatiotemporal receptive field of a number of cells in the dorsal lateral geniculate nucleus (dLGN) of the mouse. To estimate these receptive fields we will use the spike times evoked from single neurons in response to a **binary dense noise** stimulus. The stimulus is a movie, which looks something like this on each frame:\n",
    "<br>\n",
    "![stimulus frame](./res/binary_image.png)\n",
    "<br>\n",
    "We'll get spike times, from extracellular recordings with [Neuropixels](https://www.neuropixels.org/) and assigned to single neurons with [spyking-circus](https://spyking-circus.readthedocs.io/).\n",
    "<br>\n",
    "<br>\n",
    "These data are from [Denman et al., Journal of Neuroscience, 2017](http://www.jneurosci.org/cgi/pmidlookup?view=long&pmid=27986926) \n",
    "<br>\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##### Some of the code is pre-written, and you will simply execute it. In other cases, you will be prompted to write some code lines to advance towards this goal. **There are going to be 8 such prompts, with ** plus some bonus challenges if you want to keep going.\n",
    "\n",
    "## Imports of packages we will need for this notebook. \n",
    "All are standard"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import os,sys,glob, h5py\n",
    "import matplotlib.pyplot as plt\n",
    "from utils import cleanAxes"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Load the data:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "nwb_data = h5py.File('./res/M192079_noRawData.nwb')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Q1. What type of data is** ```nwb_data```**? write some code that will give the answer in the cell below. assign your answer to the variable** ```nwb_data_type```"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "To explore nwb_data, we need to know what is in it. An HDF5 file is of the [Heirarchical Data Format](https://support.hdfgroup.org/HDF5/whatishdf5.html), which means that it has data in a heirarchy. We can think of it kind of like a nested python ```dictionary```, with data grouped into different sections based on keys:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "nwb_data.keys()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The most important for our purposes will be ```processing``` and ```stimulus```. ```processing``` contains the results of spike sorring, in which the raw data are \"processed\" into spike times that are assigned to particular single neruons\n",
    "<br>\n",
    "For example, we can look in processing to see what is in there:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "nwb_data['processing'].keys()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "and look inside further, to see the elements in ```processing```:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "nwb_data['processing']['LGN'].keys()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "and look inside...further:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "nwb_data['processing']['LGN']['UnitTimes'].keys()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "These keys are the IDs of the neurons identified in the recording. So, there is a neuron #10, and a neuron #107, and so on\n",
    "<br>\n",
    "**Q2. What are the elements within one of these neurons? assign your answer to the variable** ```keys```"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "keys = "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Q3. Using the appropriate key from** ```keys``` **, get the spike times from cell #87. assign your answer to the variable** ```spike_times```"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "spike_times ="
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Q4: convert** ```spike_times``` **to a numpy ndarray**"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We've now got some spike times from an example neuron. Let's get the other critical data, the stimulus information. <br>\n",
    "There were several stimuli in this experiment. We will use these two:\n",
    "- ```binary_green```\n",
    "- ```binary_uv```\n",
    "<br>Both of these are \"binary\" spatiotemporal white noise, meaning there are two luminances (hence \"binary\"), a bright and dark (or, as plotted above, black and white)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "nwb_data['stimulus']['presentation'].keys()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Within each of these, there are two important entries:\n",
    "- ```data```\n",
    "- ```timestamps```"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print(nwb_data['stimulus']['presentation']['binary_green'].keys())\n",
    "print(nwb_data['stimulus']['presentation']['binary_uv'].keys())"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The ```timestamps``` entry is a 1D array, which contains the time of each frame\n",
    "<br>\n",
    "The ```data``` entry is a 3D array, with each frame (64 x 64 pixels) matching one of the timestamps \n",
    "<br>**Q5: How many frames are in the** ```binary_uv``` **stimulus? assign your answer to the variable** ```number_of_frames``` "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "stimulus = np.array(nwb_data['stimulus']['presentation']['binary_green']['data']).T\n",
    "stimulus_times = np.array(nwb_data['stimulus']['presentation']['binary_green']['timestamps'])+0.04 # this an adjustment for the hardware used, can ignore"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We can now estimate the receptive field of our example neuron, using the spike-triggered average method. The key concept is to make an average of the stimulus frames (the 64 X 64) images in the relevant ```nwb_data['stimulus']['presentation'][stimulus]['data']``` entries. The average we want to make is of these frames - but not just any, and certainly not all, frames.  We want to know what is the average frame on the screen **immediately preceding** a spike from our example cell. To get these frames, we need to find the frame _times_ that immediately precede a spike. \n",
    "<br>\n",
    "<br> One strategy would be to iterate over the spike times and find the closest stimulus frame. But notice that the stimulus times only cover a subset of the spike times:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print('Spike times go from '+str(spike_times[0])+' to '+str(spike_times[-1])+' seconds')\n",
    "print('Stimulus times go from '+str(stimulus_times[0])+' to '+str(stimulus_times[-1])+' seconds')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "For this analysis, it would be ineffecient to iterate over the spike times that happen during other stimuli, since we only care about the responses to the ```binary_green``` stimulus. \n",
    "<br>\n",
    "**Q6: Find the subset of spike times that occured during the stimulus. Assign these to a numpy ndarray called** ```stimulus_spike_times```"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now we can interate over these and find the frame that was being presented some amount of time before each spike. Let's first make an educated guess at what amount of time to use: 90 milliseconds, to give the stimulus some time to be tranduced by the retina and the response to ramp up in the LGN."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "time_before_spike = .090 # in seconds"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Q7: For each spike time, find the frame _number_ (this is an index, in the range** ```0``` **to** ```number_of_frames``` **. Make a** ```list``` **or** ```array``` **of these indices and assign it to** ```stimulus_frame_indices``` <br>_hint: you'll need to compare each adjusted spike time to_ ```stimulus_times``` "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now we know which frames were on the screen 90 milliseconds before each spike. We need to make an average of those frames.\n",
    "<br>\n",
    "**Q8 Make an average of the frames from the indexes we found. Call it** ```spatial_receptive_field```"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Finally, we will plot the spatial receptive field:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "plt.imshow(spatial_receptive_field,clim=(-0.2,0.2),cmap=plt.cm.Greys)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Copy (Ctrl+right click on the image) and paste the image you get into an email to Dan(daniel.denman@cuanschutz.edu) and Alon(alon.poleg-polsky@cuanschutz.edu). You're done!\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<br>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\n",
    "**Bonus 1 Calculate the spatial receptive field at many different times before the spike** ```time_before_spike``` **, from 0 to 250 milliseconds**"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Bonus 2 Calculate the temporal receptive field at the center of the receptive field. This will be the value, across the different** ```time_before_spike```**s you used in Bonus 1, at a given pixel defined as the center of the receptive field.**"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "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.6.7"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}