suite2p is a blind source separation toolbox for cell detection and signal extraction.
Here we illustrate how to use suite2p to detect cell locations, and then use FISSA to remove neuropil signals from the ROI signals.
The suite2p parts of this tutorial are based on their Jupyter notebook example.
Note that the below results are not representative of either suite2p or FISSA performance, as we are using a very small example dataset.
Reference: Pachitariu, M., Stringer, C., Dipoppa, M., Schröder, S., Rossi, L. F., Dalgleish, H., Carandini, M. & Harris, K. D. (2017). Suite2p: beyond 10,000 neurons with standard two-photon microscopy. bioRxiv: 061507; doi: 10.1101/061507.
In [1]:
# FISSA toolbox
import fissa
# suite2p toolbox
import suite2p.run_s2p
# numpy toolbox
import numpy as np
# Plotting toolbox, with notebook embedding options
import holoviews as hv
%load_ext holoviews.ipython
%output widgets='embed'
In [2]:
# Set your options for running
ops = suite2p.run_s2p.default_ops() # populates ops with the default options
# Provide an h5 path in 'h5py' or a tiff path in 'data_path'
# db overwrites any ops (allows for experiment specific settings)
db = {
'h5py': [], # a single h5 file path
'h5py_key': 'data',
'look_one_level_down': False, # whether to look in ALL subfolders when searching for tiffs
'data_path': ['exampleData/20150529'], # a list of folders with tiffs
# (or folder of folders with tiffs if look_one_level_down is True,
# or subfolders is not empty)
'save_path0': './', # save path
'subfolders': [], # choose subfolders of 'data_path' to look in (optional)
'fast_disk': './', # string which specifies where the binary file will be stored (should be an SSD)
'reg_tif': True, # save the motion corrected tiffs
'tau': 0.7, # timescale of gcamp6f
'fs': 1, # sampling rate
'spatial_scale': 4
}
# Run one experiment
opsEnd = suite2p.run_s2p.run_s2p(ops=ops, db=db)
In [3]:
# Extract the motion corrected tiffs (make sure that the reg_tif option is set to true, see above)
images = './suite2p/plane0/reg_tif'
# Load the detected regions of interest
stat = np.load('./suite2p/plane0/stat.npy', allow_pickle=True) # cell stats
ops = np.load('./suite2p/plane0/ops.npy', allow_pickle=True).item()
iscell = np.load('./suite2p/plane0/iscell.npy', allow_pickle=True)[:, 0]
# Get image size
Lx = ops['Lx']
Ly = ops['Ly']
# Get the cell ids
ncells = len(stat)
cell_ids = np.arange(ncells) # assign each cell an ID, starting from 0.
cell_ids = cell_ids[iscell==1] # only take the ROIs that are actually cells.
num_rois = len(cell_ids)
# Generate ROI masks in a format usable by FISSA (in this case, a list of masks)
rois = [np.zeros((Ly, Lx), dtype=bool) for n in range(num_rois)]
for i, n in enumerate(cell_ids):
# i is the position in cell_ids, and n is the actual cell number
ypix = stat[n]['ypix'][~stat[n]['overlap']]
xpix = stat[n]['xpix'][~stat[n]['overlap']]
rois[i][ypix, xpix] = 1
In [4]:
output_folder = 'fissa_suite2p_example'
exp = fissa.Experiment(images, [rois[:ncells]], output_folder)
exp.separate()
In [5]:
%%opts Curve {+axiswise}
def plot_cell_regions(roi_polys, plot_neuropil=False):
'''
Plot a single cell region, using holoviews.
'''
out = hv.Overlay()
if plot_neuropil:
# Plot the neuropil as well as the ROI
n_region = len(roi_polys)
else:
# Just plot the ROI, not the neuropil
n_region = 1
for i_region in range(n_region):
for part in roi_polys[i_region]:
x = part[:, 1]
y = part[:, 0]
out *= hv.Curve(zip(x, y)).opts(color='w')
return out
i_trial = 0
# Generate plots for all detected regions
region_plots = {
i_cell: plot_cell_regions(exp.roi_polys[i_cell][i_trial])
for i_cell in range(exp.nCell)
}
# Generate plots for raw extracts and neuropil removed
traces_plots = {
i_cell: hv.Curve(exp.raw[i_cell][i_trial][0, :], label='suite2p') *
hv.Curve(exp.result[i_cell][i_trial][0, :], label='FISSA')
for i_cell in range(exp.nCell)
}
# Generate average image
avg_img = hv.Raster(exp.means[i_trial])
# Generate overlay plot showing each cell location
cell_locs = hv.Overlay()
for c in range(exp.nCell):
roi_poly = exp.roi_polys[c][i_trial][0][0]
x = roi_poly[:, 1]
y = roi_poly[:, 0]
cell_locs *= hv.Curve(zip(x, y))
# Render holoviews
avg_img * cell_locs * hv.HoloMap(region_plots, kdims=['Cell']) + hv.HoloMap(traces_plots, kdims=['Cell'])
Out[5]:
Note that with the above settings for suite2p it seems to have detected more small local axon signals, instead of cells. This can possibly be improved with manual curation and suite2p setting changes, but as noted above these results should not be seen as indicative for either suite2p or FISSA due to the small dataset size.
Also note that the above Suite2P traces are done without suite2p's own neuropil removal algorithm.