Note

Click here to download the full example code

# Fit injected current

## 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)
# Set the default precision to float64, which is generally a good idea for
# optimization purposes.
jax.config.update("jax_enable_x64", True)
```

## Data Streaming

- Stream the data. Format is Neurodata Without Borders (NWB) standard

```
path = workshop_utils.data.download_data("allen_478498617.nwb", "https://osf.io/vf2nj/download",
'../data')
```

## Pynapple

### Data structures and preparation

- Open the NWB file with pynapple

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

`units`

: dictionary of neurons, holding each neuron's spike timestamps.`epochs`

: start and end times of different intervals, defining the experimental structure, specifying when each stimulation protocol began and ended.`stimulus`

: injected current, in Amperes, sampled at 20k Hz.`response`

: the neuron's intracellular voltage, sampled at 20k Hz.

```
trial_interval_set = data["epochs"]
# convert current from Ampere to pico-amperes, to match the above visualization
# and move the values to a more reasonable range.
current = data["stimulus"] * 1e12
spikes = data["units"]
```

`trial_interval_set`

: dictionnary of start and end times of different intervals, defining the experimental structure, specifying when each stimulation protocol began and ended.

```
trial_interval_set.keys()
```

`Noise 1`

: epochs of random noise

```
noise_interval = trial_interval_set["Noise 1"]
noise_interval
```

- Let's focus on the first epoch.

```
noise_interval = noise_interval.loc[[0]]
noise_interval
```

`current`

: Tsd (TimeSeriesData) : time index + data

```
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 dictionnary holding multiple`Ts`

(timeseries) objects with potentially different time indices.

```
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
count = spikes.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 and
# the full width of the window, given in bins. So std=50 corresponds to a
# standard deviation of 50*.001=.05 seconds
firing_rate = count.smooth(std=50, size=1000)
# convert from spikes per bin to spikes per second (Hz)
firing_rate = firing_rate / bin_size
```

Note that firing_rate is a `TsdFrame`

!

```
print(type(firing_rate))
```

```
# 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/nemos-workshop-feb-2024/blob/binder/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

- predictors must be 2d, spikes 1d

- predictors and spikes must be jax arrays

### Fitting the model

- GLM objects need regularizers and observation models

- call fit and retrieve parameters

- generate and examine model predictions.

- what do we see?

- examine tuning curve -- what do we see?

### Finishing up

- Finally, let's look at spiking and scoring/metrics

## Further Exercises

- what else can we do?

## 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: 01_current_injection_users.py