Skip to content

Commit 0b0888a

Browse files
authored
Merge branch 'main' into fix_merge_prop
2 parents 202c068 + 7a4e8eb commit 0b0888a

13 files changed

Lines changed: 416 additions & 11 deletions

File tree

doc/conf.py

Lines changed: 1 addition & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -129,6 +129,7 @@
129129
'../examples/tutorials/qualitymetrics',
130130
'../examples/tutorials/comparison',
131131
'../examples/tutorials/widgets',
132+
'../examples/tutorials/forhowto',
132133
]),
133134
'within_subsection_order': FileNameSortKey,
134135
'ignore_pattern': '/generate_*',

doc/how_to/index.rst

Lines changed: 1 addition & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -20,3 +20,4 @@ Guides on how to solve specific, short problems in SpikeInterface. Learn how to.
2020
physical_units
2121
unsigned_to_signed
2222
customize_a_plot
23+
../tutorials/forhowto/plot_1_working_with_tetrodes

doc/how_to/process_by_channel_group.rst

Lines changed: 2 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -1,3 +1,5 @@
1+
.. _recording-by-channel-group:
2+
13
Process a recording by channel group
24
====================================
35

27.3 KB
Loading
Lines changed: 106 additions & 8 deletions
Original file line numberDiff line numberDiff line change
@@ -1,30 +1,128 @@
1-
Noise cutoff (not currently implemented)
2-
========================================
1+
Noise cutoff (:code:`noise_cutoff`)
2+
===================================
33

44
Calculation
55
-----------
66

77

8-
Metric describing whether an amplitude distribution is cut off, similar to _amp_cutoff :ref:`amplitude cutoff <amp_cutoff>` but without a Gaussian assumption.
9-
A histogram of amplitudes is created and quantifies the distance between the low tail, mean number of spikes and high tail in terms of standard deviations.
8+
Metric describing whether an amplitude distribution is cut off as it approaches zero, similar to :ref:`amplitude cutoff <amp_cutoff>` but without a Gaussian assumption.
109

11-
A SpikeInterface implementation is not yet available.
10+
The **noise cutoff** metric assesses whether a unit’s spike‐amplitude distribution is truncated
11+
at the low-end, which may be due to the high amplitude detection threshold in the deconvolution step,
12+
i.e., if low‐amplitude spikes were missed. It does not assume a Gaussian shape;
13+
instead, it directly compares counts in the low‐amplitude bins to counts in high‐amplitude bins.
14+
15+
1. **Build a histogram**
16+
17+
For each unit, divide all amplitudes into ``n_bins`` equally spaced bins over the range of the amplitude.
18+
If the number of spikes is large, you may consider using a larger ``n_bins``. For a small number of spikes, consider a smaller ``n_bins``.
19+
Let :math:`n_i` denote the count in the :math:`i`-th bin.
20+
21+
2. **Identify the “low” region**
22+
- Compute the amplitude value at the specified ``low_quantile`` (for example, 0.10 = 10th percentile), denoted as :math:`\text{amp}_{low}`.
23+
- Find all histogram bins whose upper edge is below that quantile value. These bins form the "low‐quantile region".
24+
- Compute
25+
26+
.. math::
27+
L_{\mathrm{bins}} = \bigl\{i : \text{upper_edge}_i \le \text{amp}_{low}\bigr\}, \quad
28+
\mu_{\mathrm{low}} = \frac{1}{|L_{\mathrm{bins}}|}\sum_{i\in L_{\mathrm{bins}}} n_i.
29+
30+
3. **Identify the “high” region**
31+
32+
- Compute the amplitude value at the specified ``high_quantile`` (for example, 0.25 = top 25th percentile), denoted as :math:`\text{amp}_{high}`.
33+
- Find all histogram bins whose lower edge is greater than that quantile value. These bins form the "high‐quantile region".
34+
- Compute
35+
36+
.. math::
37+
H_{\mathrm{bins}} &= \bigl\{i : \text{lower_edge}_i \ge \text{amp}_{high}\bigr\}, \\
38+
\mu_{\mathrm{high}} &= \frac{1}{|H_{\mathrm{bins}}|}\sum_{i\in H_{\mathrm{bins}}} n_i, \quad
39+
\sigma_{\mathrm{high}} = \sqrt{\frac{1}{|H_{\mathrm{bins}}|-1}\sum_{i\in H_{\mathrm{bins}}}\bigl(n_i-\mu_{\mathrm{high}} \bigr)^2}.
40+
41+
4. **Compute cutoff**
42+
43+
The *cutoff* is given by how many standard deviations away the low-amplitude bins are from the high-amplitude bins, defined as
44+
45+
.. math::
46+
\mathrm{cutoff} = \frac{\mu_{\mathrm{low}} - \mu_{\mathrm{high}}}{\sigma_{\mathrm{high}}}.
47+
48+
49+
- If no low‐quantile bins exist, a warning is issued and ``cutoff = NaN``.
50+
- If no high‐quantile bins exist or :math:`\sigma_{\mathrm{high}} = 0`, a warning is issued and ``cutoff = NaN``.
51+
52+
5. **Compute the low-to-peak ratio**
53+
54+
- Let :math:`M = \max_i\,n_i` be the height of the largest bin in the histogram.
55+
- Define
56+
57+
.. math::
58+
\mathrm{ratio} = \frac{\mu_{\mathrm{low}}}{M}.
59+
60+
61+
- If there are no low bins, :math:`\mathrm{ratio} = NaN`.
62+
63+
64+
Together, ``(cutoff, ratio)`` quantify how suppressed the low‐end of the amplitude distribution is relative to the top quantile and to the peak.
1265

1366
Expectation and use
1467
-------------------
1568

1669
Noise cutoff attempts to describe whether an amplitude distribution is cut off.
17-
The metric is loosely based on [Hill]_'s amplitude cutoff, but is here adapted (originally by [IBL]_) to avoid making the Gaussianity assumption on spike distributions.
18-
Noise cutoff provides an estimate of false negative rate, so a lower value indicates fewer missed spikes (a more complete unit).
70+
Larger values of ``cutoff`` and ``ratio`` suggest that the distribution is cut-off.
71+
IBL uses the default value of 1 (equivalent to e.g. ``low_quantile=0.01, n_bins=100``) to choose the number of
72+
lower bins, with a suggested threshold of 5 for ``cutoff`` to determine whether a unit is cut off or not.
73+
In practice, the IBL threshold is quite conservative, and a lower threshold might work better for your data.
74+
We suggest plotting the data using the :py:func:`~spikeinterface.widgets.plot_amplitudes` widget to view your data when choosing your threshold.
75+
It is suggested to use this metric when the amplitude histogram is **unimodal**.
76+
77+
The metric is loosely based on [Hill]_'s amplitude cutoff, but is here adapted (originally by [IBL2024]_) to avoid making the Gaussian assumption on spike distributions.
78+
79+
Example code
80+
------------
81+
82+
.. code-block:: python
83+
84+
import numpy as np
85+
import matplotlib.pyplot as plt
86+
from spikeinterface.full as si
87+
88+
# Suppose `sorting_analyzer` has been computed with spike amplitudes:
89+
# Select units you are interested in visualizing
90+
unit_ids = ...
91+
92+
# Compute noise_cutoff:
93+
summary_dict = compute_noise_cutoff(
94+
sorting_analyzer=sorting_analyzer
95+
high_quantile=0.25,
96+
low_quantile=0.10,
97+
n_bins=100,
98+
unit_ids=unit_ids
99+
)
100+
101+
Reference
102+
---------
103+
104+
.. autofunction:: spikeinterface.qualitymetrics.misc_metrics.compute_noise_cutoff
105+
106+
Examples with plots
107+
-------------------
108+
109+
Here is shown the histogram of two units, with the vertical lines separating low- and high-amplitude regions.
110+
111+
- On the left, we have a unit with no truncation at the left end, and the cutoff and ratio are small.
112+
- On the right, we have a unit with truncation at -1, and the cutoff and ratio are much larger.
19113

114+
.. image:: example_cutoff.png
115+
:width: 600
20116

21117
Links to original implementations
22118
---------------------------------
23119

24120
* From `IBL implementation <https://github.com/int-brain-lab/ibllib/blob/2e1f91c622ba8dbd04fc53946c185c99451ce5d6/brainbox/metrics/single_units.py>`_
25121

122+
Note: Compared to the original implementation, we have added a comparison between the low-amplitude bins to the largest bin (``noise_ratio``).
123+
The selection of low-amplitude bins is based on the ``low_quantile`` rather than the number of bins.
26124

27125
Literature
28126
----------
29127

30-
Metric introduced by [IBL]_ (adapted from [Hill]_'s amplitude cutoff metric).
128+
Metric introduced by [IBL2024]_.

doc/references.rst

Lines changed: 2 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -123,6 +123,8 @@ References
123123
124124
.. [IBL] `Spike sorting pipeline for the International Brain Laboratory. 2022. <https://figshare.com/articles/online_resource/Spike_sorting_pipeline_for_the_International_Brain_Laboratory/19705522/3>`_
125125
126+
.. [IBL2024] `Spike sorting pipeline for the International Brain Laboratory - Version 2. 2024. <https://figshare.com/articles/online_resource/Spike_sorting_pipeline_for_the_International_Brain_Laboratory/19705522?file=49783080>`_
127+
126128
.. [Jackson] `Quantitative assessment of extracellular multichannel recording quality using measures of cluster separation. Society of Neuroscience Abstract. 2005. <https://www.sciencedirect.com/science/article/abs/pii/S0306452204008425>`_
127129
128130
.. [Jain] `UnitRefine: A Community Toolbox for Automated Spike Sorting Curation. 2025 <https://www.biorxiv.org/content/10.1101/2025.03.30.645770v1>`_

examples/how_to/working_with_tetrodes.py

Whitespace-only changes.
Lines changed: 4 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,4 @@
1+
For how to
2+
----------
3+
4+
These documents are files which we would like sphinx-gallery to run, but are linked by the how to guide.
Lines changed: 117 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,117 @@
1+
"""
2+
Working with tetrodes
3+
=====================
4+
5+
Tetrodes are a common recording method for electrophysiological data. It is also common
6+
to record from several tetrodes at the same time. In this 'how to' we'll see how to
7+
work with data from two tetrodes, each with four channels.
8+
9+
We'll start by importing some functions we'll use in this How To guide
10+
"""
11+
12+
import spikeinterface.preprocessing as spre
13+
from spikeinterface.widgets import plot_traces, plot_probe_map
14+
from spikeinterface import generate_ground_truth_recording
15+
16+
from probeinterface import generate_tetrode, ProbeGroup
17+
18+
##############################################################################
19+
# In practice, you would read in your raw data from a file. Instead, we will generate a
20+
# recording with eight channels. We can also set a duration, number of units and
21+
# sampling frequency.
22+
23+
recording, _ = generate_ground_truth_recording(
24+
durations = [60], # make the recording 60s long
25+
sampling_frequency=30_000,
26+
num_channels=8,
27+
num_units=10,
28+
)
29+
30+
##############################################################################
31+
# We now need to define the probe. This will tell the recording which channels came from
32+
# which tetrode. To do this, we will use the :code:`generate_tetrode` function from :code:`ProbeInterface`
33+
# to generate two 4-channel probes (representing one tetrode each). In our case, since we
34+
# don't know the relative distances between the tetrodes, we will move the second
35+
# tetrode away from the first by 100 microns. This is just so we can visualize the
36+
# results more easily. Eventually, we will sort each tetrode separately, so their
37+
# relative distance won't affect the results.
38+
39+
# Technically, we will add each tetrode to a :code:`ProbeGroup`. Read more in the ProbeInterface
40+
# docs.
41+
42+
# Create each individual tetrode
43+
tetrode_1 = generate_tetrode()
44+
tetrode_1.create_auto_shape()
45+
46+
tetrode_2 = generate_tetrode()
47+
tetrode_2.move([100, 0])
48+
tetrode_2.create_auto_shape()
49+
50+
# Add the two tetrodes to a ProbeGroup
51+
tetrode_group = ProbeGroup()
52+
tetrode_group.add_probe(tetrode_1)
53+
tetrode_group.add_probe(tetrode_2)
54+
55+
# Now we need to "wire" our tetrodes to ensure that each contact
56+
# can be associated with the correct channel when we attach it
57+
# to the recording. In this example we are just using `range`
58+
# but see ProbeInterface for more tutorials on wiring
59+
tetrode_group.set_global_device_channel_indices(range(8))
60+
61+
##############################################################################
62+
# We can now attach the :code:`tetrode_group` to our recording. To check if this worked, we'll
63+
# plot the probe map
64+
65+
recording_with_probe = recording.set_probegroup(tetrode_group)
66+
plot_probe_map(recording_with_probe)
67+
68+
##############################################################################
69+
# Looks good! Now that the recording is aware of the probe geometry, we can
70+
# begin a standard spike sorting pipeline. First, we can apply preprocessing.
71+
# Note that we apply this preprocessing on the entire bundle of tetrodes.
72+
73+
preprocessed_recording = spre.bandpass_filter(recording_with_probe)
74+
75+
##############################################################################
76+
# WARNING: a very common preprocessing step is to apply a common median
77+
# reference. This subtracts the median signal from all channels to help
78+
# remove noise. However, for a tetrode, a spike is often seen on all
79+
# channels. So removing the median can remove the entire spike!
80+
# This is still a danger if you have two tetrodes in a bundle, which
81+
# might pick up the same spike, but becomes less dangerous
82+
# as the number of tetrodes in your bundle increases.
83+
#
84+
# Tetrodes often have dead channels, so it is advised to try and detect
85+
# and remove these. For tetrodes, we should use a detection method which
86+
# doesn't depend on the channel locations such as std or mad:
87+
88+
recording_good_channels = spre.detect_and_remove_bad_channels(
89+
preprocessed_recording,
90+
method = "std",
91+
)
92+
93+
##############################################################################
94+
# It can be a good idea to sort your tetrode data separately for each tetrode.
95+
# When we use :code:`set_probegroup`, the channels are automatically
96+
# labelled by which probe in the probe group they belong to. We can access
97+
# this labeling using the "group" property.
98+
99+
print(recording_good_channels.get_property("group"))
100+
101+
##############################################################################
102+
# We can then use this information to split the recording by the group property:
103+
104+
grouped_recordings = recording_good_channels.split_by('group')
105+
print(grouped_recordings)
106+
107+
##############################################################################
108+
# Now that we've got preprocess, clean data. Let's take a look at a
109+
# snippet of data from the first group:
110+
111+
plot_traces(grouped_recordings[0])
112+
113+
##############################################################################
114+
# Beautiful! We are now ready to sort. To read more about sorting by group, see
115+
# :ref:`sorting-by-channel-group`. Note that many modern sorters are designed
116+
# to sort data from high-density probes and will fail for tetrodes. Please read
117+
# each spike sorter's documentation to find out if it is appropriate for tetrodes.

pyproject.toml

Lines changed: 1 addition & 3 deletions
Original file line numberDiff line numberDiff line change
@@ -202,9 +202,7 @@ docs = [
202202
"networkx",
203203
"skops", # For automated curation
204204
"scikit-learn", # For automated curation
205-
# Download data
206-
"pooch>=1.8.2",
207-
"datalad>=1.0.2",
205+
"huggingface_hub", # For automated curation
208206

209207
# for release we need pypi, so this needs to be commented
210208
"probeinterface @ git+https://github.com/SpikeInterface/probeinterface.git", # We always build from the latest version

0 commit comments

Comments
 (0)