Batch buffer overlap error causing autocorrelogram peaks & double counted spikes ( ops.ntbuff )

[czuba]

Bug in ops.ntbuff implementation

It appears the batch buffer count (ops.ntbuff) is getting improperly applied/trimmed somewhere, causing batches to overlap in time (either in fread content, or assignment to st3).

Simple test case

Adjust advanced parameters in the gui to increase .ntbuff by a large amount (e.g. ntbuff==64*300), then use Phy to examine autocorrelogram (ACG) peaks in the results to results using the default (ntbuff==64).

The following example code will retrieve a pointer to the GUI, increase the buffer length .ntbuff, and adjust the .NT accordingly (default config states special constraints on increments and value relative to .ntbuff), then apply updated settings to the ksGUI object:

% ksGUI_updatePars.m
% script to update GUI pars with recommend settings for stereo Uprobes (32 chan, 50/100 um geometry)
% 
fprintf(['\n',repmat('=-',1,20),'\n']);

%% Retrieve gui settings
ks = get(gcf, 'UserData') %#ok<NOPTS>
ops = ks.ops %#ok<NOPTS>
fprintf([repmat('--',1,20),'\n']);

%% Apply changes
ops.ntbuff              = 64*300; % (def=64) samples of symmetrical buffer for whitening and spike detection  
% TBC version define batches in seconds of data, not abstract bit chunks
    batchSec                = 5;  % define batch number in seconds of data (TBC: 5 seems good for 1-2 hr files)
    ops.NT                  = ceil(batchSec*ops.fs/32)*32 +ops.ntbuff;   % must be multiple of 32 + ntbuff. This is the batch size (try decreasing if out of memory).
 
%% Apply to updates to gui params object
fprintf('New ops settings:\n')
disp(ops)
ks.ops = ops;
fprintf(['\nDone.\tNew settings have been applied to kilosort GUI object [ks].\n',repmat('=-',1,20),'\n']);

! NOTE: ! You'll need to run this on a freshly Reset GUI instance, and be sure these settings have been applied prior to Preprocessing.
...multiple copies/sources of the ops structure exist after either Preprocess or Spikesort stages have been initiated in the GUI, and intermediate-stage changes to ks.ops do not [reliably] to propagate.

It seems like there are multiple locations along the processing stream where this hiccup might be occurring, so I don't have a ready-made fix to create a PR.

• trackAndSort.m, line 182 seems like a possible candidate. To my eye, it looks to be applying an .ntbuff offset to data that was read based on .NT-spaced samples. But that would seem to be a mismatch on the indexing side of things, rather than an overlapping readout type error (...or I may just be misreading that .ntbuff application all together)

Background

I've been working to fix/reduce the occurrence of temporally stuttered clusters in our Kilosort2 outputs using 32 channel linear array recordings (Plexon uProbes: stereo geometry, 16 stereo-pairs, 50 µm within pairs [x], 100 µm between pairs [y]).
Because our recordings have significantly fewer channels/samples (relative to neuropixels) within a given batch, I was hoping that increasing the buffer might help prevent whatever idiosyncrasies Kilosort seemed to be stumbling over between batches.

Buffer adjustments don't seem to have been a magic bullet so far (disabling batch reordering has helped the most, but there's more to be done...), but this bug is something that would affect all Kilosort users' results, so it seems like something worth flagging independent of our particular application issues.

Thanks!

[biomath]

This is likely related to my problem here:
#223

[marius10p]

hi @czuba, sorry I missed this. Why are you trying to change ntbuff? There is no telling what could happen. It serves as buffer in two different places in fact, and it's a little tricky to make sense of it.

[biomath]

Hi @marius10p,
As @czuba hasn't answered yet, do you think this is related to the issue I raised here?#223
Thank you!

[czuba]

Hey @marius10p,
Yes, I was working on the temporal stuttering issue that has been a persistent artifact on my KS2 sorting results. In addition to my hunch above (about significantly fewer samples/batch in 32ch recordings) my hope/intuition for lengthening ops.ntbuff was that if something about the signal within adjacent batches is causing templates from the same unit to appear incongruent to the template matching algorithm, then incorporating a more significant chunk of shared buffer between batches might help.

@biomath, I too fell into the head-scratching hole of trying recreate the filtered .dat temp file (both for understanding how it was being sampled, and to weave in a faster memmapped read).

More background & examples

Testing a longer ops.ntbuff with all other parameters set to default values (& batch reordering enabled), I'm still seeing lots of temporal stuttering:

...granted this is not the most rockstar stable/high-amplitude unit, but it should certainly be recoverable.

Disabling batch reordering & optimizing a few parameters (incl. 10 sec batches) helped reduce the intense stuttering, but this example unit was still getting chopped up into a set of bizarrely interdigitated units:

Both of these examples also had really strong CCG peaks, which I suspect is due to buffer segments that are not getting fully trimmed/accounted for during the actual spike detection pass. Since the default ops.ntbuff is only 64 samples, and the default waveform length is 61 samples, there is very little opportunity for this kind of double counting (...but still non-zero).

After reverting to the standard ops.ntbuff — mostly because the buffer overlap CCG peaks made that solution a nonstarter — and further tuning parameters to our recordings I've managed dial in the sorting results into usable shape:

Note: this is the output before any manual curation.

If the above interpretation of ops.ntbuff artifacts holds up, it seems like a bug worth really tracking down, but I the depths of it (extending into the structure of filtered dat files & sampling/handling of batch data w/in cuda code) were beyond me.

---

I'll split some of the above into a separate "How I tuned parameters for Uprobe recordings" post. In brief, the key parameter adjustments for Uprobe stereotrode recordings seemed to be:

• disabling batch reordering
• increasing batch duration from the 2 second default to 5-10 seconds
• increasing sigma to 70 (to better allow units to span the 100µm distance between stereopairs)
• refined code to sweep through entire timecourse during template fitting (see also Batch ordering in initial pass misses half the sorted batches? #84 )

I still see a handful of nonsensical unit dropouts & splits, but they're far fewer and considerably more manageable than the extreme batch-wise stuttering apparent in my first example (see also #60 & #63 ).

Overall, I still don't understand why I'm seeing sooo much of this temporal stuttering artifact...it obviously isn't affecting everyone (at least not to this degree), otherwise it'd be 'issue no.1' across-the-board, but the root of it does seem to crop up in various other issues/forms (e.g. #175 "Sudden change of the template over consecutive batches").