Plot_unit_locations error?

[JazlynnXM-Tan]

Hi,

I have a recording on neuropixel 2.0 using openephys and performed sorting on kilosort 4 (outside of spikeinterface). I then loaded both the recording and sorting into spikeinterface to create a sorting_analyzer object. After computing random spikes, templates, spike and unit locations, I tried plot_unit_locations which gave a plot where most of the units seemed to be on the rightmost shank, while there are also many units on the third from left shank, where there aren't even recording channels selected. Finally, some units are plotted in between shanks. These locations do not match the spike locations visualization from kilosort either.

Plot_unit_locations output:

I plotted the spike_positions.npy file from kilosort 4:

To look at a single unit in detail, here is unit 52 highlighted in red on the second from left shank (for ease of viewing I did not plot all the spikes from other units). The plot_unit_locations widget however highlights unit 52 on the 4th shank, at a different depth as well.

May I know if I am doing something wrong?

This is my code:

recording = si.read_openephys(recording_file_path)
kilosort_results = si.read_kilosort(kilosort_file_path)
analyzer = si.create_sorting_analyzer(sorting=kilosort_results, recording=recording,format="binary_folder",return_scaled=True, folder='my_save_path')

si.compute_random_spikes(analyzer)
si.compute_templates(analyzer)
si.compute_spike_locations(analyzer)
si.compute_unit_locations(analyzer)

si.plot_unit_locations(analyzer,backend='ipywidgets')

[chrishalcrow]

Hi @JazlynnXM-Tan ,

The unit_locations are likely computed in different ways. By default, spikeinterface uses monopolar_triangulation (when I googled, the first results was page 4 here: https://www.eneuro.org/content/eneuro/11/2/ENEURO.0229-23.2023.full.pdf ). This assumes the cell is monopolar current source, and will give some “interesting” results if the unit template doesn’t look like a "standard" action potential.

I have a feeling you’ll get closer to kilosort results if you compute again using center_of_mass. You can do this using

si.compute_spike_locations(analyzer, method=center_of_mass")

Could you recompute the spike locations, and let us know what the results look like? Thanks!

We’ve actually documented this wrong - our docs claim that center_of_mass is the default, but it isn’t! I'll make a fix now.

[JazlynnXM-Tan]

Hi, I reran it as you said. The units no longer appear on regions without recording channels. However, they do still seem to be clustered on the rightmost shank. The position of the unit 52 has changed but is also not in the expected location from kilosort results.

From plot_unit_locations

From plot_spike_locations

I did get this error from computing either spike and unit locations though:
spikeinterface/sortingcomponents/peak_localization.py:226: RuntimeWarning: invalid value encountered in divide coms = np.dot(wf_data, local_contact_locations) / (np.sum(wf_data, axis=1)[:, np.newaxis])

[JazlynnXM-Tan]

I ran about 15 min of the same recording in spikeinterface with kilosort4 in run_sorter_by_property, grouping by shank. The output from plot_unit_locations look normal, using the center of mass method of computing locations. Could the problem I had when using kilosort externally be related to how the kilosort results were loaded into spikeinterface (perhaps at the read_kilosort step)? Or maybe it was due to the improper separation of the shanks in kilosort? I am using the latest version of kilosort (both when calling from spikeinterface and running it separately) and my probe map includes kcoords information which the kilosort documentation says should restrict template placement to within the shank.

[alejoe91]

HI @JazlynnXM-Tan

This is helpful! Can you ask on the Kilosort GH repo as well? It might be indeed some problem with kcoords! Spike sorting by shank in SpikeInterface is bullet-proof :)

[JazlynnXM-Tan]

Sure, have asked on kilosort: MouseLand/Kilosort#826

[alejoe91]

@JazlynnXM-Tan actually the problem is probably in the channel map from SpikeGLX: see this issue int-brain-lab/ibl-sorter#24

Thansk @oliche for the hint!

[JazlynnXM-Tan]

Hi @alejoe91, do you mean that the raw data itself has to be ordered by channel number (ie. CH1, CH2, etc.)? I was under the impression that the channel map file contains the information for the sorters to reorder the data. The channel maps I used did have the x, y and k coordinate information.

Do you have any suggestions on what we can do?

Thanks!

[alejoe91]

I think that performing sorting within spikeinterface might fix the problem. Otherwise I'll ask @oliche what recommendation he has!

[oliche]

I would be curious to look at the spikeglx meta data file and if possible a few seconds of raw data to understand the issue. @grg2rsr would also be interested !