mda extractors failed when sorting concatenated files

[yuhang21]

Hi I want to sort my concatenate .mda files which loaded using concatenate_recordings. I used mountainsort 5 and load the sorter using run_sorter My concatenated recording is like the following

ConcatenateSegmentRecording: 128 channels - 30.0kHz - 1 segments - 78,947,056 samples
2,631.57s (43.86 minutes) - int16 dtype - -1264390144.00 B
C:\Users\xiela\anaconda3\Lib\site-packages\spikeinterface\core\baserecording.py:222: RuntimeWarning: overflow encountered in scalar multiply
memory_bytes = num_samples * num_channels * dtype_size_bytes

After running the sorter, there is an error from the mdaextractor.py

File ~\anaconda3\Lib\site-packages\spikeinterface\sorters\basesorter.py:289 in run_from_folder
raise SpikeSortingError(

SpikeSortingError: Spike sorting error trace:
Traceback (most recent call last):
File "C:\Users\xiela\anaconda3\Lib\site-packages\spikeinterface\sorters\basesorter.py", line 254, in run_from_folder
SorterClass._run_from_folder(sorter_output_folder, sorter_params, verbose)
File "C:\Users\xiela\anaconda3\Lib\site-packages\spikeinterface\sorters\external\mountainsort5.py", line 176, in _run_from_folder
sorting = ms5.sorting_scheme2(recording=recording, sorting_parameters=scheme2_sorting_parameters)
^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
File "C:\Users\xiela\anaconda3\Lib\site-packages\mountainsort5\schemes\sorting_scheme2.py", line 74, in sorting_scheme2
training_recording = get_sampled_recording_for_training(
^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
File "C:\Users\xiela\anaconda3\Lib\site-packages\mountainsort5\core\get_sampled_recording_for_training.py", line 44, in get_sampled_recording_for_training
traces_list.append(recording.get_traces(start_frame=start_frame, end_frame=end_frame))
^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
File "C:\Users\xiela\anaconda3\Lib\site-packages\spikeinterface\core\baserecording.py", line 282, in get_traces
traces = rs.get_traces(start_frame=start_frame, end_frame=end_frame, channel_indices=channel_indices)
^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
File "C:\Users\xiela\anaconda3\Lib\site-packages\spikeinterface\preprocessing\whiten.py", line 107, in get_traces
traces = self.parent_recording_segment.get_traces(start_frame, end_frame, slice(None))
^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
File "C:\Users\xiela\anaconda3\Lib\site-packages\spikeinterface\preprocessing\common_reference.py", line 146, in get_traces
all_traces = self.parent_recording_segment.get_traces(start_frame, end_frame, slice(None))
^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
File "C:\Users\xiela\anaconda3\Lib\site-packages\spikeinterface\preprocessing\filter.py", line 130, in get_traces
traces_chunk, left_margin, right_margin = get_chunk_with_margin(
^^^^^^^^^^^^^^^^^^^^^^
File "C:\Users\xiela\anaconda3\Lib\site-packages\spikeinterface\core\recording_tools.py", line 232, in get_chunk_with_margin
traces_chunk = rec_segment.get_traces(
^^^^^^^^^^^^^^^^^^^^^^^
File "C:\Users\xiela\anaconda3\Lib\site-packages\spikeinterface\core\segmentutils.py", line 187, in get_traces
traces = rec_seg.get_traces(start_frame - seg_start, end_frame - seg_start, channel_indices)
^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
File "C:\Users\xiela\anaconda3\Lib\site-packages\spikeinterface\extractors\mdaextractors.py", line 174, in get_traces
recordings = recordings[channel_indices, :].T
~~~~~~~~~~^^^^^^^^^^^^^^^^^^^^
TypeError: 'NoneType' object is not subscriptable

Spike sorting failed. You can inspect the runtime trace in E:\EBL data\tracking_me\Processed_01-Mar-2024_221810\sorting/spikeinterface_log.json.

I used the same workflow with other binary format of data and works fine. Probably the issue relates to the mdaextractor which need some modifications? Thank you so much for your help.

[yuhang21]

print some readchunk. The single . mda works fine if sorted itself, while the concatenated mda doesn't proceed, no matter how large/ small the files are.

Using training recording of duration 300 sec with the sampling mode uniform
Start frame: 0, End frame: 300150, Channel indices: slice(None, None, None)
Recordings after readChunk: [[-314 -242 -269 ... 187 197 138]
[1371 1357 1308 ... 242 199 86]
[1310 1427 1584 ... 159 165 63]
...
[ -70 -85 -66 ... -46 -67 -132]
[ 634 607 629 ... 70 58 30]
[ 193 208 143 ... -62 55 2]]
Start frame: 2711817, End frame: 3012117, Channel indices: slice(None, None, None)
Recordings after readChunk: [[-155 -144 -135 ... -263 -277 -258]
[ 27 103 136 ... 66 23 6]
[-316 -302 -343 ... -190 -216 -160]
...
[ 36 28 55 ... -413 -471 -650]
[-350 -274 -273 ... 50 134 117]
[ 18 166 195 ... -562 -637 -605]]
Start frame: 5423784, End frame: 5724084, Channel indices: slice(None, None, None)
Recordings after readChunk: [[ 1104 1161 1152 ... -443 -559 -487]
[ 1051 1117 1048 ... 498 410 604]
[ 1036 1074 946 ... 679 589 582]
...
[ 1545 1525 1550 ... -1155 -1054 -1017]
[ 925 777 746 ... -122 -10 -50]
[ 1615 1488 1412 ... -1309 -1252 -1211]]
Start frame: 8135751, End frame: 8436051, Channel indices: slice(None, None, None)
Recordings after readChunk: [[-1309 -1288 -1301 ... 1773 1829 1791]
[ -463 -409 -457 ... 1603 1654 1557]
[ -645 -644 -651 ... 1486 1440 1472]
...
[-1768 -1729 -1683 ... 2476 2458 2381]
[ -842 -828 -892 ... 1513 1496 1576]
[-1743 -1753 -1749 ... 2488 2332 2342]]
Start frame: 10847718, End frame: 11148018, Channel indices: slice(None, None, None)
[Errno 22] Invalid argument
Problem reading chunk from file: 'my_mda_path'
Recordings after readChunk: None
Error running mountainsort5
C:\Users\xiela\anaconda3\Lib\site-packages\spikeinterface\extractors\mdaextractors.py:422: RuntimeWarning: overflow encountered in scalar multiply

[zm711]

@yuhang21,

which version of spikeinterface do you have? Your errors aren't lining up with the current code base. Could you:

• Try updating to latest (or installing from source)
• Give us the script for how you are concatenating your recordings

Once we see this then maybe this will help. This could be that mda does lazy reading so the concatenation process doesn't correctly propagate the memmap, but I'm not sure. The immediate error it gives is the overflow in scalar multiply, but in the current code everything should be in python ints and not numpy ints so it shouldn't overflow. I want to make sure that the int problem hasn't already been fixed in the current code and then if this persists we could explore the memmap avenue.

[yuhang21]

@zm711
Thanks for getting back. Please see my codes below. It might be the .mda data structure issue since I tried to concatenate some other formats like .dat/ .rec etc. and all of them work fine.

file_folder='E:/data/20240301_152922.rec/20240301_152922.mountainsort/'
file_name_1 = '20240301_152922.group0.mda'
file_name_2 = '20240302_181122.group0.mda'
rec1=se.read_mda_recording(folder_path=file_folder, raw_fname=file_name_1, params_fname='params_1.json', geom_fname='geom_1.csv')
rec2=se.read_mda_recording(folder_path=file_folder, raw_fname=file_name_2, params_fname='params_1.json', geom_fname='geom_1.csv')
pi = read_probeinterface('util/NET-1-128_mda.json')
probe = pi.probes[0]
rec1=rec1.set_probe(probe, in_place=True)
rec2=rec2.set_probe(probe, in_place=True)
recordings_list=[]
recordings_list.append(rec1)
recordings_list.append(rec2)
rec=si.concatenate_recordings(recordings_list)

rec.set_channel_gains(0.195) # Need to manually set for now
rec.set_channel_offsets(0)

channel_ids = rec.get_channel_ids()
fs = rec.get_sampling_frequency()

After loading I used mountainsort5 as the sorter after the postprocessing

rec_filt = sp.bandpass_filter(
    rec, freq_min=300, freq_max=6000, dtype='int32')
rec_cr = sp.common_reference(
    rec_filt, operator="median", reference="global")
rec_preprocessed = sp.whiten(rec_cr, dtype='float32')
rec_for_wvf_extraction = rec_cr 
rec_for_wvf_extraction.dump(Path(save_folder) / "recording.json")


ms5_params = {
        'scheme': '2',
        'detect_sign': -1,
        'detect_threshold': 5.5,
        'npca_per_channel': 3,  # default 3
        'npca_per_subdivision': 10,  # default 10
        'snippet_mask_radius': 100, #seems necessary for the template alignment 
        'scheme2_detect_channel_radius': 200. 
        'scheme2_training_duration_sec': 300,
        'filter': False,
        'whiten': False,
        #"temporary_base_dir": "E:/temp",
        
}

#Mountainsort5 

sorting = si.run_sorter(
    sorter_name='mountainsort5',
    recording=rec_preprocessed,
    output_folder=Path(save_folder) / 'sorting',
    remove_existing_folder=True, 
    verbose=True,
    **ms5_params
)

The code will stop here by showing the above error.
Please let me know. Thank you!

[zm711]

@yuhang21, sorry one last thing. Which version of spikeinterface? I think @h-mayorquin worked on some of this stuff so maybe it is fixed in the most recent version. But if you are on the most recent version then we might have a problem.

C:\Users\xiela\anaconda3\Lib\site-packages\spikeinterface\extractors\mdaextractors.py:422: RuntimeWarning: overflow encountered in scalar multiply

This line can occur due to overflow when doing multiplication in numpy. We made some recent patches to try to prevent these overflows. The only issue is that MS5 in the most recent version has a bug on Windows that we are actively trying to patch right now.

[yuhang21]

@zm711
I used 0.99.1. I tried to update to 0.100.0 but exactly I found the bug in ms5 so I switched back. I would prefer if you can help me locate the updated mda concatenate source files and all the files in use during the concatenation. I can replace them with the most recent one and try again. Thank you

[zm711]

@yuhang21, are you doing at least python 3.10. To fix the ms5 bug is actually super easy but only works if you are past 3.10. It will actually be way easier than trying to fix any overflow errors (also I don't remember exactly where they are--because they behave differently between linux and windows).

[yuhang21]

@zm711 Thanks. Seems we're closer to the solutions. I am using 3.11 and a windows10. If you have a record on where's the error is please let me know. On my end when I used ms5 in 0.100.0 it didn't start sorting (maybe it started but super slow) so that I cannot see where're the errors. Typically my 20-min recording needs 1.2~1.5k seconds to finish the sorting but in 0.100.0 it cannot finish in 3hr. I also see some comments on other threads but want to make sure I get the correct solutions and will try it out. Thank you.

[zm711]

So in the src/spikeinterface/sorters/external/mountainsort5.py you need to edit the

with TemporaryDirectory(dir=p["temporary_base_dir"]) as tmpdir:

to

with TemporaryDirectory(dir=p["temporary_base_dir"], ignore_cleanup_errors=True) as tmpdir:

This will prevent the cached recording from failing. The change in speed is actually due to the creating the cached recording. I think part of the problem is that you need to set n_jobs_for_preprocessing otherwise it will switch to mono-job for creating the cache which will be really slow. This wasn't the case in 0.99.0.

Could you try those two changes? And let us know what happens. (Also update to 0.100.0 before making those edits).

[alejoe91]

I'm working on a PR for this :) Will push it soon

[zm711]

I'm working on a PR for this :) Will push it soon

Excellent, I've been working from a branch for all my ms5 runs for ages :)

[yuhang21]

@zm711 Hi sorry to let you know that the proposed strategies didn't work. My total recordings are about 40min in total with 30k. The console displaying the following:

runcell('start soring from here', 'E:/31-Jan-2024/scriptsprocessing_pipeline_mult_SG.py')
WARNING: Not whitening (MountainSort5 expects whitened data)

As I did the whitening and filtering in the preprocessing, we should neglect the warning here.
However, the console stuck here and didn't use any CPU. It cannot be manually stopped anyway. I believe we didn't solve ms5 bug, so I don't know whether the .mda concatenation can be solved or not. If I switch back to ms4 then I need to downgrade the isosplit. I haven't read thoroughly but is there any technical challenges if I want to edit from 0.99.1 with a new concatenation code from 0.100.5? Or it would be easier to let you guys finish the modification first in 0.100.5? Thank you (most recent is 0.100.5 as today instead of 0.100.0)

[zm711]

You're right 0.100.5 is the most recent!

I think we need to figure out what your problem is. Right now I see a few things

1. Extractor was not working (I thought that would be fixed by 0.100.5)
2. MS5 doesn't work on Windows right now (for 3.10+ the temp dir fix helps, but doesn't solve the problem completely)

It might be helpful to show us the part of the script you ran, the version of mountainsort5 you have and to rerun with verbose = True in run_sorter. It seems like you didn't even make it to the fix I suggested (for tempdir) so it would be helpful to know if this is still an extractor problem or a mountainsort problem or something else.