Skip to content

Note

Click here to download the full example code

Fit injected current

Data for this notebook is a patch clamp experiment with a mouse V1 neuron, from the Allen Brain Atlas

Allen Brain Atlas view of the data we will analyze.

Learning objectives

  • Learn how to explore spiking data and do basic analyses using pynapple
  • Learn how to structure data for NeMoS
  • Learn how to fit a basic Generalized Linear Model using NeMoS
  • Learn how to retrieve the parameters and predictions from a fit GLM for intrepetation.
# Import everything
import jax
import os
import matplotlib.pyplot as plt
import nemos as nmo
import nemos.glm
import numpy as np
import pynapple as nap
import workshop_utils

# configure plots some
plt.style.use(workshop_utils.STYLE_FILE)

Data Streaming

path = workshop_utils.fetch_data("allen_478498617.nwb")

Pynapple

Data structures and preparation

data = nap.load_file(path)
print(data)

Annotated view of the data we will analyze.

  • stimulus: Tsd containing injected current, in Amperes, sampled at 20k Hz.
  • response: Tsd containing the neuron's intracellular voltage, sampled at 20k Hz.
  • units: Tsgroup, dictionary of neurons, holding each neuron's spike timestamps.
  • epochs: IntervalSet, dictionary with start and end times of different intervals, defining the experimental structure, specifying when each stimulation protocol began and ended.
epochs = data["epochs"]
epochs.keys()
  • Noise 1: epochs of random noise
noise_interval = epochs["Noise 1"]
noise_interval
  • Let's focus on the first epoch.
noise_interval = noise_interval[0]
noise_interval
  • current : Tsd (TimeSeriesData) : time index + data
# convert current from Ampere to pico-amperes, to match the Allen Institute figures and
# move the values to a more reasonable range.
current = data["stimulus"] * 1e12
current
  • restrict : restricts a time series object to a set of time intervals delimited by an IntervalSet object
current = current.restrict(noise_interval)
current
  • TsGroup : a custom dictionary holding multiple Ts (timeseries) objects with potentially different time indices.
spikes = data["units"]
spikes

We can index into the TsGroup to see the timestamps for this neuron's spikes:

spikes[0]

Let's restrict to the same epoch noise_interval:

spikes = spikes.restrict(noise_interval)
print(spikes)
spikes[0]

Let's visualize the data from this trial:

fig, ax = plt.subplots(1, 1, figsize=(8, 2))
ax.plot(current, "grey")
ax.plot(spikes.to_tsd([-5]), "|", color="k", ms = 10)
ax.set_ylabel("Current (pA)")
ax.set_xlabel("Time (s)")

Basic analyses

The Generalized Linear Model gives a predicted firing rate. First we can use pynapple to visualize this firing rate for a single trial.

  • count : count the number of events within bin_size
# bin size in seconds
bin_size = 0.001
# Get spikes for neuron 0
count = spikes[0].count(bin_size)
count

Let's convert the spike counts to firing rate :

  • smooth : convolve with a Gaussian kernel
# the inputs to this function are the standard deviation of the gaussian in seconds and
# the full width of the window, in standard deviations. So std=.05 and size_factor=20
# gives a total filter size of 0.05 sec * 20 = 1 sec.
firing_rate = count.smooth(std=.05, size_factor=20)
# convert from spikes per bin to spikes per second (Hz)
firing_rate = firing_rate / bin_size

# we're hiding the details of the plotting function for the purposes of this
# tutorial, but you can find it in the associated github repo if you're
# interested:
# https://github.com/flatironinstitute/ccn-workshop-fens-2024/blob/main/src/workshop_utils/plotting.py
workshop_utils.plotting.current_injection_plot(current, spikes, firing_rate)

What is the relationship between the current and the spiking activity? compute_1d_tuning_curves : compute the firing rate as a function of a 1-dimensional feature.

tuning_curve = nap.compute_1d_tuning_curves(spikes, current, nb_bins=15)
tuning_curve

Let's plot the tuning curve of the neuron.

workshop_utils.plotting.tuning_curve_plot(tuning_curve)

NeMoS

Preparing data

Get data from pynapple to NeMoS-ready format:

  • predictors and spikes must have same number of time points
# enter code here
  • predictors must be 2d, spikes 1d
# enter code here

Fitting the model

  • GLM objects need regularizers and observation models
# enter code here
  • call fit and retrieve parameters
# enter code here
  • generate and examine model predictions.
# enter code here
  • what do we see?
# enter code here
  • examine tuning curve — what do we see?
# enter code here

Extending the model

  • choose a length of time over which the neuron integrates the input current
# enter code here
  • define a basis object
# enter code here
  • create the design matrix
  • examine the features it contains
# enter code here
  • create and fit the GLM
  • examine the parameters
# enter code here
  • compare the predicted firing rate to the data and the old model
  • what do we see?
# enter code here
  • examine the predicted average firing rate and tuning curve
  • what do we see?
# enter code here
  • use log-likelihood to compare models
# enter code here

Finishing up

  • what if you want to compare models across datasets?
# enter code here
  • what about spiking?
# enter code here

Further Exercises

  • what else can we do?

Data citation

The data used in this tutorial is from the Allen Brain Map, with the following citation:

Contributors: Agata Budzillo, Bosiljka Tasic, Brian R. Lee, Fahimeh Baftizadeh, Gabe Murphy, Hongkui Zeng, Jim Berg, Nathan Gouwens, Rachel Dalley, Staci A. Sorensen, Tim Jarsky, Uygar Sümbül Zizhen Yao

Dataset: Allen Institute for Brain Science (2020). Allen Cell Types Database -- Mouse Patch-seq [dataset]. Available from brain-map.org/explore/classes/multimodal-characterization.

Primary publication: Gouwens, N.W., Sorensen, S.A., et al. (2020). Integrated morphoelectric and transcriptomic classification of cortical GABAergic cells. Cell, 183(4), 935-953.E19. https://doi.org/10.1016/j.cell.2020.09.057

Patch-seq protocol: Lee, B. R., Budzillo, A., et al. (2021). Scaled, high fidelity electrophysiological, morphological, and transcriptomic cell characterization. eLife, 2021;10:e65482. https://doi.org/10.7554/eLife.65482

Mouse VISp L2/3 glutamatergic neurons: Berg, J., Sorensen, S. A., Miller, J., Ting, J., et al. (2021) Human neocortical expansion involves glutamatergic neuron diversification. Nature, 598(7879):151-158. doi: 10.1038/s41586-021-03813-8

Total running time of the script: ( 0 minutes 0.000 seconds)

Download Python source code: 05_current_injection_users.py

Download Jupyter notebook: 05_current_injection_users.ipynb

Gallery generated by mkdocs-gallery