Merge pull request #1514 from flatironinstitute/dev

Dev to main for release 1.12.2

Author: pgunn

README.md

@@ -168,7 +168,6 @@ Special thanks to the following people for letting us use their datasets in demo
 
 Also a special thanks to:
 * Eric Thompson, for various strong contributions to code and demos, both before and during his employment at the Flatiron Institute.
-* Cai Changjia, for Volpy
 * Ethan Blackwood, for several contributions in various areas
 
 # License

VERSION

@@ -1 +1 @@
-1.12.1
+1.12.2

bin/caiman_gui.py

@@ -16,6 +16,23 @@
 from caiman.source_extraction.cnmf.cnmf import load_CNMF
 from caiman.source_extraction.cnmf.params import CNMFParams
 
+##############################
+# This is a tool that can help you visualise dumps of the CNMF object,
+# generated as demonstrated in demo_pipeline.
+#
+# If you're already using it, that's fine, but there are better tools
+# out there (such as the Mesmerize-vis package). To use, point it at the
+# hdf5 file that CLI demo or jupyter notebook made, or have your own code
+# call the save method on the CNMF object, and then feed that (and then the
+# mmap file) into this script when it asks, and you'll see your ROIs and
+# be able to experiment with thresholding (and re-save with those different
+# thresholds).
+#
+# mesmerize-vis offers a better approach to these things, and is much more
+# recent and better maintained.
+
+
+
 # Always start by initializing Qt (only once per application)
 app = QtWidgets.QApplication([])
 
@@ -66,7 +83,7 @@ def make_color_img(img, gain=255, min_max=None, out_type=np.uint8):
     M = FileDialog()
     d, f = os.path.split(cnm_obj.mmap_file)
     cnm_obj.mmap_file = M.getOpenFileName(caption='Load memory mapped file',
-                                          directory=d, filter=f + ';;*.mmap')[0]
+                                          filter=f + ';;*.mmap')[0]
 
 if fpath[-3:] == 'nwb':
     mov = caiman.load(cnm_obj.mmap_file,
@@ -405,23 +422,23 @@ def move(event):
 
 
 ## PARAMS
-params = [{'name': 'min_cnn_thr', 'type': 'float', 'value': 0.99, 'limits': (0, 1),'step':0.01},
-            {'name': 'cnn_lowest', 'type': 'float', 'value': 0.1, 'limits': (0, 1),'step':0.01},
-            {'name': 'rval_thr', 'type': 'float', 'value': 0.85, 'limits': (-1, 1),'step':0.01},
-            {'name': 'rval_lowest', 'type': 'float', 'value': -1, 'limits': (-1, 1),'step':0.01},
-            {'name': 'min_SNR', 'type': 'float', 'value': 2, 'limits': (0, 20),'step':0.1},
-            {'name': 'SNR_lowest', 'type': 'float', 'value': 0, 'limits': (0, 20),'step':0.1},
-            {'name': 'RESET', 'type': 'action'},
+params = [{'name': 'min_cnn_thr',       'type': 'float', 'value': 0.99, 'limits': (0, 1),  'step':0.01},
+            {'name': 'cnn_lowest',      'type': 'float', 'value': 0.1,  'limits': (0, 1),  'step':0.01},
+            {'name': 'rval_thr',        'type': 'float', 'value': 0.85, 'limits': (-1, 1), 'step':0.01},
+            {'name': 'rval_lowest',     'type': 'float', 'value': -1,   'limits': (-1, 1), 'step':0.01},
+            {'name': 'min_SNR',         'type': 'float', 'value': 2,    'limits': (0, 20), 'step':0.1},
+            {'name': 'SNR_lowest',      'type': 'float', 'value': 0,    'limits': (0, 20), 'step':0.1},
+            {'name': 'RESET',           'type': 'action'},
             {'name': 'SHOW BACKGROUND', 'type': 'action'},
-            {'name': 'SHOW NEURONS', 'type': 'action'}
+            {'name': 'SHOW NEURONS',    'type': 'action'}
             ]
 
 ## Create tree of Parameter objects
 pars = Parameter.create(name='params', type='group', children=params) 
 
 
 params_action = [{'name': 'Filter components', 'type': 'bool', 'value': True, 'tip': "Filter components"},          
-                 {'name': 'View components', 'type': 'list', 'values': ['All','Accepted',
+                 {'name': 'View components',   'type': 'list', 'limits': ['All','Accepted',
                                                        'Rejected', 'Unassigned'], 'value': 'All'},
                  {'name': 'ADD GROUP', 'type': 'action'},
                  {'name': 'REMOVE GROUP', 'type': 'action'},

caiman/__init__.py

@@ -1,10 +1,11 @@
 #!/usr/bin/env python
 
-import pkg_resources
+import importlib.metadata
+
 from caiman.base.movies import movie, load, load_movie_chain, _load_behavior, play_movie
 from caiman.base.timeseries import concatenate
 from caiman.cluster import start_server, stop_server
 from caiman.mmapping import load_memmap, save_memmap, save_memmap_each, save_memmap_join
 from caiman.summary_images import local_correlations
 
-__version__ = pkg_resources.get_distribution('caiman').version
+__version__ = importlib.metadata.version('caiman')

caiman/base/movies.py

@@ -187,7 +187,7 @@ def bin_median(self, window: int = 10) -> np.ndarray:
                 median image
 
         """
-        T, d1, d2 = np.shape(self)
+        T, d1, d2 = self.shape
         num_windows = int(T // window)
         num_frames = num_windows * window
         return np.nanmedian(np.nanmean(np.reshape(self[:num_frames], (window, num_windows, d1, d2)), axis=0), axis=0)
@@ -207,7 +207,7 @@ def bin_median_3d(self, window=10):
                 median image
 
         """
-        T, d1, d2, d3 = np.shape(self)
+        T, d1, d2, d3 = self.shape
         num_windows = int(T // window)
         num_frames = num_windows * window
         return np.nanmedian(np.nanmean(np.reshape(self[:num_frames], (window, num_windows, d1, d2, d3)), axis=0),
@@ -465,18 +465,18 @@ def computeDFF(self, secsWindow: int = 5, quantilMin: int = 8, method: str = 'on
         if np.min(self) <= 0 and method != 'only_baseline':
             raise ValueError("All pixels must be positive")
 
-        numFrames, linePerFrame, pixPerLine = np.shape(self)
+        numFrames, linePerFrame, pixPerLine = self.shape
         downsampfact = int(secsWindow * self.fr)
         logger.debug(f"Downsample factor: {downsampfact}")
         elm_missing = int(np.ceil(numFrames * 1.0 / downsampfact) * downsampfact - numFrames)
         padbefore = int(np.floor(elm_missing / 2.))
         padafter = int(np.ceil(elm_missing / 2.))
 
-        logger.debug(f'Initial Size Image: {np.shape(self)}')
+        logger.debug(f'Initial Size Image: {self.shape}')
         sys.stdout.flush()
         mov_out = movie(np.pad(self.astype(np.float32), ((padbefore, padafter), (0, 0), (0, 0)), mode='reflect'),
                         **self.__dict__)
-        numFramesNew, linePerFrame, pixPerLine = np.shape(mov_out)
+        numFramesNew, linePerFrame, pixPerLine = mov_out.shape
 
         # compute baseline quickly
         logger.debug("binning data ...")
@@ -540,7 +540,7 @@ def IPCA(self, components: int = 50, batch: int = 1000) -> tuple[np.ndarray, np.
         """
 
         # vectorize the images
-        num_frames, h, w = np.shape(self)
+        num_frames, h, w = self.shape
         frame_size = h * w
         frame_samples = np.reshape(self, (num_frames, frame_size)).T
 
@@ -611,7 +611,7 @@ def IPCA_stICA(self,
         joint_ics = ica.fit_transform(eigenstuff)
 
         # extract the independent frames
-        _, h, w = np.shape(self)
+        _, h, w = self.shape
         frame_size = h * w
         ind_frames = joint_ics[:frame_size, :]
         ind_frames = np.reshape(ind_frames.T, (componentsICA, h, w))
@@ -623,7 +623,7 @@ def IPCA_denoise(self, components: int = 50, batch: int = 1000):
         Create a denoised version of the movie using only the first 'components' components
         """
         _, _, clean_vectors = self.IPCA(components, batch)
-        self = self.__class__(np.reshape(np.float32(clean_vectors.T), np.shape(self)), **self.__dict__)
+        self = self.__class__(np.reshape(np.float32(clean_vectors.T), self.shape), **self.__dict__)
         return self
 
     def local_correlations(self,
@@ -1214,7 +1214,7 @@ def load(file_name: Union[str, list[str]],
 
             if input_arr.ndim == 2:
                 if shape is not None:
-                    _, T = np.shape(input_arr)
+                    _, T = input_arr.shape
                     d1, d2 = shape
                     input_arr = np.transpose(np.reshape(input_arr, (d1, d2, T), order='F'), (2, 0, 1))
                 else:
@@ -1361,13 +1361,13 @@ def load_movie_chain(file_list: list[str],
             if m.ndim == 2:
                 m = m[np.newaxis, :, :]
 
-            _, h, w = np.shape(m)
+            _, h, w = m.shape
             m = m[:, top:h - bottom, left:w - right]
         else:
             if m.ndim == 3:
                 m = m[np.newaxis, :, :, :]
 
-            _, h, w, d = np.shape(m)
+            _, h, w, d = m.shape
             m = m[:, top:h - bottom, left:w - right, z_top:d - z_bottom]
 
         mov.append(m)

caiman/base/rois.py

@@ -123,14 +123,14 @@ def extract_binary_masks_from_structural_channel(Y,
 def mask_to_2d(mask):
     # todo todocument
     if mask.ndim > 2:
-        _, d1, d2 = np.shape(mask)
+        _, d1, d2 = mask.shape
         dims = d1, d2
         return scipy.sparse.coo_matrix(np.reshape(mask[:].transpose([1, 2, 0]), (
             np.prod(dims),
             -1,
         ), order='F'))
     else:
-        dims = np.shape(mask)
+        dims = mask.shape
         return scipy.sparse.coo_matrix(np.reshape(mask, (
             np.prod(dims),
             -1,
@@ -140,7 +140,7 @@ def mask_to_2d(mask):
 def get_distance_from_A(masks_gt, masks_comp, min_dist=10) -> list:
     # todo todocument
 
-    _, d1, d2 = np.shape(masks_gt)
+    _, d1, d2 = masks_gt.shape
     dims = d1, d2
     A_ben = scipy.sparse.csc_matrix(np.reshape(masks_gt[:].transpose([1, 2, 0]), (
         np.prod(dims),
@@ -214,7 +214,7 @@ def nf_match_neurons_in_binary_masks(masks_gt,
     """
     logger = logging.getLogger("caiman")
 
-    _, d1, d2 = np.shape(masks_gt)
+    _, d1, d2 = masks_gt.shape
     dims = d1, d2
 
     # transpose to have a sparse list of components, then reshaping it to have a 1D matrix red in the Fortran style
@@ -244,8 +244,8 @@ def nf_match_neurons_in_binary_masks(masks_gt,
 
     # compute precision and recall
     TP = np.sum(np.array(costs) < thresh_cost) * 1.
-    FN = np.shape(masks_gt)[0] - TP
-    FP = np.shape(masks_comp)[0] - TP
+    FN = masks_gt.shape[0] - TP
+    FP = masks_comp.shape[0] - TP
     TN = 0
 
     performance = dict()
@@ -259,9 +259,9 @@ def nf_match_neurons_in_binary_masks(masks_gt,
     idx_tp_ben = matches[0][idx_tp]    # ground truth
     idx_tp_cnmf = matches[1][idx_tp]   # algorithm - comp
 
-    idx_fn = np.setdiff1d(list(range(np.shape(masks_gt)[0])), matches[0][idx_tp])
+    idx_fn = np.setdiff1d(list(range(masks_gt.shape[0])), matches[0][idx_tp])
 
-    idx_fp = np.setdiff1d(list(range(np.shape(masks_comp)[0])), matches[1][idx_tp])
+    idx_fp = np.setdiff1d(list(range(masks_comp.shape[0])), matches[1][idx_tp])
 
     idx_fp_cnmf = idx_fp
 
@@ -736,8 +736,8 @@ def distance_masks(M_s:list, cm_s: list[list], max_dist: float, enclosed_thr: Op
         test_comp = test_comp.copy()[:, :]
 
         # the number of components for each
-        nb_gt = np.shape(gt_comp)[-1]
-        nb_test = np.shape(test_comp)[-1]
+        nb_gt   = gt_comp.shape[-1]
+        nb_test = test_comp.shape[-1]
         D = np.ones((nb_gt, nb_test))
 
         cmgt_comp = np.array(cmgt_comp)

caiman/behavior/behavior.py

@@ -46,7 +46,7 @@ def select_roi(img: np.ndarray, n_rois: int = 1) -> list:
         fig = plt.figure()
         plt.imshow(img, cmap=matplotlib.cm.gray)
         pts = fig.ginput(0, timeout=0)
-        mask = np.zeros(np.shape(img), dtype=np.int32)
+        mask = np.zeros(img.shape, dtype=np.int32)
         pts = np.asarray(pts, dtype=np.int32)
         cv2.fillConvexPoly(mask, pts, (1, 1, 1), lineType=cv2.LINE_AA)
         masks.append(mask)
@@ -317,11 +317,11 @@ def extract_components(mov_tot,
             mov_tot = mov_tot / norm_fact[:, np.newaxis, np.newaxis, np.newaxis]
         else:
             norm_fact = np.array([1., 1.])
-        c, T, d1, d2 = np.shape(mov_tot)
+        c, T, d1, d2 = mov_tot.shape
 
     else:
         norm_fact = 1
-        T, d1, d2 = np.shape(mov_tot)
+        T, d1, d2 = mov_tot.shape
         c = 1
 
     tt = time.time()

caiman/caimanmanager.py

@@ -127,7 +127,7 @@ def do_check_install(targdir: str, inplace: bool = False) -> None:
 
 
 def do_run_pytest(targdir: str) -> None:
-    out, err, ret = runcmd(["pytest", "--verbose", "--pyargs", "caiman"])
+    out, err, ret = runcmd(["pytest", "--verbose", "-W", "ignore::sklearn.exceptions.ConvergenceWarning", "--pyargs", "caiman"])
     if ret != 0:
         print(f"pytest failed with return code {ret}")
         sys.exit(ret)

caiman/cluster.py

@@ -222,12 +222,12 @@ def setup_cluster(backend:str = 'multiprocessing',
 
     if backend == 'multiprocessing' or backend == 'local':
         if backend == 'local':
-            logger.warn('The local backend is an alias for the multiprocessing backend, and the alias may be removed in some future version of Caiman')
+            logger.warning('The local backend is an alias for the multiprocessing backend, and the alias may be removed in some future version of Caiman')
         if len(multiprocessing.active_children()) > 0:
             if ignore_preexisting:
-                logger.warn('Found an existing multiprocessing pool. '
-                            'This is often indicative of an already-running CaImAn cluster. '
-                            'You have configured the cluster setup to not raise an exception.')
+                logger.warning('Found an existing multiprocessing pool. '
+                               'This is often indicative of an already-running CaImAn cluster. '
+                               'You have configured the cluster setup to not raise an exception.')
             else:
                 raise Exception(
                     'A cluster is already running. Terminate with dview.terminate() if you want to restart.')
@@ -256,7 +256,7 @@ def setup_cluster(backend:str = 'multiprocessing',
 
     elif backend == "single" or single_thread:
         if single_thread:
-            logger.warn('The single_thread flag to setup_cluster() is deprecated and may be removed in the future')
+            logger.warning('The single_thread flag to setup_cluster() is deprecated and may be removed in the future')
         dview = None
         c = None
         n_processes = 1

caiman/components_evaluation.py

@@ -67,14 +67,17 @@ def compute_event_exceptionality(traces: np.ndarray,
         erfc: ndarray
             probability at each time step of observing the N consecutive actual trace values given the distribution of noise
 
-        noise_est: ndarray
-            the components ordered according to the fitness
+        std_r:
+            Standard deviation of r
+
+        mode:
+            Mode of the traces
     """
     if N == 0:
         # Without this, numpy ranged syntax does not work correctly, and also N=0 is conceptually incoherent
         raise Exception("FATAL: N=0 is not a valid value for compute_event_exceptionality()")
 
-    T = np.shape(traces)[-1]
+    T = traces.shape[-1]
     if use_mode_fast:
         md = caiman.utils.stats.mode_robust_fast(traces, axis=1)
     else:
@@ -85,7 +88,6 @@ def compute_event_exceptionality(traces: np.ndarray,
     # only consider values under the mode to determine the noise standard deviation
     ff1 = -ff1 * (ff1 < 0)
     if robust_std:
-
         # compute 25 percentile
         ff1 = np.sort(ff1, axis=1)
         ff1[ff1 == 0] = np.nan
@@ -138,7 +140,7 @@ def find_activity_intervals(C, Npeaks: int = 5, tB=-3, tA=10, thres: float = 0.3
     # todo todocument
     logger = logging.getLogger("caiman")
 
-    K, T = np.shape(C)
+    K, T = C.shape
     L:list = []
     for i in range(K):
         if np.sum(np.abs(np.diff(C[i, :]))) == 0:
@@ -208,7 +210,7 @@ def classify_components_ep(Y, A, C, b, f, Athresh=0.1, Npeaks=5, tB=-3, tA=10, t
     """
     logger = logging.getLogger("caiman")
 
-    K, _ = np.shape(C)
+    K, _ = C.shape
     A = csc_matrix(A)
     AA = (A.T * A).toarray()
     nA = np.sqrt(np.array(A.power(2).sum(0)))
@@ -415,7 +417,7 @@ def evaluate_components(Y: np.ndarray,
     tB = np.minimum(-2, np.floor(-5. / 30 * final_frate))
     tA = np.maximum(5, np.ceil(25. / 30 * final_frate))
     logger.info(f'{tB=},{tA=}')
-    dims, T = np.shape(Y)[:-1], np.shape(Y)[-1]
+    dims, T = Y.shape[:-1], Y.shape[-1]
 
     Yr = np.reshape(Y, (np.prod(dims), T), order='F')
 
@@ -427,7 +429,7 @@ def evaluate_components(Y: np.ndarray,
 
     logger.debug('Removing Baseline')
     if remove_baseline:
-        num_samps_bl = np.minimum(np.shape(traces)[-1]// 5, 800)
+        num_samps_bl = np.minimum(traces.shape[-1]// 5, 800)
         slow_baseline = False
         if slow_baseline:
 
@@ -441,7 +443,7 @@ def evaluate_components(Y: np.ndarray,
             padbefore = int(np.floor(elm_missing / 2.))
             padafter = int(np.ceil(elm_missing / 2.))
             tr_tmp = np.pad(traces.T, ((padbefore, padafter), (0, 0)), mode='reflect')
-            numFramesNew, num_traces = np.shape(tr_tmp)
+            numFramesNew, num_traces = tr_tmp.shape
                                                                                              # compute baseline quickly
             logger.debug("binning data ...")
             tr_BL = np.reshape(tr_tmp, (downsampfact, numFramesNew // downsampfact, num_traces), order='F')