From 3131926a5d26d3965eb6a8c689ab9b4decef1170 Mon Sep 17 00:00:00 2001 From: Matthias H Hennig Date: Fri, 7 Jun 2024 23:53:19 +0100 Subject: [PATCH 1/9] fixed brackets in warning --- src/spikeinterface/comparison/multicomparisons.py | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/src/spikeinterface/comparison/multicomparisons.py b/src/spikeinterface/comparison/multicomparisons.py index 77adcaa8ca..247cc82f63 100644 --- a/src/spikeinterface/comparison/multicomparisons.py +++ b/src/spikeinterface/comparison/multicomparisons.py @@ -186,7 +186,7 @@ def save_to_folder(self, save_folder): warnings.warn( "save_to_folder() is deprecated. " "You should save and load the multi sorting comparison object using pickle." - "\n>>> pickle.dump(mcmp, open('mcmp.pkl', 'wb')))))\n>>> mcmp_loaded = pickle.load(open('mcmp.pkl', 'rb'))", + "\n>>> pickle.dump(mcmp, open('mcmp.pkl', 'wb'))\n>>> mcmp_loaded = pickle.load(open('mcmp.pkl', 'rb'))", DeprecationWarning, stacklevel=2, ) @@ -218,7 +218,7 @@ def load_from_folder(folder_path): warnings.warn( "load_from_folder() is deprecated. " "You should save and load the multi sorting comparison object using pickle." - "\n>>> pickle.dump(mcmp, open('mcmp.pkl', 'wb')))))\n>>> mcmp_loaded = pickle.load(open('mcmp.pkl', 'rb'))", + "\n>>> pickle.dump(mcmp, open('mcmp.pkl', 'wb'))\n>>> mcmp_loaded = pickle.load(open('mcmp.pkl', 'rb'))", DeprecationWarning, stacklevel=2, ) From 81904bafe75c436487576c59d033a32ff3abc4ca Mon Sep 17 00:00:00 2001 From: mhhennig Date: Mon, 15 Jul 2024 21:58:07 +0100 Subject: [PATCH 2/9] Now exclusive support for HS v0.4 (Lightning) --- .../sorters/external/herdingspikes.py | 225 ++++++------------ 1 file changed, 77 insertions(+), 148 deletions(-) diff --git a/src/spikeinterface/sorters/external/herdingspikes.py b/src/spikeinterface/sorters/external/herdingspikes.py index 9b9b751534..e4b9c18d9d 100644 --- a/src/spikeinterface/sorters/external/herdingspikes.py +++ b/src/spikeinterface/sorters/external/herdingspikes.py @@ -1,13 +1,10 @@ from __future__ import annotations from pathlib import Path -import copy from packaging import version from ..basesorter import BaseSorter -from spikeinterface.core.old_api_utils import NewToOldRecording -from spikeinterface.core import load_extractor from spikeinterface.extractors import HerdingspikesSortingExtractor @@ -19,90 +16,72 @@ class HerdingspikesSorter(BaseSorter): requires_locations = True compatible_with_parallel = {"loky": True, "multiprocessing": True, "threading": False} _default_params = { - # core params - "clustering_bandwidth": 5.5, # 5.0, - "clustering_alpha": 5.5, # 5.0, + "chunk_size": None, + "rescale": True, + "rescale_value": -1280.0, + "common_reference": "median", + "spike_duration": 1.0, + "amp_avg_duration": 0.4, + "threshold": 8.0, + "min_avg_amp": 1.0, + "AHP_thr": 0.0, + "neighbor_radius": 90.0, + "inner_radius": 70.0, + "peak_jitter": 0.25, + "rise_duration": 0.26, + "decay_filtering": False, + "decay_ratio": 1.0, + "localize": True, + "save_shape": True, + "out_file": "HS2_detected", + "left_cutout_time": 0.3, + "right_cutout_time": 1.8, + "verbose": True, + "clustering_bandwidth": 4.0, + "clustering_alpha": 4.5, "clustering_n_jobs": -1, "clustering_bin_seeding": True, - "clustering_min_bin_freq": 16, # 10, + "clustering_min_bin_freq": 4, "clustering_subset": None, - "left_cutout_time": 0.3, # 0.2, - "right_cutout_time": 1.8, # 0.8, - "detect_threshold": 20, # 24, #15, - # extra probe params - "probe_masked_channels": [], - "probe_inner_radius": 70, - "probe_neighbor_radius": 90, - "probe_event_length": 0.26, - "probe_peak_jitter": 0.2, - # extra detection params - "t_inc": 100000, - "num_com_centers": 1, - "maa": 12, - "ahpthr": 11, - "out_file_name": "HS2_detected", - "decay_filtering": False, - "save_all": False, - "amp_evaluation_time": 0.4, # 0.14, - "spk_evaluation_time": 1.0, - # extra pca params "pca_ncomponents": 2, "pca_whiten": True, - # bandpass filter - "freq_min": 300.0, - "freq_max": 6000.0, - "filter": True, - # rescale traces - "pre_scale": True, - "pre_scale_value": 20.0, - # remove duplicates (based on spk_evaluation_time) - "filter_duplicates": True, } _params_description = { - # core params - "clustering_bandwidth": "Meanshift bandwidth, average spatial extent of spike clusters (um)", - "clustering_alpha": "Scalar for the waveform PC features when clustering.", - "clustering_n_jobs": "Number of cores to use for clustering.", - "clustering_bin_seeding": "Enable clustering bin seeding.", - "clustering_min_bin_freq": "Minimum spikes per bin for bin seeding.", - "clustering_subset": "Number of spikes used to build clusters. All by default.", - "left_cutout_time": "Cutout size before peak (ms).", - "right_cutout_time": "Cutout size after peak (ms).", - "detect_threshold": "Detection threshold", - # extra probe params - "probe_masked_channels": "Masked channels", - "probe_inner_radius": "Radius of area around probe channel for localization", - "probe_neighbor_radius": "Radius of area around probe channel for neighbor classification.", - "probe_event_length": "Duration of a spike event (ms)", - "probe_peak_jitter": "Maximum peak misalignment for synchronous spike (ms)", - # extra detection params - "t_inc": "Number of samples per chunk during detection.", - "num_com_centers": "Number of centroids to average when localizing.", - "maa": "Minimum summed spike amplitude for spike acceptance.", - "ahpthr": "Requires magnitude of spike rebound for acceptance", - "out_file_name": "File name for storage of unclustered detected spikes", - "decay_filtering": "Experimental: Set to True at your risk", - "save_all": "Save all working files after sorting (slow)", - "amp_evaluation_time": "Amplitude evaluation time (ms)", - "spk_evaluation_time": "Spike evaluation time (ms)", - # extra pca params - "pca_ncomponents": "Number of principal components to use when clustering", - "pca_whiten": "If true, whiten data for pca", - # bandpass filter - "freq_min": "High-pass filter cutoff frequency", - "freq_max": "Low-pass filter cutoff frequency", - "filter": "Enable or disable filter", - # rescale traces - "pre_scale": "Scales recording traces to optimize HerdingSpikes performance", - "pre_scale_value": "Scale to apply in case of pre-scaling of traces", - # remove duplicates (based on spk_evaluation_time) - "filter_duplicates": "Remove spike duplicates (based on spk_evaluation_time)", + "localize": "Perform spike localization. (`bool`, `True`)", + "save_shape": "Save spike shape. (`bool`, `True`)", + "out_file": "Path and filename to store detection and clustering results. (`str`, `HS2_detected`)", + "verbose": "Print progress information. (`bool`, `True`)", + "chunk_size": " Number of samples per chunk during detection. If `None`, a suitable value will be estimated. (`int`, `None`)", + "common_reference": "Method for common reference filtering, can be `average` or `median` (`str`, `median`)", + "rescale": "Automatically re-scale the data. (`bool`, `True`)", + "rescale_value": "Factor by which data is re-scaled. (`float`, `-1280.0`)", + "threshold": "Spike detection threshold. (`float`, `8.0`)", + "spike_duration": "Maximum duration over which a spike is evaluated (ms). (`float`, `1.0`)", + "amp_avg_duration": "Maximum duration over which the spike amplitude is evaluated (ms). (`float`, `0.4`)", + "min_avg_amp": "Minimum integrated spike amplitude for a true spike. (`float`, `1.0`)", + "AHP_thr": "Minimum value of the spike repolarisation for a true spike. (`float`, `1.0`)", + "neighbor_radius": "Radius of area around probe channel for neighbor classification (microns). (`float`, `90.0`)", + "inner_radius": "Radius of area around probe channel for spike localisation (microns). (`float`, `70.0`)", + "peak_jitter": "Maximum peak misalignment for synchronous spike (ms). (`float`, `0.25`)", + "rise_duration": "Maximum spike rise time, in milliseconds. (`float`, `0.26`)", + "decay_filtering": "Exclude duplicate spikes based on spatial decay pattern, experimental. (`bool`,`False`)", + "decay_ratio": "Spatial decay rate for `decay_filtering`. (`float`,`1.0`)", + "left_cutout_time": "Length of cutout before peak (ms). (`float`, `0.3`)", + "right_cutout_time": "Length of cutout after peak (ms). (`float`, `1.8`)", + "pca_ncomponents": "Number of principal components to use when clustering. (`int`, `2`)", + "pca_whiten": "If `True`, whiten data for PCA. (`bool`, `True`)", + "clustering_bandwidth": "Meanshift bandwidth, average spatial extent of spike clusters (microns). (`float`, `4.0`)", + "clustering_alpha": "Scalar for the waveform PC features when clustering. (`float`, `4.5`)", + "clustering_n_jobs": "Number of cores to use for clustering, use `-1` for all available cores. (`int`, `-1`)", + "clustering_bin_seeding": "Enable clustering bin seeding. (`bool`, `True`)", + "clustering_min_bin_freq": "Minimum spikes per bin for bin seeding. (`int`, `4`)", + "clustering_subset": "Number of spikes used to build clusters. All by default. (`int`, `None`)", } - sorter_description = """Herding Spikes is a density-based spike sorter designed for high-density retinal recordings. + sorter_description = """Herding Spikes is a density-based spike sorter designed for large-scale high-density recordings. It uses both PCA features and an estimate of the spike location to cluster different units. - For more information see https://doi.org/10.1016/j.jneumeth.2016.06.006""" + For more information see https://www.sciencedirect.com/science/article/pii/S221112471730236X""" installation_mesg = """\nTo use HerdingSpikes run:\n >>> pip install herdingspikes @@ -134,96 +113,46 @@ def _check_apply_filter_in_params(cls, params): @classmethod def _setup_recording(cls, recording, sorter_output_folder, params, verbose): - # nothing to copy inside the folder : Herdingspikes used natively spikeinterface + # nothing to copy inside the folder : Herdingspikes uses spikeinterface natively pass @classmethod def _run_from_folder(cls, sorter_output_folder, params, verbose): import herdingspikes as hs - from spikeinterface.preprocessing import bandpass_filter, normalize_by_quantile hs_version = version.parse(hs.__version__) - if hs_version >= version.parse("0.3.99"): - new_api = True + if hs_version >= version.parse("0.4"): + lightning_api = True else: - new_api = False + lightning_api = False + + assert ( + lightning_api + ), "HerdingSpikes version <0.4 is no longer supported. run:\n>>> pip install --upgrade herdingspikes" recording = cls.load_recording_from_folder(sorter_output_folder.parent, with_warnings=False) + sorted_file = str(sorter_output_folder / "HS2_sorted") + params["out_file"] = str(sorter_output_folder / "HS2_detected") p = params - # Bandpass filter - if p["filter"] and p["freq_min"] is not None and p["freq_max"] is not None: - recording = bandpass_filter(recording=recording, freq_min=p["freq_min"], freq_max=p["freq_max"]) - - if p["pre_scale"]: - recording = normalize_by_quantile( - recording=recording, scale=p["pre_scale_value"], median=0.0, q1=0.05, q2=0.95 - ) - - if new_api: - recording_to_hs = recording - else: - print( - "herdingspikes version<0.3.99 uses the OLD spikeextractors with NewToOldRecording.\n" - "Consider updating herdingspikes (pip install herdingspikes>=0.3.99)" - ) - recording_to_hs = NewToOldRecording(recording) - - # this should have its name changed - Probe = hs.probe.RecordingExtractor( - recording_to_hs, - masked_channels=p["probe_masked_channels"], - inner_radius=p["probe_inner_radius"], - neighbor_radius=p["probe_neighbor_radius"], - event_length=p["probe_event_length"], - peak_jitter=p["probe_peak_jitter"], + det = hs.HSDetectionLightning(recording, p) + det.DetectFromRaw() + C = hs.HSClustering(det) + C.ShapePCA() + C.CombinedClustering( + alpha=p["clustering_alpha"], + cluster_subset=p["clustering_subset"], + bandwidth=p["clustering_bandwidth"], + bin_seeding=p["clustering_bin_seeding"], + min_bin_freq=p["clustering_min_bin_freq"], + n_jobs=p["clustering_n_jobs"], ) - H = hs.HSDetection( - Probe, - file_directory_name=str(sorter_output_folder), - left_cutout_time=p["left_cutout_time"], - right_cutout_time=p["right_cutout_time"], - threshold=p["detect_threshold"], - to_localize=True, - num_com_centers=p["num_com_centers"], - maa=p["maa"], - ahpthr=p["ahpthr"], - out_file_name=p["out_file_name"], - decay_filtering=p["decay_filtering"], - save_all=p["save_all"], - amp_evaluation_time=p["amp_evaluation_time"], - spk_evaluation_time=p["spk_evaluation_time"], - ) - - H.DetectFromRaw(load=True, tInc=int(p["t_inc"])) - - sorted_file = str(sorter_output_folder / "HS2_sorted.hdf5") - if not H.spikes.empty: - C = hs.HSClustering(H) - C.ShapePCA(pca_ncomponents=p["pca_ncomponents"], pca_whiten=p["pca_whiten"]) - C.CombinedClustering( - alpha=p["clustering_alpha"], - cluster_subset=p["clustering_subset"], - bandwidth=p["clustering_bandwidth"], - bin_seeding=p["clustering_bin_seeding"], - n_jobs=p["clustering_n_jobs"], - min_bin_freq=p["clustering_min_bin_freq"], - ) - else: - C = hs.HSClustering(H) - - if p["filter_duplicates"]: - uids = C.spikes.cl.unique() - for u in uids: - s = C.spikes[C.spikes.cl == u].t.diff() < p["spk_evaluation_time"] / 1000 * Probe.fps - C.spikes = C.spikes.drop(s.index[s]) - if verbose: print("Saving to", sorted_file) - C.SaveHDF5(sorted_file, sampling=Probe.fps) + C.SaveHDF5(sorted_file, sampling=recording.get_sampling_frequency()) @classmethod def _get_result_from_folder(cls, sorter_output_folder): From 023614d6891478779cf5849a0ad1ddf336ae0bf5 Mon Sep 17 00:00:00 2001 From: mhhennig Date: Tue, 16 Jul 2024 09:05:54 +0100 Subject: [PATCH 3/9] updated documentation --- doc/get_started/install_sorters.rst | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) diff --git a/doc/get_started/install_sorters.rst b/doc/get_started/install_sorters.rst index 12233784fd..61abc0f129 100644 --- a/doc/get_started/install_sorters.rst +++ b/doc/get_started/install_sorters.rst @@ -43,9 +43,10 @@ Herdingspikes2 * Python + C++ * Url: https://github.com/mhhennig/hs2 -* Authors: Matthias Hennig, Jano Horvath,Cole Hurwitz, Oliver Muthmann, Albert Puente Encinas, Martino Sorbaro, Cesar Juarez Ramirez, Raimon Wintzer: GUI and visualisation +* Authors: Matthias Hennig, Jano Horvath, Cole Hurwitz, Rickey K. Liang, Oliver Muthmann, Albert Puente Encinas, Martino Sorbaro, Cesar Juarez Ramirez, Raimon Wintzer * Installation:: + pip install cython numpy pip install herdingspikes From bfcfabe0120a29690faa2280f7271453385123cf Mon Sep 17 00:00:00 2001 From: mhhennig Date: Tue, 16 Jul 2024 09:10:42 +0100 Subject: [PATCH 4/9] updated sorter doc and fix in confusion matrix --- src/spikeinterface/comparison/comparisontools.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/spikeinterface/comparison/comparisontools.py b/src/spikeinterface/comparison/comparisontools.py index 87d0bf512b..809f512591 100644 --- a/src/spikeinterface/comparison/comparisontools.py +++ b/src/spikeinterface/comparison/comparisontools.py @@ -707,7 +707,7 @@ def do_confusion_matrix(event_counts1, event_counts2, match_12, match_event_coun N2 = len(unit2_ids) matched_units1 = match_12[match_12 != -1].index - matched_units2 = match_12[match_12 != -1].values + matched_units2 = match_12[match_12 != -1].values.astype(int) unmatched_units1 = match_12[match_12 == -1].index unmatched_units2 = unit2_ids[~np.isin(unit2_ids, matched_units2)] From f92c790df18caf53266c38ef1b76f5fb334ee406 Mon Sep 17 00:00:00 2001 From: "Matthias H. Hennig" Date: Tue, 16 Jul 2024 14:23:44 +0100 Subject: [PATCH 5/9] revert change to comparisontools.py Correct, the table has `sorting2.get_unit_ids()` - however `pd.DataFrame` casts these to `float` even when all channel inds are `int`. Needs type checking in `make_match_count_matrix`... --- src/spikeinterface/comparison/comparisontools.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/spikeinterface/comparison/comparisontools.py b/src/spikeinterface/comparison/comparisontools.py index 809f512591..87d0bf512b 100644 --- a/src/spikeinterface/comparison/comparisontools.py +++ b/src/spikeinterface/comparison/comparisontools.py @@ -707,7 +707,7 @@ def do_confusion_matrix(event_counts1, event_counts2, match_12, match_event_coun N2 = len(unit2_ids) matched_units1 = match_12[match_12 != -1].index - matched_units2 = match_12[match_12 != -1].values.astype(int) + matched_units2 = match_12[match_12 != -1].values unmatched_units1 = match_12[match_12 == -1].index unmatched_units2 = unit2_ids[~np.isin(unit2_ids, matched_units2)] From 05c62613cd091525a220bdd4a09e4041e9dc33a1 Mon Sep 17 00:00:00 2001 From: "Matthias H. Hennig" Date: Tue, 16 Jul 2024 17:36:36 +0100 Subject: [PATCH 6/9] data is no longer filtered --- src/spikeinterface/sorters/external/herdingspikes.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/spikeinterface/sorters/external/herdingspikes.py b/src/spikeinterface/sorters/external/herdingspikes.py index e4b9c18d9d..3c22b8ba02 100644 --- a/src/spikeinterface/sorters/external/herdingspikes.py +++ b/src/spikeinterface/sorters/external/herdingspikes.py @@ -109,7 +109,7 @@ def get_sorter_version(cls): @classmethod def _check_apply_filter_in_params(cls, params): - return params["filter"] + return False @classmethod def _setup_recording(cls, recording, sorter_output_folder, params, verbose): From aa72e8d30d6b64cda9e9cef7acd9e4d8a9e397ec Mon Sep 17 00:00:00 2001 From: mhhennig Date: Wed, 17 Jul 2024 10:41:27 +0100 Subject: [PATCH 7/9] fix file name --- src/spikeinterface/sorters/external/herdingspikes.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/spikeinterface/sorters/external/herdingspikes.py b/src/spikeinterface/sorters/external/herdingspikes.py index e4b9c18d9d..81b5264be2 100644 --- a/src/spikeinterface/sorters/external/herdingspikes.py +++ b/src/spikeinterface/sorters/external/herdingspikes.py @@ -133,7 +133,7 @@ def _run_from_folder(cls, sorter_output_folder, params, verbose): recording = cls.load_recording_from_folder(sorter_output_folder.parent, with_warnings=False) - sorted_file = str(sorter_output_folder / "HS2_sorted") + sorted_file = str(sorter_output_folder / "HS2_sorted.hdf5") params["out_file"] = str(sorter_output_folder / "HS2_detected") p = params From 2e141a925e52c4a3592c0ba702cbc1e2a24cfd25 Mon Sep 17 00:00:00 2001 From: Matthias H Hennig Date: Fri, 19 Jul 2024 23:08:54 +0100 Subject: [PATCH 8/9] bump version with msvc fix --- src/spikeinterface/sorters/external/herdingspikes.py | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/src/spikeinterface/sorters/external/herdingspikes.py b/src/spikeinterface/sorters/external/herdingspikes.py index e09f579191..16ad51fb5b 100644 --- a/src/spikeinterface/sorters/external/herdingspikes.py +++ b/src/spikeinterface/sorters/external/herdingspikes.py @@ -122,14 +122,14 @@ def _run_from_folder(cls, sorter_output_folder, params, verbose): hs_version = version.parse(hs.__version__) - if hs_version >= version.parse("0.4"): + if hs_version >= version.parse("0.4.001"): lightning_api = True else: lightning_api = False assert ( lightning_api - ), "HerdingSpikes version <0.4 is no longer supported. run:\n>>> pip install --upgrade herdingspikes" + ), "HerdingSpikes version <0.4.001 is no longer supported. run:\n>>> pip install --upgrade herdingspikes" recording = cls.load_recording_from_folder(sorter_output_folder.parent, with_warnings=False) From a20c5c309205c625e7a156a503ff2f1f3c7f4e1f Mon Sep 17 00:00:00 2001 From: "Matthias H. Hennig" Date: Fri, 19 Jul 2024 23:23:12 +0100 Subject: [PATCH 9/9] Update src/spikeinterface/sorters/external/herdingspikes.py Co-authored-by: Chris Halcrow <57948917+chrishalcrow@users.noreply.github.com> --- src/spikeinterface/sorters/external/herdingspikes.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/spikeinterface/sorters/external/herdingspikes.py b/src/spikeinterface/sorters/external/herdingspikes.py index 16ad51fb5b..a84c05c240 100644 --- a/src/spikeinterface/sorters/external/herdingspikes.py +++ b/src/spikeinterface/sorters/external/herdingspikes.py @@ -60,7 +60,7 @@ class HerdingspikesSorter(BaseSorter): "spike_duration": "Maximum duration over which a spike is evaluated (ms). (`float`, `1.0`)", "amp_avg_duration": "Maximum duration over which the spike amplitude is evaluated (ms). (`float`, `0.4`)", "min_avg_amp": "Minimum integrated spike amplitude for a true spike. (`float`, `1.0`)", - "AHP_thr": "Minimum value of the spike repolarisation for a true spike. (`float`, `1.0`)", + "AHP_thr": "Minimum value of the spike repolarisation for a true spike. (`float`, `0.0`)", "neighbor_radius": "Radius of area around probe channel for neighbor classification (microns). (`float`, `90.0`)", "inner_radius": "Radius of area around probe channel for spike localisation (microns). (`float`, `70.0`)", "peak_jitter": "Maximum peak misalignment for synchronous spike (ms). (`float`, `0.25`)",