Getting data from a session

Getting data from a session#

We’re going to examine the data available for a single session. We load this using pynwb. This loads all of the data available for this session.

nwb_path = r'/data/222426_2016-02-04_10-25-24_nwb_2026-08-19_17-50-08/222426_2016-02-04_10-25-24_nwb_2026-08-19_17-50-08.nwb.zarr'

nwbfile = pynwb.read_nwb(nwb_path)

---------------------------------------------------------------------------
ValueError                                Traceback (most recent call last)
Cell In[2], line 3
      1 nwb_path = r'/data/222426_2016-02-04_10-25-24_nwb_2026-08-19_17-50-08/222426_2016-02-04_10-25-24_nwb_2026-08-19_17-50-08.nwb.zarr'
----> 3 nwbfile = pynwb.read_nwb(nwb_path)

File /opt/envs/ctlut/lib/python3.10/site-packages/hdmf/utils.py:596, in docval.<locals>.dec.<locals>.func_call(*args, **kwargs)
    594 def func_call(*args, **kwargs):
    595     pargs = _check_args(args, kwargs)
--> 596     return func(**pargs)

File /opt/envs/ctlut/lib/python3.10/site-packages/pynwb/__init__.py:622, in read_nwb(**kwargs)
    620         return NWBZarrIO.read_nwb(path=path)
    621     else:
--> 622         raise ValueError(
    623             f"Unable to read file: '{path}'. The file is not recognized as "
    624             "either a valid HDF5 or Zarr NWB file. Please ensure the file exists and contains valid NWB data."
    625         )
    626 except ImportError:
    627     raise ValueError(
    628         f"Unable to read file: '{path}'. The file is not recognized as an HDF5 NWB file. "
    629         "If you are trying to read a Zarr file, please install hdmf-zarr using: pip install hdmf-zarr"
    630     )

ValueError: Unable to read file: '/data/222426_2016-02-04_10-25-24_nwb_2026-08-19_17-50-08/222426_2016-02-04_10-25-24_nwb_2026-08-19_17-50-08.nwb.zarr'. The file is not recognized as either a valid HDF5 or Zarr NWB file. Please ensure the file exists and contains valid NWB data.

nwbfile

Let’s explore:

Maximum projection#

This is the projection of the full motion corrected movie. It shows all of the cells imaged during the session.

max_projection = nwbfile.processing['ophys'].data_interfaces['SummaryImages'].images['maximum_intensity_projection'][:]

fig = plt.figure(figsize=(6,6))
plt.imshow(max_projection, cmap='gray')
plt.axis('off')

ROI Masks#

ROIs are all of the segmented masks for cell bodies identified in this session. These are stored in the PlaneSegmentation table using a sparse array.

seg = nwbfile.processing["ophys"]["ImageSegmentation"]["PlaneSegmentation"].to_dataframe()
seg.head()

Let’s look at how this is represented. Each mask is a list of (x,y,weight) for only the pixels where the ROI mask is located.

seg['pixel_mask'][517473350]

Plot the masks for all the ROIs.

rois = np.zeros((512,512))
for index,row in seg.iterrows():
    for x, y, weight in row.pixel_mask:
        rois[int(x), int(y)] = weight

plt.imshow(rois)

Knowing the location of a neuron is valuable if you want to examine the spatial relationships between neurons. For instance, you can calculate the center of an ROI (take the mean of the x and y pixel locations) and use that to measure the distance between two neurons.

Fluorescence and DF/F traces#

The NWB file contains a number of traces reflecting the processing that is done to the extracted fluorescence before we analyze it. The fluorescence traces are the mean fluorescence of all the pixels contained within a ROI mask. In addition to the raw fluorescence, there are also neuropil corrected traces, demixed traces, and DF/F traces.

The signal we are most interested in is the DFF - the change in fluorescence normalized by the baseline fluorescence. The baseline fluorescence was computed as the median fluorescence in a 180s window centered on each time point. The result is the dff trace:

dff_series = nwbfile.processing["ophys"]["DfOverF"]["DfOverF"]
dff = dff_series.data[:].T  # Transpose to get (n_cells, n_timepoints)
ts = dff_series.timestamps[:]

fig = plt.figure(figsize=(8,3))
plt.plot(ts, dff[122,:], color='gray')
plt.xlabel("Time (s)")
plt.xlim(1900,2200)
plt.ylabel("DFF")

Extracted events#

In addition to these traces, we also provide events extracted from the DF/F traces using the L0 method developed by Sean Jewell and Daniella Witten.

# Get DfOverF events in a RoiResponseSeries
dff_events_series = nwbfile.processing["ophys"]["DfOverF"]["DfOverFEvents"]
dff_events = dff_events_series.data[:].T  # Transpose to get (n_cells, n_timepoints)
ts = dff_events_series.timestamps[:]

fig = plt.figure(figsize=(8,3))
plt.plot(ts, dff[122,:], color='gray')
plt.plot(ts, 2*dff_events[122,:]+5, color='black')
plt.xlabel("Time (s)")
plt.xlim(1900,2200)
plt.ylabel("DFF")

Stimulus epochs#

Several stimuli are shown during each imaging session, interleaved with each other. The stimulus epoch table provides information of these interleaved stimulus epochs, revealing when each epoch starts and ends. .

stim_epoch = nwbfile.intervals['epochs'].to_dataframe()
stim_epoch

Column	Description
start_time	The time at the start of the epoch
stop_time	The time at the end of the epoch
stimulus_type	The name of the stimulus for the epoch

Let’s plot the DFF traces of a number of cells and overlay stimulus epochs.

fig = plt.figure(figsize=(14,8))

#here we plot the first 50 neurons in the session
for i in range(50):
    plt.plot(ts, dff[i,:]+(i*2), color='gray')

#here we shade the plot when each stimulus is presented
colors = ['blue','orange','green','red']
for c, stim_name in enumerate(stim_epoch.stimulus_type.unique()):
    stim = stim_epoch[stim_epoch.stimulus_type==stim_name]
    for j in range(len(stim)):
        plt.axvspan(xmin=stim.start_time.iloc[j], xmax=stim.stop_time.iloc[j], color=colors[c], alpha=0.1)

Running speed#

The running speed of the animal on the rotating disk during the entire session.

running_speed_series = nwbfile.processing["behavior"]["BehavioralTimeSeries"]['running_speed']
dxcm = running_speed_series.data[:]
running_ts = running_speed_series.timestamps[:]

Plot the running speed.

plt.plot(running_ts, dxcm)
plt.ylabel("Running speed (cm/s)", fontsize=18)
plt.xlabel("Time (s)", fontsize=18)

Add the running speed to the neural activity and stimulus epoch figure we made above

fig = plt.figure(figsize=(14,8))

#here we plot the first 50 neurons in the session
for i in range(50):
    plt.plot(ts, dff[i,:]+(i*2), color='gray')

#here we shade the plot when each stimulus is presented
colors = ['blue','orange','green','red']
for c, stim_name in enumerate(stim_epoch.stimulus_type.unique()):
    stim = stim_epoch[stim_epoch.stimulus_type==stim_name]
    for j in range(len(stim)):
        plt.axvspan(xmin=stim.start_time.iloc[j], xmax=stim.stop_time.iloc[j], color=colors[c], alpha=0.1)

#here we add the running speed (scaled and offset)
plt.plot(running_ts, (0.2*dxcm)-20)

Stimulus Table and Template#

Each stimulus that is shown has a stimulus table that details what each trial is and when it is presented. Additionally, the natural scenes, natural movies, and locally sparse noise stimuli have a stimulus template that shows the exact image that is presented to the mouse. We detail how to access and use these items in Visual stimuli.

Cell ids and indices#

Each neuron in the dataset has a unique id. These IDs are stored in the PlaneSegmentation table we looked at before.

# Cell IDs are stored in the PlaneSegmentation table
plane_seg = nwbfile.processing["ophys"]["ImageSegmentation"]["PlaneSegmentation"]
cell_ids = plane_seg.id[:]
cell_ids

Within each individual session, a cell id is associated with an index. This index maps into the dff or event arrays. Pick one cell id from the list above and find the index for that neuron.

target_cell_id = cell_ids[0]  # Use first cell as example
cell_index = np.where(cell_ids == target_cell_id)[0]
print(f"Cell ID {target_cell_id} is at index {cell_index}")

During data processing, we matched identified ROIs across each of the sessions within experiment containers. Approximately one third of the neurons in the dataset were matched across all three sessions, one third were matched in two of the three session, and one third were only found in one session. When neurons are matched across sessions, that neuron will have the same cell id in all said sessions. This is explored in Cross session data.

How come we don’t always match ROIs across all three session for all neurons?

There are a few factors that could explain why we don’t always match ROIs across all sessions that include biological, experimental, and analytical reasons. Biologically, a neuron must be active within a session to be identifiable during segmentation. For various reasons, a neuron might not be active during some sessions while it is active during others. Experimentally, there are challenges to returning to the precise same field of view. Being at a slightly different depth, or having just a bit of tilt in the imaging plane, might result in some neurons that were in view during one session not being in view during another. Analytically, the method for identifying ROIs as well as for matching ROIs from multiple sessions can make mistakes.