Kilosort
Kilosort copied to clipboard
MouseLand
Getting following error when spike sorting tetrodes: The size of tensor a must match the size of tensor b at non-singleton dimension 0
Exception has occurred: RuntimeError
The size of tensor a (16) must match the size of tensor b (2561540) at non-singleton dimension 0
File "C:\Users\m329786\code_repositories\nwb2spikesort\run_kilosort4.py", line 57, in <module>
run_kilosort(settings=settings, probe_name=probe_file, save_preprocessed_copy=True)
RuntimeError: The size of tensor a (16) must match the size of tensor b (2561540) at non-singleton dimension 0
The full output is below:
kilosort.run_kilosort: Kilosort version 4.0.36
kilosort.run_kilosort: Python version 3.9.21
kilosort.run_kilosort: ----------------------------------------
kilosort.run_kilosort: System information:
kilosort.run_kilosort: Windows-10-10.0.26100-SP0 AMD64
kilosort.run_kilosort: Intel64 Family 6 Model 170 Stepping 4, GenuineIntel
kilosort.run_kilosort: Using GPU for PyTorch computations. Specify `device` to change this.
kilosort.run_kilosort: Using CUDA device: NVIDIA RTX 2000 Ada Generation Laptop GPU 8.00GB
kilosort.run_kilosort: ----------------------------------------
kilosort.run_kilosort: Sorting [WindowsPath('C:/Users/m329786/Data/spikeinterface_output/kilosort4/sorter_output/sub-MSEL02698_ses-stim01_task-none_ieeg/sorter_output/recording.dat')]
kilosort.run_kilosort: Interpreting binary file as default dtype='int16'. If data was saved in a different format, specify `data_dtype`.
kilosort.run_kilosort:
kilosort.run_kilosort: Resource usage before sorting
kilosort.run_kilosort: ********************************************************
kilosort.run_kilosort: CPU usage: 6.70 %
kilosort.run_kilosort: Mem used: 27.20 % | 17.28 GB
kilosort.run_kilosort: Mem avail: 46.16 / 63.43 GB
kilosort.run_kilosort: ------------------------------------------------------
kilosort.run_kilosort: GPU usage: `conda install pynvml` for GPU usage
kilosort.run_kilosort: GPU memory: 13.29 % | 1.06 / 8.00 GB
kilosort.run_kilosort: Allocated: 0.00 % | 0.00 / 8.00 GB
kilosort.run_kilosort: Max alloc: 0.00 % | 0.00 / 8.00 GB
kilosort.run_kilosort: ********************************************************
kilosort.run_kilosort:
kilosort.run_kilosort: Computing preprocessing variables.
kilosort.run_kilosort: ----------------------------------------
kilosort.run_kilosort: N samples: 227219712
kilosort.run_kilosort: N seconds: 7096.715515136719
kilosort.run_kilosort: N batches: 89
c:\Users\m329786\.conda\envs\kilosort\lib\site-packages\kilosort\preprocessing.py:110: UserWarning: The use of `x.T` on tensors of dimension other than 2 to reverse their shape is deprecated and it will throw an error in a future release. Consider `x.mT` to transpose batches of matrices or `x.permute(*torch.arange(x.ndim - 1, -1, -1))` to reverse the dimensions of a tensor. (Triggered internally at C:\actions-runner\_work\pytorch\pytorch\pytorch\aten\src\ATen\native\TensorShape.cpp:4416.)
CC = CC + (X @ X.T)/X.shape[1]
kilosort.run_kilosort: Encountered error in `run_kilosort`:
Traceback (most recent call last):
File "c:\Users\m329786\.conda\envs\kilosort\lib\site-packages\kilosort\run_kilosort.py", line 269, in _sort
ops = compute_preprocessing(ops, device, tic0=tic0, file_object=file_object)
File "c:\Users\m329786\.conda\envs\kilosort\lib\site-packages\kilosort\run_kilosort.py", line 591, in compute_preprocessing
whiten_mat = preprocessing.get_whitening_matrix(bfile, xc, yc, nskip=nskip,
File "c:\Users\m329786\.conda\envs\kilosort\lib\site-packages\kilosort\preprocessing.py", line 110, in get_whitening_matrix
CC = CC + (X @ X.T)/X.shape[1]
RuntimeError: The size of tensor a (16) must match the size of tensor b (2561540) at non-singleton dimension 0
Backend qtagg is interactive backend. Turning interactive mode on.
The code that produces errors is shown below:
import json
import numpy as np
from kilosort import run_kilosort
#with open("settings.json", "r") as file:
# settings_file_data = json.load(file)
settings = dict()
settings['filename'] = r'C:\Users\m329786\Data\spikeinterface_cache\sub-MSEL02698_ses-stim01_task-none_ieeg\traces_cached_seg0.bin'
#settings['filename'] = r'C:\Users\m329786\Data\spikeinterface_output\kilosort4\sorter_output\sub-MSEL02698_ses-stim01_task-none_ieeg\sorter_output\recording.dat'
settings['n_chan_bin'] = 16 #settings_file_data['main']['n_chan_bin']
settings['Th_learned'] = 8 # Reduce this to increase the number of detected units (default=8)
settings['Th_single_ch'] = 6 # Reduce this to increase the number of detected units (default=6)
settings['Th_universal'] = 9 # Used in universal template calculations (default=9)
settings['acg_threshold'] = 0.2 # Fraction of refractory period violations used to assign "good" units (default=0.2)
settings['artifact_threshold'] = np.inf # Absolute values to zero out (default=np.inf)
settings['batch_size'] = 2561410 #settings_file_data['main']['batch_size'] # Number of samples per batch of data (default=60000)
settings['binning_depth'] = 5 # The size of the drift correction 2D histogram bin in microns (default=5)
settings['ccg_threshold'] = 0.25 # Fraction of refractory period violations that are allowed in the CCG compared to baseline; used to perform splits and merges (default=0.25)
settings['clear_cache'] = False
settings['cluster_downsampling'] = 20 # Inverse fraction of nodes used as landmarks during clustering (default=20)
settings['dmin'] = 35.0 #settings_file_data['extra']['dmin'] # Vertical spacing between channels in microns (default=None). This is adjusted automatically by Kilosort but may have to be set manually (e.g., median distance) for irregualrly spaced contacts.
settings['dminx'] = 35.0 #settings_file_data['extra']['dminx'] # Horizontal spacing between channels in microns (default=32). This is adjusted automatically by Kilosort but may have to be set manually (e.g., median distance) for irregualrly spaced contacts.
settings['do_CAR'] = True # Whether to apply common average referencing (default=True)
settings['drift_smoothing'] = [0.5, 0.5, 0.5] # Amount of gaussian smoothing to apply to the spatiotemporal drift estimation, for correlation, time (units of registration blocks), and y (units of batches) axes (default=[0.5, 0.5, 0.5])
settings['duplicate_spike_ms'] = 0.25 # Time in ms for which subsequent spikes from the same cluster are assumed to be artifacts (default=0.25). A value of 0 disables this step.
settings['fs'] = 32017.58778625954 #settings_file_data['main']['fs'] # Sampling frequency of probe (default=30000)
settings['highpass_cutoff'] = 300 # Critical frequency (Hz) for highpass Butterworth filter applied to data (default=300).
settings['invert_sign'] = False
settings['max_channel_distance'] = 35.0 #settings_file_data['extra']['max_channel_distance'] # Templates farther away than this from their nearest channel will not be used. Also limits distance between compared channels during clustering. This should be set based on the probe geometry (default=32).
settings['max_peels'] = 100 # Number of iterations to do over each batch of data in the matching pursuit step. More iterations may detect more overlapping spikes (default=100).
settings['min_template_size'] = 10 # Standard deviation of the smallest, spatial envelope Gaussian used for universal templates (default=10). This is presumably in sample points.
settings['n_pcs'] = 6 # Number of single-channel PCs to use for extracting spike features (only used if templates_from_data is True; default=6).
settings['n_templates'] = 6 # Number of single-channel templates to use for the universal templates (only used if templates_from_data is True; default=6).
settings['nblocks'] = 0 # Number of non-overlapping blocks for drift correction (default=1). For probes with fewer channels (around 64 or less) or with sparser spacing (around 50um or more between contacts), drift estimates are not likely to be accurate, so drift correction should be skipped by setting nblocks = 0.
settings['nearest_chans'] = 4 #settings_file_data['extra']['nearest_chans'] # Number of nearest channels to consider when finding local maxima during spike detection (default=10).
settings['nearest_templates'] = 16 #settings_file_data['extra']['nearest_templates'] # Number of nearest spike template locations to consider when finding local maxima during spike detection (default=100).
settings['nskip'] = 25 # Batch stride for computing whitening matrix (default=25).
settings['nt'] = 65 #settings_file_data['extra']['nt'] # The number of time samples used to represent spike waveforms (use an odd number; default=61).
settings['nt0min'] = None # Sample index for aligning waveforms, so that their minimum or maximum value happens here. Defaults to int(20 * settings['nt']/61).
settings['position_limit'] = 100 # Maximum distance (in microns) between channels that can be used to estimate spike positions in `postprocessing.compute_spike_positions` (default=100). This does not affect spike sorting, only how positions are estimated after sorting is complete.
settings['save_preprocessed_copy'] = False
settings['scale'] = None # Scaling factor to apply to data before all other operations. In most cases this should be left as None, but may be necessary for float32 data for example. If needed, `shift` and `scale` should be set such that data is roughly in the range -100 to +100.
settings['shift'] = None # Scalar shift to apply to data before all other operations. In most cases this should be left as None, but may be necessary for float32 data for example. If needed, `shift` and `scale` should be set such that data is roughly in the range -100 to +100.
settings['sig_interp'] = 20 # Approximate spatial smoothness scale in units of microns (default=20).
settings['template_sizes'] = 5 # Number of sizes for universal spike templates (multiples of the min_template_size; default=5).
settings['templates_from_data'] = True # Indicates whether spike shapes used in universal templates should be estimated from the data or loaded from the predefined templates (default=True).
settings['whitening_range'] = 16 #settings_file_data['extra']['whitening_range'] # Number of nearby channels used to estimate the whitening matrix (default=32).
settings['x_centers'] = None # Number of x-positions to use when determining center points for template groupings. If None, this will be determined automatically by finding peaks in channel density (default).
settings['data_dtype'] = 'int32'
probe_file = r'C:\Users\m329786\Data\spikeinterface_cache\sub-MSEL02698_ses-stim01_task-none_ieeg\probe_ks4.json'
ops, st, clu, tF, Wall, similar_templates, is_ref, est_contam_rate, kept_spikes = \
run_kilosort(settings=settings, probe_name=probe_file, save_preprocessed_copy=True)
Please let me know if you have any suggestions of what is causing this error. I can also provide download lonk for the binary file. The probe file is attached.
Please start by uploading kilosort4.log from the results directory.
Your probe file is improperly formatted. Each of the array entries ("xc", "yc", etc) should have shape (N,), yours have shape (N,1).
I reformatted the probe file
I managed to get a little bit further but ran into another issue
Exception has occurred: IndexError
too many indices for array: array is 1-dimensional, but 2 were indexed
File "C:\Users\m329786\code_repositories\nwb2spikesort\run_kilosort4.py", line 56, in <module>
run_kilosort(settings=settings, probe_name=probe_file, save_preprocessed_copy=True)
IndexError: too many indices for array: array is 1-dimensional, but 2 were indexed
The ks log file is also attached
kcoords is still formatted incorrectly in your probe file (same issue as before).
I fixed kcoords field and got a new error. The files are attached.
Exception has occurred: RuntimeError
The expanded size of the tensor (639) must match the existing size (1371) at non-singleton dimension 0. Target sizes: [639]. Tensor sizes: [1371]
File "C:\Users\m329786\code_repositories\nwb2spikesort\run_kilosort4.py", line 56, in <module>
run_kilosort(settings=settings, probe_name=probe_file, save_preprocessed_copy=True)
RuntimeError: The expanded size of the tensor (639) must match the existing size (1371) at non-singleton dimension 0. Target sizes: [639]. Tensor sizes: [1371]
That is strange... can you provide the link for the binary file please? If you don't want to post it here, you can e-mail it to me at [email protected]
In your e-mail you indicated the data is saved in int32 format, you need to specify that for run_kilosort as well with run_kilosort(..., data_dtype='int32'). The source of the error seems to be the large batch size, but I'm still trying to figure out why. I was able to sort it with batch_size = 60000 (the default) or batch_size = 120000 if you want to try it with those in the meantime. Anything larger resulted in the same error.
I also want to point out that I noticed a lot of brief artifacts showing up across all channels like the screenshot below, that were not removed by our default preprocessing steps. I'm assuming they're not neural because of the tight timing across all channels and the large difference in amplitude compared to the other activity. Those will interfere with the clustering step, so you will get better results if you do something to remove those from your data prior to sorting.
I do specify data_dtype='int32'. Kilosort completed with batch_size = 120000. But it seems that all detected clusters are some kind of unaligned noise.
Those artifacts you see in the data across channels are caused by electric stimulation.
The log files you attached say you're using int16, that's why I mentioned it. Make sure you're setting it as an argument to run_kilosort, it won't do anything if you're including it in the settings dictionary.
Artifacts from electric stimulation are very likely to interfere with clustering, so removing those would be my first suggestion for improving the sorting results.