MatNWB Optical Physiology Tutorial
Introduction
In this tutorial, we will create fake data for a hypothetical optical physiology experiment with a freely moving animal. The types of data we will convert are:
- Acquired two-photon images
- Image segmentation (ROIs)
- Fluorescence and dF/F response
It is recommended to first work through the Introduction to MatNWB tutorial, which demonstrates installing MatNWB and creating an NWB file with subject information, animal position, and trials, as well as writing and reading NWB files in MATLAB. Set up the NWB file
An NWB file represents a single session of an experiment. Each file must have a session_description, identifier, and session start time. Create a new NWBFile object with those and additional metadata. For all MatNWB functions, we use the Matlab method of entering keyword argument pairs, where arguments are entered as name followed by value.
'session_description', 'mouse in open exploration',...
'identifier', 'Mouse5_Day3', ...
'session_start_time', datetime(2018, 4, 25, 2, 30, 3, 'TimeZone', 'local'), ...
'timestamps_reference_time', datetime(2018, 4, 25, 3, 0, 45, 'TimeZone', 'local'), ...
'general_experimenter', 'LastName, FirstName', ... % optional
'general_session_id', 'session_1234', ... % optional
'general_institution', 'University of My Institution', ... % optional
'general_related_publications', 'DOI:10.1016/j.neuron.2016.12.011'); % optional
Optical Physiology
Optical physiology results are written in four steps:
- Create imaging plane
- Acquired two-photon images
- Image segmentation
- Fluorescence and dF/F responses
Imaging Plane
First, you must create an ImagingPlane object, which will hold information about the area and method used to collect the optical imaging data. This requires creation of a Device object for the microscope and an OpticalChannel object. Then you can create an ImagingPlane. optical_channel = types.core.OpticalChannel( ...
'description', 'description', ...
'emission_lambda', 500.);
device = types.core.Device();
nwb.general_devices.set('Device', device);
imaging_plane_name = 'imaging_plane';
imaging_plane = types.core.ImagingPlane( ...
'optical_channel', optical_channel, ...
'description', 'a very interesting part of the brain', ...
'device', types.untyped.SoftLink(device), ...
'excitation_lambda', 600., ...
'location', 'my favorite brain location');
nwb.general_optophysiology.set(imaging_plane_name, imaging_plane);
Storing One-Photon Data
Now that we have our ImagingPlane, we can create a OnePhotonSeries object to store raw one-photon imaging data. Here, we have two options. The first option is to supply the raw image data using the data argument. The second option is to provide a path to the image files. These two options have trade-offs, so it is worth considering how you want to store this data. Direct Storage
% using internal data. this data will be stored inside the NWB file
InternalOnePhoton = types.core.OnePhotonSeries( ...
'data', ones(100, 100, 1000), ...
'imaging_plane', types.untyped.SoftLink(imaging_plane), ...
'starting_time_rate', 1.0, ...
'data_unit', 'normalized amplitude' ...
nwb.acquisition.set('1pInternal', InternalOnePhoton);
External Linkage
% using external data. only the file paths will be stored inside the NWB file
ExternalOnePhoton = types.core.OnePhotonSeries( ...
'dimension', [100, 100], ...
'external_file', 'images.tiff', ...
'external_file_starting_frame', 0, ...
'imaging_plane', types.untyped.SoftLink(imaging_plane), ...
'format', 'external', ...
'starting_time', 0.0, ...
'starting_time_rate', 1.0 ...
nwb.acquisition.set('1pExternal', ExternalOnePhoton);
Storing Two-Photon Data
Direct Storage
InternalTwoPhoton = types.core.TwoPhotonSeries( ...
'imaging_plane', types.untyped.SoftLink(imaging_plane), ...
'starting_time', 0.0, ...
'starting_time_rate', 3.0, ...
'data', ones(200, 100, 1000), ...
nwb.acquisition.set('2pInternal', InternalTwoPhoton);
External Linkage
ExternalTwoPhoton = types.core.TwoPhotonSeries( ...
'external_file', 'images.tiff', ...
'imaging_plane', types.untyped.SoftLink(imaging_plane), ...
'external_file_starting_frame', 0, ...
'starting_time_rate', 3.0, ...
'starting_time', 0.0, ...
'data', NaN(2, 2, 2), ...
nwb.acquisition.set('2pExternal', ExternalTwoPhoton);
Plane Segmentation
Regions of interest (ROIs)
Add ROIs using an image mask. An image mask is an array that is the same size as a single frame of the TwoPhotonSeries, and indicates where a single region of interest is. This image mask may be boolean or continuous between 0 and 1. % generate fake image_mask data
imaging_shape = [100, 100];
image_mask = zeros(y, x, n_rois);
image_mask(start(1):start(1)+10, start(2):start(2)+10, 1) = 1;
% add data to NWB structures
plane_segmentation = types.core.PlaneSegmentation( ...
'colnames', {'image_mask'}, ...
'description', 'output from segmenting my favorite imaging plane', ...
'id', types.hdmf_common.ElementIdentifiers('data', int64(0:19)), ...
'imaging_plane', types.untyped.SoftLink(imaging_plane), ...
'image_mask', types.hdmf_common.VectorData('data', image_mask, 'description', 'image masks') ...
img_seg = types.core.ImageSegmentation();
img_seg.planesegmentation.set('PlaneSegmentation', plane_segmentation)
ans =
Set with properties:
PlaneSegmentation: [types.core.PlaneSegmentation]
ophys_module = types.core.ProcessingModule( ...
'description', 'contains optical physiology data')
ophys_module =
ProcessingModule with properties:
description: 'contains optical physiology data'
dynamictable: [0×1 types.untyped.Set]
nwbdatainterface: [0×1 types.untyped.Set]
ophys_module.nwbdatainterface.set('ImageSegmentation', img_seg);
nwb.processing.set('ophys', ophys_module);
Storing fluorescence of ROIs over time
Now that ROIs are stored, you can store fluorescence dF/F data for these regions of interest. This type of data is stored using the RoiResponseSeries class. You will not need to instantiate this class directly to create objects of this type, but it is worth noting that this is the class you will work with after you read data back in. First, create a data interface to store this data in
roi_table_region = types.hdmf_common.DynamicTableRegion( ...
'table', types.untyped.ObjectView(plane_segmentation), ...
'description', 'all_rois', ...
roi_response_series = types.core.RoiResponseSeries( ...
'rois', roi_table_region, ...
'data', NaN(n_rois, 100), ...
'data_unit', 'lumens', ...
'starting_time_rate', 3.0, ...
fluorescence = types.core.Fluorescence();
fluorescence.roiresponseseries.set('RoiResponseSeries', roi_response_series);
ophys_module.nwbdatainterface.set('Fluorescence', fluorescence);
nwb.processing.set('ophys', ophys_module);
Write the NWB file
nwbExport(nwb, 'ophys_tutorial.nwb');
Reading the NWB file
read_nwb = nwbRead('ophys_tutorial.nwb', 'ignorecache');
Data arrays are read passively from the file. Calling TimeSeries.data does not read the data values, but presents an HDF5 object that can be indexed to read data.
read_nwb.processing.get('ophys').nwbdatainterface.get('Fluorescence')...
.roiresponseseries.get('RoiResponseSeries').data
ans =
DataStub with properties:
filename: 'ophys_tutorial.nwb'
path: '/processing/ophys/Fluorescence/RoiResponseSeries/data'
dims: [20 100]
ndims: 2
dataType: 'double'
This allows you to conveniently work with datasets that are too large to fit in RAM all at once. Access the data in the matrix using the load method.
load with no input arguments reads the entire dataset:
read_nwb.processing.get('ophys').nwbdatainterface.get('Fluorescence'). ...
roiresponseseries.get('RoiResponseSeries').data.load
NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN
NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN
NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN
NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN
NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN
NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN
NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN
NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN
NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN
NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN
If all you need is a section of the data, you can read only that section by indexing the DataStub object like a normal array in MATLAB. This will just read the selected region from disk into RAM. This technique is particularly useful if you are dealing with a large dataset that is too big to fit entirely into your available RAM.
read_nwb.processing.get('ophys'). ...
nwbdatainterface.get('Fluorescence'). ...
roiresponseseries.get('RoiResponseSeries'). ...
data(1:5, 1:10)
NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN
NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN
NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN
NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN
NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN
% read back the image masks
read_nwb.processing.get('ophys'). ...
nwbdatainterface.get('ImageSegmentation'). ...
planesegmentation.get('PlaneSegmentation'). ...
Learn more!
See the API documentation to learn what data types are available.
Other MatNWB tutorials
Python tutorials
See our tutorials for more details about your data type:
Check out other tutorials that teach advanced NWB topics: