diff --git a/caiman/source_extraction/cnmf/estimates.py b/caiman/source_extraction/cnmf/estimates.py index 7a26c8fa0..80c1688bf 100644 --- a/caiman/source_extraction/cnmf/estimates.py +++ b/caiman/source_extraction/cnmf/estimates.py @@ -381,6 +381,55 @@ def nb_view_components(self, Yr=None, img=None, idx=None, thr=thr, denoised_color=denoised_color, cmap=cmap) return self + def hv_view_components(self, Yr=None, img=None, idx=None, + denoised_color=None, cmap='viridis'): + """view spatial and temporal components interactively in a notebook + + Args: + Yr : np.ndarray + movie in format pixels (d) x frames (T) + + img : np.ndarray + background image for contour plotting. Default is the mean + image of all spatial components (d1 x d2) + + idx : list + list of components to be plotted + + denoised_color: string or None + color name (e.g. 'red') or hex color code (e.g. '#F0027F') + + cmap: string + name of colormap (e.g. 'viridis') used to plot image_neurons + """ + if 'csc_matrix' not in str(type(self.A)): + self.A = scipy.sparse.csc_matrix(self.A) + + plt.ion() + nr, T = self.C.shape + if self.R is None: + self.R = self.YrA + if self.R.shape != [nr, T]: + if self.YrA is None: + self.compute_residuals(Yr) + else: + self.R = self.YrA + + if img is None: + img = np.reshape(np.array(self.A.mean(axis=1)), self.dims, order='F') + + if idx is None: + hv_plot = caiman.utils.visualization.hv_view_patches( + Yr, self.A, self.C, self.b, self.f, self.dims[0], self.dims[1], + YrA=self.R, image_neurons=img, denoised_color=denoised_color, + cmap=cmap) + else: + hv_plot = caiman.utils.visualization.hv_view_patches( + Yr, self.A.tocsc()[:, idx], self.C[idx], self.b, self.f, + self.dims[0], self.dims[1], YrA=self.R[idx], image_neurons=img, + denoised_color=denoised_color, cmap=cmap) + return hv_plot + def nb_view_components_3d(self, Yr=None, image_type='mean', dims=None, max_projection=False, axis=0, denoised_color=None, cmap='jet', thr=0.9): diff --git a/caiman/utils/visualization.py b/caiman/utils/visualization.py index 8af20bb31..d90b13d3a 100644 --- a/caiman/utils/visualization.py +++ b/caiman/utils/visualization.py @@ -30,6 +30,8 @@ from tempfile import NamedTemporaryFile from typing import Dict from warnings import warn +import holoviews as hv +import functools as fct from ..base.rois import com from ..summary_images import local_correlations @@ -248,6 +250,76 @@ def nb_view_patches(Yr, A, C, b, f, d1, d2, YrA=None, image_neurons=None, thr=0. return Y_r +def hv_view_patches(Yr, A, C, b, f, d1, d2, YrA=None, image_neurons=None, denoised_color=None, cmap='viridis'): + """ + Interactive plotting utility for ipython notebook + + Args: + Yr: np.ndarray + movie + + A,C,b,f: np.ndarrays + outputs of matrix factorization algorithm + + d1,d2: floats + dimensions of movie (x and y) + + YrA: np.ndarray + ROI filtered residual as it is given from update_temporal_components + If not given, then it is computed (K x T) + + image_neurons: np.ndarray + image to be overlaid to neurons (for instance the average) + + denoised_color: string or None + color name (e.g. 'red') or hex color code (e.g. '#F0027F') + + cmap: string + name of colormap (e.g. 'viridis') used to plot image_neurons + """ + + nr, T = C.shape + nA2 = np.ravel(np.power(A, 2).sum(0)) if type( + A) == np.ndarray else np.ravel(A.power(2).sum(0)) + b = np.squeeze(b) + f = np.squeeze(f) + if YrA is None: + Y_r = np.array( + spdiags(old_div(1, nA2), 0, nr, nr) * + (A.T * np.matrix(Yr) - + (A.T * np.matrix(b[:, np.newaxis])) * np.matrix(f[np.newaxis]) - + A.T.dot(A) * np.matrix(C)) + C) + else: + Y_r = C + YrA + if image_neurons is None: + image_neurons = A.mean(1).reshape((d1, d2), order='F') + smp = cm.ScalarMappable(cmap=cmap) + im_rgb = smp.to_rgba(image_neurons)[:, :, :3] + im_hsv = mpl.colors.rgb_to_hsv(im_rgb) + + def norm(a, rg=(0, 1)): + a_norm = (a - a.min()) / (a.max() - a.min()) + return a_norm * (rg[1] - rg[0]) + rg[0] + + Ad = np.asarray(A.todense()).reshape((d1, d2, -1), order='F') + + def plot_unit(uid, scl): + trace = ( + hv.Curve(Y_r[uid, :], kdims='time').opts(framewise=True) + * (hv.Curve(C[uid, :], kdims='time') + .opts(color=denoised_color, framewise=True)) + ).opts(aspect=2, frame_height=200) + A_scl = norm(Ad[:, :, uid], (scl, 1)) + im_hsv_scl = im_hsv.copy() + im_hsv_scl[:, :, 2] = im_hsv[:, :, 2] * A_scl + im_u = (hv.HSV(im_hsv_scl, kdims=['height', 'width']) + .opts(aspect='equal', frame_height=200)) + return im_u + trace + + return (hv.DynamicMap(plot_unit, kdims=['unit_id', 'scale']) + .redim.range(unit_id=(0, nr), scale=(0.0, 1.0))) + + def get_contours(A, dims, thr=0.9, thr_method='nrg', swap_dim=False): """Gets contour of spatial components and returns their coordinates @@ -992,6 +1064,36 @@ def GetBox(centers, R, dims): ax.axis('off') pl.subplots_adjust(0, 0, 1, 1, .06, .06) + +def nb_inspect_correlation_pnr(corr, pnr): + """ + inspect correlation and pnr images to infer the min_corr, min_pnr + + Args: + corr: ndarray + correlation image created with caiman.summary_images.correlation_pnr + + pnr: ndarray + peak-to-noise image created with caiman.summary_images.correlation_pnr + """ + hv_corr = hv.Image(corr, vdims='corr', label='correlation') + hv_pnr = hv.Image(pnr, vdims='pnr', label='pnr') + + def hist(im, rx, ry): + obj = im.select(x=rx, y=ry) if rx and ry else im + return hv.operation.histogram(obj) + + str_corr = (hv.streams.RangeXY(source=hv_corr) + .rename(x_range='rx', y_range='ry')) + str_pnr = (hv.streams.RangeXY(source=hv_pnr) + .rename(x_range='rx', y_range='ry')) + hist_corr = hv.DynamicMap( + fct.partial(hist, im=hv_corr), streams=[str_corr]) + hist_pnr = hv.DynamicMap( + fct.partial(hist, im=hv_pnr), streams=[str_pnr]) + return (hv_corr << hist_corr) + (hv_pnr << hist_pnr) + + def inspect_correlation_pnr(correlation_image_pnr, pnr_image): """ inspect correlation and pnr images to infer the min_corr, min_pnr diff --git a/demos/notebooks/demo_pipeline_cnmfE.ipynb b/demos/notebooks/demo_pipeline_cnmfE.ipynb index 0313fd0c0..258d81478 100644 --- a/demos/notebooks/demo_pipeline_cnmfE.ipynb +++ b/demos/notebooks/demo_pipeline_cnmfE.ipynb @@ -43,7 +43,7 @@ "import caiman as cm\n", "from caiman.source_extraction import cnmf\n", "from caiman.utils.utils import download_demo\n", - "from caiman.utils.visualization import inspect_correlation_pnr\n", + "from caiman.utils.visualization import inspect_correlation_pnr, nb_inspect_correlation_pnr\n", "from caiman.motion_correction import MotionCorrect\n", "from caiman.source_extraction.cnmf import params as params\n", "from caiman.utils.visualization import plot_contours, nb_view_patches, nb_plot_contour\n", @@ -54,7 +54,9 @@ "except:\n", " pass\n", "import bokeh.plotting as bpl\n", - "bpl.output_notebook()" + "import holoviews as hv\n", + "bpl.output_notebook()\n", + "hv.notebook_extension('bokeh')" ] }, { @@ -306,7 +308,7 @@ "# compute some summary images (correlation and peak to noise)\n", "cn_filter, pnr = cm.summary_images.correlation_pnr(images[::1], gSig=gSig[0], swap_dim=False) # change swap dim if output looks weird, it is a problem with tiffile\n", "# inspect the summary images and set the parameters\n", - "inspect_correlation_pnr(cn_filter, pnr)" + "nb_inspect_correlation_pnr(cn_filter, pnr)" ] }, { @@ -440,7 +442,7 @@ "outputs": [], "source": [ "# accepted components\n", - "cnm.estimates.nb_view_components(img=cn_filter, idx=cnm.estimates.idx_components,\n", + "cnm.estimates.hv_view_components(img=cn_filter, idx=cnm.estimates.idx_components,\n", " denoised_color='red', cmap='gray')" ] }, @@ -453,7 +455,7 @@ "outputs": [], "source": [ "# rejected components\n", - "cnm.estimates.nb_view_components(img=cn_filter, idx=cnm.estimates.idx_components_bad,\n", + "cnm.estimates.hv_view_components(img=cn_filter, idx=cnm.estimates.idx_components_bad,\n", " denoised_color='red', cmap='gray')" ] },