Closed DaohanZhang closed 2 months ago
import os
os.system('newgrp docker')
import numpy as np
import matplotlib.pyplot as plt
import spikeinterface.full as si # import core only
import spikeinterface.extractors as se
import spikeinterface.preprocessing as spre
import spikeinterface.sorters as ss
import spikeinterface.postprocessing as spost
import spikeinterface.qualitymetrics as sqm
import spikeinterface.comparison as sc
import spikeinterface.exporters as sexp
import spikeinterface.curation as scur
import spikeinterface.widgets as sw
from spikeinterface.sortingcomponents.peak_detection import detect_peaks
from spikeinterface.sortingcomponents.peak_localization import localize_peaks
print('#')
/home/zhangdaohan20h/.conda/envs/kilosort4/lib/python3.9/site-packages/tqdm/auto.py:21: TqdmWarning: IProgress not found. Please update jupyter and ipywidgets. See https://ipywidgets.readthedocs.io/en/stable/user_install.html
from .autonotebook import tqdm as notebook_tqdm
#
#ss.available_sorters()
print(si.__version__)
0.101.0
lfpraw = se.read_spikeglx("/home/zhangdaohan20h/public_data/NPX_examples/Pt01", load_sync_channel=False,
stream_id="imec0.lf")
lfpraw
/home/zhangdaohan20h/public_data/NPX_examples/Pt01/Pt01_g0_t0.imec0.ap.meta
/home/zhangdaohan20h/public_data/NPX_examples/Pt01/Pt01_g0_t0.imec0.lf.meta
@samuelgarcia @cwindolf. This seems like an issue with the Motion object. Any ideas? I haven't read the stack trace too closely, but I haven't looked at the Motion code at all, so better just to let you guys take a look.
Hi @zm711 , thanks for flagging this! I think that you may be the first user of LFP-based dredge in SpikeInterface @DaohanZhang , thanks for trying it out. At first glance this seems like it could relate to un-pickling of the motion object? I can try to dig into that when I get a chance but Sam probably knows better than me if that's the right idea...
Hi thanks for reporting. I can also have look. @cwindolf : did you already fixed it ?
No, haven't had a chance to look into this (or https://github.com/SpikeInterface/spikeinterface/issues/3328) in detail yet
Hello, everyone. It seems that it is a problem related to the motion object in general, I tried to spikesort after using AP based dregde and I got the same message:
SpikeSortingError Traceback (most recent call last)
Cell In[11], line 8
6 sorting = si.read_kilosort(Analysis / 'KS4' / 'sorter_output') # load sorting object
7 else:
----> 8 sorting = si.run_sorter(sorter_name = 'kilosort4', recording = recording,
9 folder = Analysis_folder / 'KS4',
10 do_correction = False, verbose = True, docker_image = False)
File [~\spikeinterface\src\spikeinterface\sorters\runsorter.py:216](http://localhost:8888/lab/tree/OneDrive%20-%20The%20University%20of%20Western%20Ontario/Documents/code/Poster%20-%20SFN%202024/Spikesorting/~/spikeinterface/src/spikeinterface/sorters/runsorter.py#line=215), in run_sorter(sorter_name, recording, folder, remove_existing_folder, delete_output_folder, verbose, raise_error, docker_image, singularity_image, delete_container_files, with_output, output_folder, **sorter_params)
205 raise RuntimeError(
206 "The python `spython` package must be installed to "
207 "run singularity. Install with `pip install spython`"
208 )
210 return run_sorter_container(
211 container_image=container_image,
212 mode=mode,
213 **common_kwargs,
214 )
--> 216 return run_sorter_local(**common_kwargs)
File [~\spikeinterface\src\spikeinterface\sorters\runsorter.py:276](http://localhost:8888/lab/tree/OneDrive%20-%20The%20University%20of%20Western%20Ontario/Documents/code/Poster%20-%20SFN%202024/Spikesorting/~/spikeinterface/src/spikeinterface/sorters/runsorter.py#line=275), in run_sorter_local(sorter_name, recording, folder, remove_existing_folder, delete_output_folder, verbose, raise_error, with_output, output_folder, **sorter_params)
274 SorterClass.set_params_to_folder(recording, folder, sorter_params, verbose)
275 SorterClass.setup_recording(recording, folder, verbose=verbose)
--> 276 SorterClass.run_from_folder(folder, raise_error, verbose)
277 if with_output:
278 sorting = SorterClass.get_result_from_folder(folder, register_recording=True, sorting_info=True)
File [~\spikeinterface\src\spikeinterface\sorters\basesorter.py:301](http://localhost:8888/lab/tree/OneDrive%20-%20The%20University%20of%20Western%20Ontario/Documents/code/Poster%20-%20SFN%202024/Spikesorting/~/spikeinterface/src/spikeinterface/sorters/basesorter.py#line=300), in BaseSorter.run_from_folder(cls, output_folder, raise_error, verbose)
298 print(f"{sorter_name} run time {run_time:0.2f}s")
300 if has_error and raise_error:
--> 301 raise SpikeSortingError(
302 f"Spike sorting error trac[e:\n](file:///E:/n){error_log_to_display}\n"
303 f"Spike sorting failed. You can inspect the runtime trace in {output_folder}/spikeinterface_log.json."
304 )
306 return run_time
SpikeSortingError: Spike sorting error trace:
Traceback (most recent call last):
File "[C:\Users\Juan\spikeinterface\src\spikeinterface\sorters\basesorter.py", line 261](file:///C:/Users/Juan/spikeinterface/src/spikeinterface/sorters/basesorter.py#line=260), in run_from_folder
SorterClass._run_from_folder(sorter_output_folder, sorter_params, verbose)
File "[C:\Users\Juan\spikeinterface\src\spikeinterface\sorters\external\kilosort4.py", line 180](file:///C:/Users/Juan/spikeinterface/src/spikeinterface/sorters/external/kilosort4.py#line=179), in _run_from_folder
recording = cls.load_recording_from_folder(sorter_output_folder.parent, with_warnings=False)
File "[C:\Users\Juan\spikeinterface\src\spikeinterface\sorters\basesorter.py", line 212](file:///C:/Users/Juan/spikeinterface/src/spikeinterface/sorters/basesorter.py#line=211), in load_recording_from_folder
recording = load_extractor(json_file, base_folder=output_folder)
File "[C:\Users\Juan\spikeinterface\src\spikeinterface\core\base.py", line 1184](file:///C:/Users/Juan/spikeinterface/src/spikeinterface/core/base.py#line=1183), in load_extractor
return BaseExtractor.load(file_or_folder_or_dict, base_folder=base_folder)
File "[C:\Users\Juan\spikeinterface\src\spikeinterface\core\base.py", line 769](file:///C:/Users/Juan/spikeinterface/src/spikeinterface/core/base.py#line=768), in load
extractor = BaseExtractor.from_dict(d, base_folder=base_folder)
File "[C:\Users\Juan\spikeinterface\src\spikeinterface\core\base.py", line 515](file:///C:/Users/Juan/spikeinterface/src/spikeinterface/core/base.py#line=514), in from_dict
extractor = _load_extractor_from_dict(dictionary)
File "[C:\Users\Juan\spikeinterface\src\spikeinterface\core\base.py", line 1123](file:///C:/Users/Juan/spikeinterface/src/spikeinterface/core/base.py#line=1122), in _load_extractor_from_dict
extractor = extractor_class(**new_kwargs)
File "[C:\Users\Juan\spikeinterface\src\spikeinterface\sortingcomponents\motion\motion_interpolation.py", line 302](file:///C:/Users/Juan/spikeinterface/src/spikeinterface/sortingcomponents/motion/motion_interpolation.py#line=301), in __init__
assert channel_locations.ndim >= motion.dim, (
AttributeError: 'dict' object has no attribute 'dim'
As @cwindolf suggested is probably something off with the pickling/unpickling. @JuanPimientoCaicedo do you have a spikeinterface_recording.json
or .pkl
in the sorting output folder?
Yes, I have it, let me share it to you, this is the code that I use to correct the motion and it seems the dim attribute is there.
# Apply drift correction
motion_info = si.load_motion_info(Analysis_folder / 'motion_correction')
recording = interpolate_motion(recording = recording,
motion = motion_info['motion'],
**motion_info['parameters']['interpolate_motion_kwargs'])
print(motion_info['motion'].dim)
1
I will keep digging into it.
Maybe during the json loading _load_extractor_from_dict
is confused by Motion objects, since they're not BaseExtractors, so they just get left as dicts? Would that mean that InterpolateMotionRecordings are only pickleable and not json-serializable? This is just my guess -- maybe there is also a way to make them work with json
I think I misunderstood your question, in the sorter_output folder I only have the probe.prb file
Ok, I guess that is the issue. When we dump to JSON, we call the Motion.to_dict()
here. The Motion dict is then reloaded as a dict, and not instantiated as a Motion object anywhere..
I think I misunderstood your question, in the sorter_output folder I only have the probe.prb file
No no you perfectly understood! That's what I wanted to know :)
A quick fix would be to make the InterpolateMotionRecording
not JSON serializable. This will generate a pickle file instead of a JSON file, and probably will fix this issue.
Let me make a test PR for you guys to try!
Can you test this out? https://github.com/SpikeInterface/spikeinterface/pull/3341
@alejoe91 yes, it works. I was able to pass that point.
Thank you!
Thanks a lot! But actually I do not think this issue got entirely fixed, despite testing results using interpolation method dredge_ap seem okay. I still encountered the similar issue reported in https://github.com/SpikeInterface/spikeinterface/issues/3328
import numpy as np
import matplotlib.pyplot as plt
import spikeinterface.full as si # import core only
import spikeinterface.extractors as se
import spikeinterface.preprocessing as spre
import spikeinterface.sorters as ss
import spikeinterface.postprocessing as spost
import spikeinterface.qualitymetrics as sqm
import spikeinterface.comparison as sc
import spikeinterface.exporters as sexp
import spikeinterface.curation as scur
import spikeinterface.widgets as sw
from spikeinterface.sortingcomponents.peak_detection import detect_peaks
from spikeinterface.sortingcomponents.peak_localization import localize_peaks
from spikeinterface.sortingcomponents.motion import estimate_motion
from spikeinterface import set_global_job_kwargs
global_job_kwargs = dict(n_jobs=20, chunk_duration="1s", progress_bar=True)
set_global_job_kwargs(**global_job_kwargs)
lfpraw = se.read_spikeglx("/home/zhangdaohan20h/public_data/NPX_examples/Pt01/", load_sync_channel=False,
stream_id="imec0.lf")
lfprec = si.astype(lfpraw, np.float32)
lfprec = si.depth_order(lfprec)
cutoff_um = 8000
if cutoff_um is not None:
geom = lfprec.get_channel_locations()
lfprec = lfprec.remove_channels(lfprec.channel_ids[geom[:, 1] > cutoff_um])
lfprec = si.bandpass_filter(
lfprec,
freq_min=0.5,
freq_max=250,
margin_ms=1000,
filter_order=3,
dtype="float32",
add_reflect_padding=True,
)
bad_chans, labels = si.detect_bad_channels(lfprec, psd_hf_threshold=1.4, num_random_chunks=100, seed=0)
print("Found bad channels", bad_chans)
lfprec = lfprec.remove_channels(bad_chans)
lfprec = si.phase_shift(lfprec)
lfprec = si.common_reference(lfprec)
lfprec = si.resample(lfprec, 250, margin_ms=1000)
lfprec = si.directional_derivative(lfprec, order=2, edge_order=1)
lfprec = si.average_across_direction(lfprec)
from spikeinterface.sortingcomponents.motion import estimate_motion
motion = estimate_motion(lfprec, method='dredge_lfp', rigid=False, progress_bar=True, max_disp_um=1000)
rec = se.read_spikeglx("/home/zhangdaohan20h/public_data/NPX_examples/Pt01",
load_sync_channel=False, stream_id="imec0.ap")
rec = si.astype(rec, np.float32)
rec = si.bandpass_filter(rec)
rec = si.common_reference(rec)
#interpolate 30kHz AP with motion of 250Hz downsampled LFP
rec = interpolate_motion(rec, motion, border_mode='remove_channels',
spatial_interpolation_method='kriging', sigma_um=20.0, p=1,
num_closest=3, interpolation_time_bin_centers_s=None,
interpolation_time_bin_size_s=None, dtype=None)
sorting = ss.run_sorter(sorter_name="kilosort4", recording=rec, docker_image=False, remove_existing_folder=True, delete_output_folder = False,
torch_device=("cuda:0"), verbose=True)
sorting
========================================
Loading recording with SpikeInterface...
number of samples: 25014692
number of channels: 70
numbef of segments: 1
sampling rate: 30000.0
dtype: float32
========================================
Preprocessing filters computed in 355.39s; total 355.39s
computing drift
Re-computing universal templates from data.
100%|█████████▉| 416/417 [2:20:33<00:20, 20.27s/it]
Error running kilosort4
---------------------------------------------------------------------------
SpikeSortingError Traceback (most recent call last)
Cell In[9], line 1
----> 1 sorting = ss.run_sorter(sorter_name="kilosort4", recording=rec, docker_image=False, remove_existing_folder=True, delete_output_folder = False,
2 torch_device=("cuda:0"), verbose=True)
4 sorting
File ~/.conda/envs/kilosort4/lib/python3.9/site-packages/spikeinterface/sorters/runsorter.py:216, in run_sorter(sorter_name, recording, folder, remove_existing_folder, delete_output_folder, verbose, raise_error, docker_image, singularity_image, delete_container_files, with_output, output_folder, **sorter_params)
205 raise RuntimeError(
206 "The python `spython` package must be installed to "
207 "run singularity. Install with `pip install spython`"
208 )
210 return run_sorter_container(
211 container_image=container_image,
212 mode=mode,
213 **common_kwargs,
214 )
--> 216 return run_sorter_local(**common_kwargs)
File ~/.conda/envs/kilosort4/lib/python3.9/site-packages/spikeinterface/sorters/runsorter.py:276, in run_sorter_local(sorter_name, recording, folder, remove_existing_folder, delete_output_folder, verbose, raise_error, with_output, output_folder, **sorter_params)
274 SorterClass.set_params_to_folder(recording, folder, sorter_params, verbose)
275 SorterClass.setup_recording(recording, folder, verbose=verbose)
--> 276 SorterClass.run_from_folder(folder, raise_error, verbose)
277 if with_output:
278 sorting = SorterClass.get_result_from_folder(folder, register_recording=True, sorting_info=True)
File ~/.conda/envs/kilosort4/lib/python3.9/site-packages/spikeinterface/sorters/basesorter.py:301, in BaseSorter.run_from_folder(cls, output_folder, raise_error, verbose)
298 print(f"{sorter_name} run time {run_time:0.2f}s")
300 if has_error and raise_error:
--> 301 raise SpikeSortingError(
302 f"Spike sorting error trace:\n{error_log_to_display}\n"
303 f"Spike sorting failed. You can inspect the runtime trace in {output_folder}/spikeinterface_log.json."
304 )
306 return run_time
SpikeSortingError: Spike sorting error trace:
Traceback (most recent call last):
File "/home/zhangdaohan20h/.conda/envs/kilosort4/lib/python3.9/site-packages/spikeinterface/sorters/basesorter.py", line 261, in run_from_folder
SorterClass._run_from_folder(sorter_output_folder, sorter_params, verbose)
File "/home/zhangdaohan20h/.conda/envs/kilosort4/lib/python3.9/site-packages/spikeinterface/sorters/external/kilosort4.py", line 256, in _run_from_folder
ops, bfile, st0 = compute_drift_correction(
File "/home/zhangdaohan20h/.conda/envs/kilosort4/lib/python3.9/site-packages/kilosort/run_kilosort.py", line 344, in compute_drift_correction
ops, st = datashift.run(ops, bfile, device=device, progress_bar=progress_bar)
File "/home/zhangdaohan20h/.conda/envs/kilosort4/lib/python3.9/site-packages/kilosort/datashift.py", line 192, in run
st, _, ops = spikedetect.run(ops, bfile, device=device, progress_bar=progress_bar)
File "/home/zhangdaohan20h/.conda/envs/kilosort4/lib/python3.9/site-packages/kilosort/spikedetect.py", line 231, in run
X = bfile.padded_batch_to_torch(ibatch, ops)
File "/home/zhangdaohan20h/.conda/envs/kilosort4/lib/python3.9/site-packages/kilosort/io.py", line 666, in padded_batch_to_torch
X = super().padded_batch_to_torch(ibatch)
File "/home/zhangdaohan20h/.conda/envs/kilosort4/lib/python3.9/site-packages/kilosort/io.py", line 486, in padded_batch_to_torch
data = self.file[bstart : bend]
File "/home/zhangdaohan20h/.conda/envs/kilosort4/lib/python3.9/site-packages/kilosort/io.py", line 899, in __getitem__
samples = self.recording.get_traces(start_frame=i, end_frame=j,
File "/home/zhangdaohan20h/.conda/envs/kilosort4/lib/python3.9/site-packages/spikeinterface/core/baserecording.py", line 342, in get_traces
traces = rs.get_traces(start_frame=start_frame, end_frame=end_frame, channel_indices=channel_indices)
File "/home/zhangdaohan20h/.conda/envs/kilosort4/lib/python3.9/site-packages/spikeinterface/sortingcomponents/motion/motion_interpolation.py", line 467, in get_traces
traces = interpolate_motion_on_traces(
File "/home/zhangdaohan20h/.conda/envs/kilosort4/lib/python3.9/site-packages/spikeinterface/sortingcomponents/motion/motion_interpolation.py", line 139, in interpolate_motion_on_traces
bin_time = time_bins[bin_ind]
IndexError: index 208455 is out of bounds for axis 0 with size 208455
Spike sorting failed. You can inspect the runtime trace in /share/home/zhangdaohan20h/CODES/NPX/preanalysis/MGH01/kilosort4_output/spikeinterface_log.json.
The index 208455 is the number of timestep of motion.displacement calculated by dredge_lfp. I' m a little bit doubtful about @cwindolf "interpolate_motion should work regardless of what time bins the Motion object has." https://github.com/SpikeInterface/spikeinterface/issues/3328, which should not lead to such error in the following sorting procedure. Besides, the same issue did not raise directly in ks3 and ks2.5, maybe it is a hidden error which has not been recognized.
Hi! Using ks4 without docker is terribly convenient dealing with Neuropixels recording from sipkeGLX! But after interpolating my action potential (AP) with motion calculated from local field potentials (LFP) using both
dredge_lfp
in spikeinterface or rigid_registration inDREDge, I failed to run Kilosort4. There raised an error that variance "motion", which should not be a dict, happened to be a dict and lost its attribute "motion.dim". I was using the lattest version (V0.101.0) of spikeinterface, and my environment meet the requirements, since kilosort4 can be applied to raw ap data smoothly. My notebook is attached below Thanks!