{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "# NWB Tutorial - Extracellular Electrophysiology" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Introduction\n", "In this tutorial, we will create an NWB file for a hypothetical extracellular electrophysiology experiment with a freely moving animal. The types of data we will convert are:\n", "- Subject information (species, strain, age, etc.) \n", "- Animal position\n", "- Trials\n", "- Raw data\n", "- LFP\n", "- Spike times" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Installing PyNWB\n", "First, install PyNWB using pip or conda. You will need Python 3.5+ installed.\n", "- `pip install pynwb`\n", "- `conda install -c conda-forge pynwb`" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Set up the NWB file\n", "An NWB file represents a single session of an experiment. Each file must have a session description, identifier, and session start time. Importantly, the session start time is the reference time for all timestamps in the file. For example, an event with a timestamp of 0 in the file means the event occurred exactly at the session start time. \n", "\n", "Create a new `NWBFile` object with those and additional metadata. For all constructors in PyNWB, we recommend using keyword arguments for clarity." ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "from pynwb import NWBFile\n", "from datetime import datetime\n", "from dateutil import tz\n", "\n", "session_start_time = datetime(2018, 4, 25, 2, 30, 3, tzinfo=tz.gettz('US/Pacific'))\n", "\n", "nwbfile = NWBFile(\n", " session_description='Mouse exploring an open field',\n", " identifier='Mouse5_Day3',\n", " session_start_time=session_start_time,\n", " session_id='session_1234', # optional\n", " experimenter='My Name', # optional\n", " lab='My Lab Name', # optional\n", " institution='University of My Institution', # optional\n", " related_publications='DOI:10.1016/j.neuron.2016.12.011' # optional\n", ")\n", "print(nwbfile)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Subject information\n", "Create a `Subject` object to store information about the experimental subject, such as age, species, genotype, sex, and a description. Then set `nwbfile.subject` to the new `Subject` object.\n", "\n", "\n", "\n", "Each of these fields is free-form text, so any values will be valid, but we recommend these values follow particular conventions to help software tools interpret the data:\n", "- For age, we recommend using the [ISO 8601 Duration format](https://en.wikipedia.org/wiki/ISO_8601#Durations), e.g., \"P90D\" for 90 days old\n", "- For species, we recommend using the formal latin binomal name, e.g., \"*Mus musculus*\", \"*Homo sapiens*\"\n", "- For sex, we recommend using \"F\" (female), \"M\" (male), \"U\" (unknown), and \"O\" (other)" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "from pynwb.file import Subject\n", "\n", "nwbfile.subject = Subject(\n", " subject_id='001',\n", " age='P90D', \n", " description='mouse 5',\n", " species='Mus musculus', \n", " sex='M'\n", ")" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## SpatialSeries and Position\n", "PyNWB contains classes that are specialized for different types of data. To store the spatial position of a subject, we will use the `SpatialSeries` and `Position` classes. \n", "\n", "`SpatialSeries` is a subclass of `TimeSeries`. `TimeSeries` is a common base class for measurements sampled over time, and provides fields for data and time (regularly or irregularly sampled).\n", "\n", "\n", "\n", "Create a `SpatialSeries` object named `'SpatialSeries'` with some fake data." ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "import numpy as np\n", "from pynwb.behavior import SpatialSeries\n", "\n", "# create fake data with shape (50, 2)\n", "# the first dimension should always represent time\n", "position_data = np.array([np.linspace(0, 10, 50),\n", " np.linspace(0, 8, 50)]).T\n", "position_timestamps = np.linspace(0, 50) / 200\n", "\n", "spatial_series_obj = SpatialSeries(\n", " name='SpatialSeries', \n", " description='(x,y) position in open field',\n", " data=position_data,\n", " timestamps=position_timestamps,\n", " reference_frame='(0,0) is bottom left corner'\n", ")" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "You can print the `SpatialSeries` object to view its contents." ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "print(spatial_series_obj)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "To help data analysis and visualization tools know that this `SpatialSeries` object represents the position of the subject, store the `SpatialSeries` object inside of a `Position` object, which can hold one or more `SpatialSeries` objects.\n", "\n", "" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "from pynwb.behavior import Position\n", "\n", "position_obj = Position(spatial_series=spatial_series_obj) # name is set to 'Position' by default" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Behavior Processing Module\n", "\n", "NWB differentiates between raw, *acquired data*, which should never change, and *processed data*, which are the results of preprocessing algorithms and could change. Let's assume that the subject's position was computed from a video tracking algorithm, so it would be classified as processed data. Since processed data can be diverse, NWB allows us to create processing modules, which are like folders, to store related processed data. \n", "\n", "Create a processing module called \"behavior\" for storing behavioral data in the `NWBFile` and add the `Position` object to the processing module." ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "behavior_module = nwbfile.create_processing_module(\n", " name='behavior', \n", " description='processed behavioral data'\n", ")\n", "behavior_module.add(position_obj)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Write to file" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Now, write the NWB file that we have built so far." ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "from pynwb import NWBHDF5IO\n", "\n", "with NWBHDF5IO('ecephys_tutorial.nwb', 'w') as io:\n", " io.write(nwbfile)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Read from an NWB file\n", "\n", "We can now read the file and print it to inspect its contents. \n", "\n", "We can also print the SpatialSeries data that we created by referencing the names of the objects in the hierarchy that contain it. We can access a processing module by indexing `nwbfile.processing` with the name of the processing module, which in our case is \"behavior\". \n", "\n", "Then, we can access the `Position` object inside of the \"behavior\" processing module by indexing it with the name of the `Position` object. The default name of `Position` objects is \"Position\". \n", "\n", "Finally, we can access the `SpatialSeries` object inside of the `Position` object by indexing it with the name of the `SpatialSeries` object, which we named `'SpatialSeries'`." ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "with NWBHDF5IO('ecephys_tutorial.nwb', 'r') as io:\n", " read_nwbfile = io.read()\n", " print(read_nwbfile.processing['behavior'])\n", " print(read_nwbfile.processing['behavior']['Position'])\n", " print(read_nwbfile.processing['behavior']['Position']['SpatialSeries'])" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "We can also use the [HDFView](https://www.hdfgroup.org/downloads/hdfview/) tool to inspect the resulting NWB file.\n", "\n", "" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Trials\n", "\n", "Trials are stored in a `TimeIntervals` object which is a subclass of `DynamicTable`. `DynamicTable` objects are used to store tabular metadata throughout NWB, including for trials, electrodes, and sorted units. They offer flexibility for tabular data by allowing required columns, optional columns, and custom columns not defined in the standard.\n", "\n", "\n", "\n", "The trials DynamicTable can be thought of as a table with this structure:\n", "\n", "\n", "\n", "We can add custom, user-defined columns to the trials table to hold data and metadata specific to this experiment or session. Continue adding to our `NWBFile` by creating a new column for the trials table named `'correct'`, which will be a boolean array." ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "nwbfile.add_trial_column(name='correct', description='whether the trial was correct')\n", "nwbfile.add_trial(start_time=1.0, stop_time=5.0, correct=True)\n", "nwbfile.add_trial(start_time=6.0, stop_time=10.0, correct=False)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "We can view and inspect the trials table in tabular form by converting it to a pandas dataframe." ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "df = nwbfile.trials.to_dataframe()\n", "df" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Extracellular electrophysiology\n", "\n", "## Electrodes table\n", "In order to store extracellular electrophysiology data, you first must create an electrodes table describing the electrodes that generated this data. Extracellular electrodes are stored in an `electrodes` table, which is also a `DynamicTable`. The `electrodes` table has several required fields: x, y, z, impedence, location, filtering, and electrode group.\n", "\n", "\n", "\n", "Since this is a DynamicTable, we can add additional metadata fields. We will be adding a \"label\" column to the table.\n", "Use the following code to add electrodes for an array with 4 shanks and 3 channels per shank." ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "nwbfile.add_electrode_column(name='label', description='label of electrode')\n", "\n", "nshanks = 4;\n", "nchannels_per_shank = 3;\n", "electrode_counter = 0\n", "device = nwbfile.create_device(\n", " name='array',\n", " description='the best array',\n", " manufacturer='Probe Company 9000'\n", ")\n", "for ishank in range(nshanks):\n", " # create an electrode group for this shank\n", " electrode_group = nwbfile.create_electrode_group(\n", " name='shank{}'.format(ishank),\n", " description='electrode group for shank {}'.format(ishank),\n", " device=device,\n", " location='brain area'\n", " )\n", " # add electrodes to the electrode table\n", " for ielec in range(nchannels_per_shank):\n", " nwbfile.add_electrode(\n", " x=5.3, y=1.5, z=8.5, imp=np.nan,\n", " location='unknown', \n", " filtering='unknown',\n", " group=electrode_group,\n", " label='shank{}elec{}'.format(ishank, ielec)\n", " )\n", " electrode_counter += 1" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Like for the trials table, we can view the electrodes table in tabular form by converting it to a pandas dataframe." ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "df = nwbfile.electrodes.to_dataframe()\n", "df" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Links\n", "In the above loop, we created `ElectrodeGroup` objects in the `NWBFile`, and when we added an electrode to the `NWBFile`, we passed in the `ElectrodeGroup` object for the required `'group'` argument. This creates a reference from the `electrodes` table to individual `ElectrodeGroup` objects, one per row (electrode).\n", "\n", "## ElectricalSeries and DynamicTableRegion\n", "Raw voltage data and LFP data are stored in `ElectricalSeries` objects. `ElectricalSeries` is a subclass of `TimeSeries` specialized for voltage data. In order to create our `ElectricalSeries` objects, we will need to reference a set of rows in the `electrodes` table to indicate which electrodes were recorded. We will do this by creating a `DynamicTableRegion`, which is a type of link that allows you to reference specific rows of a `DynamicTable`, such as the `electrodes` table, by row indices.\n", "\n", "Create a `DynamicTableRegion` that references all rows of the `electrodes` table." ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "all_table_region = nwbfile.create_electrode_table_region(\n", " region=list(range(electrode_counter)), # reference row indices 0 to N-1\n", " description='all electrodes'\n", ")" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Raw data\n", "\n", "Now create an `ElectricalSeries` object to hold raw data collected during the experiment, passing in this `DynamicTableRegion` reference to all rows of the `electrodes` table.\n", "\n", "" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "from pynwb.ecephys import ElectricalSeries\n", "\n", "raw_data = np.random.randn(50, 4)\n", "raw_elec_series = ElectricalSeries(\n", " name='ElectricalSeries', \n", " data=raw_data, \n", " electrodes=all_table_region, \n", " starting_time=0.,\n", " rate=20000.\n", ")" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "NWB organizes data into different groups depending on the type of data. Groups can be thought of as folders within the file. Here are some of the groups within an NWB file and the types of data they are intended to store:\n", "- **acquisition**: raw, acquired data that should never change\n", "- **processing**: processed data, typically the results of preprocessing algorithms and could change\n", "- **analysis**: results of data analysis\n", "- **stimuli**: stimuli used in the experiment (e.g., images, videos, light pulses)\n", "\n", "Since this `ElectricalSeries` represents raw data from the data acquisition system, add it to the acquisition group of the NWB file." ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "nwbfile.add_acquisition(raw_elec_series)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## LFP\n", "\n", "Now create an `ElectricalSeries` object to hold LFP data collected during the experiment, again passing in the `DynamicTableRegion` reference to all rows of the `electrodes` table." ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "lfp_data = np.random.randn(50, 4)\n", "lfp_elec_series = ElectricalSeries(\n", " name='ElectricalSeries', \n", " data=lfp_data, \n", " electrodes=all_table_region, \n", " starting_time=0.,\n", " rate=200.\n", ")" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "To help data analysis and visualization tools know that this `ElectricalSeries` object represents LFP data, store the `ElectricalSeries` object inside of an `LFP` object. This is analogous to how we stored the `SpatialSeries` object inside of a `Position` object earlier.\n", "\n", "" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "from pynwb.ecephys import LFP\n", "\n", "lfp = LFP(electrical_series=lfp_elec_series)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Unlike the raw data, which we put into the acquisition group of the NWB file, LFP data is typically considered processed data because the raw data was filtered and downsampled to generate the LFP. \n", "\n", "Create a processing module named `'ecephys'` and add the `LFP` object to it. This is analogous to how we stored the `Position` object in a processing module named `'behavior'` earlier." ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "ecephys_module = nwbfile.create_processing_module(\n", " name='ecephys', \n", " description='processed extracellular electrophysiology data'\n", ")\n", "ecephys_module.add(lfp)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Spike Times\n", "Spike times are stored in the `Units` table, which is another subclass of `DynamicTable`. We can add columns to the `Units` table just like we did for the electrodes and trials tables. \n", "\n", "Generate some random spike data and populate the `Units` table using `nwbfile.add_unit`. Then display the `Units` table as a pandas dataframe." ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "nwbfile.add_unit_column(name='quality', description='sorting quality')\n", "\n", "poisson_lambda = 20\n", "firing_rate = 20\n", "n_units = 10\n", "for n_units_per_shank in range(n_units):\n", " n_spikes = np.random.poisson(lam=poisson_lambda)\n", " spike_times = np.round(np.cumsum(np.random.exponential(1/firing_rate, n_spikes)), 5)\n", " nwbfile.add_unit(spike_times=spike_times, quality='good', waveform_mean=[1., 2., 3., 4., 5.])\n", "\n", "df = nwbfile.units.to_dataframe()\n", "df" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Write the file" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "with NWBHDF5IO('ecephys_tutorial.nwb', 'w') as io:\n", " io.write(nwbfile)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Read the NWB file\n", "\n", "We can access the raw data by indexing `nwbfile.acquisition` with the name of the `ElectricalSeries`, which we named `'ElectricalSeries'`. \n", "\n", "We can also access the LFP data by indexing `nwbfile.processing` with the name of the processing module, \"ecephys\". Then, we can access the `LFP` object inside of the \"ecephys\" processing module by indexing it with the name of the `LFP` object. The default name of `LFP` objects is \"LFP\". \n", "\n", "Finally, we can access the `ElectricalSeries` object inside of the `LFP` object by indexing it with the name of the `ElectricalSeries` object, which we named `'ElectricalSeries'`." ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "with NWBHDF5IO('ecephys_tutorial.nwb', 'r') as io:\n", " read_nwbfile = io.read()\n", " print(read_nwbfile.acquisition['ElectricalSeries'])\n", " print(read_nwbfile.processing['ecephys'])\n", " print(read_nwbfile.processing['ecephys']['LFP'])\n", " print(read_nwbfile.processing['ecephys']['LFP']['ElectricalSeries'])" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Reading NWB data\n", "\n", "Data arrays are read passively from the file. Calling the `data` attribute on a `TimeSeries` such as an `ElectricalSeries` or `SpatialSeries` does not read the data values, but presents an `h5py` object that can be indexed to read data. You can use the `[:]` operator to read the entire data array into memory.\n", "\n", "Load and print all the `data` values of the `ElectricalSeries` object representing the LFP data." ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "with NWBHDF5IO('ecephys_tutorial.nwb', 'r') as io:\n", " read_nwbfile = io.read()\n", " print(read_nwbfile.processing['ecephys']['LFP']['ElectricalSeries'].data[:])" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Accessing data regions\n", "It is often preferable to read only a portion of the data. To do this, index or slice into the `data` attribute just like if you were indexing or slicing a numpy array.\n", "\n", "The following code prints elements 0:10 in the first dimension (time) and 0:3 in the second dimension (electrodes) from the LFP data we have written.\n", "\n", "Accessing data from a `DynamicTable` is similar: `read_nwbfile.units['spike_times'][0]` reads only the spike times from the 0th unit (the value at the 0th row and column named 'spike_times' of the `Units` table)." ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "with NWBHDF5IO('ecephys_tutorial.nwb', 'r') as io:\n", " read_nwbfile = io.read()\n", "\n", " print('section of lfp:')\n", " print(read_nwbfile.processing['ecephys']['LFP']['ElectricalSeries'].data[:10,:3])\n", " print('')\n", " print('spike times from 0th unit:')\n", " print(read_nwbfile.units['spike_times'][0])" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "# Learn more!\n", "\n", "## Python tutorials\n", "### See our tutorials for more details about your data type:\n", "* [Extracellular electrophysiology](https://pynwb.readthedocs.io/en/stable/tutorials/domain/ecephys.html#sphx-glr-tutorials-domain-ecephys-py)\n", "* [Calcium imaging](https://pynwb.readthedocs.io/en/stable/tutorials/domain/ophys.html#sphx-glr-tutorials-domain-ophys-py)\n", "* [Intracellular electrophysiology](https://pynwb.readthedocs.io/en/stable/tutorials/domain/icephys.html#sphx-glr-tutorials-domain-icephys-py)\n", "\n", "### Check out other tutorials that teach advanced NWB topics:\n", "* [Iterative data write](https://pynwb.readthedocs.io/en/stable/tutorials/general/iterative_write.html#sphx-glr-tutorials-general-iterative-write-py)\n", "* [Extensions](https://pynwb.readthedocs.io/en/stable/tutorials/general/extensions.html#sphx-glr-tutorials-general-extensions-py)\n", "* [Advanced HDF5 I/O](https://pynwb.readthedocs.io/en/stable/tutorials/general/advanced_hdf5_io.html#sphx-glr-tutorials-general-advanced-hdf5-io-py)\n", "\n", "\n", "## MATLAB tutorials\n", "* [Extracellular electrophysiology](https://neurodatawithoutborders.github.io/matnwb/tutorials/html/ecephys.html)\n", "* [Calcium imaging](https://neurodatawithoutborders.github.io/matnwb/tutorials/html/ophys.html)\n", "* [Intracellular electrophysiology](https://neurodatawithoutborders.github.io/matnwb/tutorials/html/icephys.html)\n" ] } ], "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.7" } }, "nbformat": 4, "nbformat_minor": 4 }