BUG: clusters with <=10 spikes export ContamPct = 0.0, the best possible score

[cindykrafft]

Describe the issue:

CCG.refract initializes the per-cluster contamination estimate to zeros and only computes it for clusters with more than 10 spikes and a nonzero time span:

R12 = np.zeros(Nfilt, )
for kk in range(Nfilt):
    ...
    if (len(st1) > 10) and ((st1.max() - st1.min()) != 0):
        is_refractory[kk], cross_refractory[kk], R12[kk] = check_CCG(...)

io.save_to_phy then writes est_contam_rate * 100 to cluster_ContamPct.tsv. Any cluster with ≤10 spikes is therefore exported with ContamPct = 0.0 — the best possible score — indistinguishable from a genuinely clean unit. (Its KSLabel is mua, but ContamPct is commonly used as an independent selection criterion.)

Kilosort 2.5/3 handled this case the opposite way: set_cutoff.m defaulted est_contam_rate to 1, so unevaluable units exported ContamPct = 100.

Why it matters:

Selection criteria of the form "ContamPct < x" (with x = 10 or 20, both common in published Neuropixels analyses) silently admit every tiny junk cluster unless the pipeline also applies a minimum-spike-count criterion.

Reproduction (kilosort 4.1.7):

import numpy as np
from kilosort import CCG

rng = np.random.default_rng(1)
st_good = np.sort(rng.uniform(0, 600, 5000))  # Poisson unit, no refractoriness
st_junk = np.sort(rng.uniform(0, 600, 8))     # 8-spike junk cluster
st = np.concatenate([st_good, st_junk])
clu = np.concatenate([np.zeros(5000, int), np.ones(8, int)]).astype(np.int32)
o = np.argsort(st)
is_ref, est = CCG.refract(clu[o], st[o])
print(est * 100)   # -> [93.9  0.0]

The 5000-spike contaminated unit is correctly reported at 93.9; the 8-spike junk cluster gets 0.0.

Suggested fix:

Initialize R12 = np.ones(Nfilt,) so unevaluable clusters export as ContamPct = 100, restoring the KS2.5/3 convention (fix PR incoming from this account). Exporting NaN instead would also work if you prefer "not evaluated" to be visibly distinct from "fully contaminated" — happy to adjust the PR either way.

Version: reproduced on kilosort 4.1.7 (pip), torch 2.14 CPU; the code path is unchanged on current main.

---

Generated by Claude Code