Quality Metrics Tutorial

After spike sorting, you might want to validate the ‘goodness’ of the sorted units. This can be done using the qualitymetrics submodule, which computes several quality metrics of the sorted units.

import spikeinterface.core as si
import spikeinterface.extractors as se
from spikeinterface.postprocessing import compute_principal_components
from spikeinterface.qualitymetrics import (compute_snrs, compute_firing_rates,
    compute_isi_violations, calculate_pc_metrics, compute_quality_metrics)

First, let’s download a simulated dataset from the repo ‘https://gin.g-node.org/NeuralEnsemble/ephy_testing_data’

local_path = si.download_dataset(remote_path='mearec/mearec_test_10s.h5')
recording, sorting = se.read_mearec(local_path)
print(recording)
print(sorting)

MEArecRecordingExtractor: 32 channels - 32.0kHz - 1 segments - 320,000 samples - 10.00s
                          float32 dtype - 39.06 MiB
  file_path: /home/docs/spikeinterface_datasets/ephy_testing_data/mearec/mearec_test_10s.h5
MEArecSortingExtractor: 10 units - 1 segments - 32.0kHz
  file_path: /home/docs/spikeinterface_datasets/ephy_testing_data/mearec/mearec_test_10s.h5

Extract spike waveforms

For convenience, metrics are computed on the WaveformExtractor object, because it contains a reference to the “Recording” and the “Sorting” objects:

we = si.extract_waveforms(recording=recording,
                          sorting=sorting,
                          folder='waveforms_mearec',
                          sparse=False,
                          ms_before=1,
                          ms_after=2.,
                          max_spikes_per_unit=500,
                          n_jobs=1,
                          chunk_durations='1s')
print(we)

extract waveforms memmap multi buffer:   0%|          | 0/10 [00:00<?, ?it/s]
extract waveforms memmap multi buffer: 100%|##########| 10/10 [00:00<00:00, 150.01it/s]
WaveformExtractor: 32 channels - 10 units - 1 segments
  before:32 after:64 n_per_units:500

The spikeinterface.qualitymetrics submodule has a set of functions that allow users to compute metrics in a compact and easy way. To compute a single metric, one can simply run one of the quality metric functions as shown below. Each function has a variety of adjustable parameters that can be tuned.

firing_rates = compute_firing_rates(we)
print(firing_rates)
isi_violation_ratio, isi_violations_count = compute_isi_violations(we)
print(isi_violation_ratio)
snrs = compute_snrs(we)
print(snrs)

{'#0': 5.3, '#1': 5.0, '#2': 4.3, '#3': 3.0, '#4': 4.8, '#5': 3.7, '#6': 5.1, '#7': 11.1, '#8': 19.5, '#9': 12.9}
{'#0': 0.0, '#1': 0.0, '#2': 0.0, '#3': 0.0, '#4': 0.0, '#5': 0.0, '#6': 0.0, '#7': 0.0, '#8': 0.0, '#9': 0.0}
{'#0': 23.631605, '#1': 25.505018, '#2': 13.756625, '#3': 21.838345, '#4': 7.4441104, '#5': 7.4596915, '#6': 20.874245, '#7': 7.4103203, '#8': 8.090734, '#9': 8.939655}

Some metrics are based on the principal component scores, so they require a WaveformsPrincipalComponent object as input:

pc = compute_principal_components(waveform_extractor=we, load_if_exists=True,
                                     n_components=3, mode='by_channel_local')
print(pc)

pc_metrics = calculate_pc_metrics(pc, metric_names=['nearest_neighbor'])
print(pc_metrics)

Fitting PCA:   0%|          | 0/10 [00:00<?, ?it/s]
Fitting PCA:  20%|██        | 2/10 [00:00<00:00, 17.16it/s]
Fitting PCA:  40%|████      | 4/10 [00:00<00:00, 18.60it/s]
Fitting PCA:  60%|██████    | 6/10 [00:00<00:00, 18.32it/s]
Fitting PCA:  80%|████████  | 8/10 [00:00<00:00, 14.98it/s]
Fitting PCA: 100%|██████████| 10/10 [00:00<00:00, 10.93it/s]
Fitting PCA: 100%|██████████| 10/10 [00:00<00:00, 12.97it/s]

Projecting waveforms:   0%|          | 0/10 [00:00<?, ?it/s]
Projecting waveforms: 100%|██████████| 10/10 [00:00<00:00, 116.11it/s]
WaveformPrincipalComponent: 32 channels - 1 segments
  mode: by_channel_local n_components: 3
{'nn_hit_rate': {'#0': 0.9952830188679245, '#1': 0.96, '#2': 0.9244186046511628, '#3': 0.9916666666666667, '#4': 0.9739583333333334, '#5': 0.972972972972973, '#6': 0.9803921568627451, '#7': 0.9363636363636364, '#8': 0.9690721649484536, '#9': 0.9437984496124031}, 'nn_miss_rate': {'#0': 0.0036127167630057803, '#1': 0.002158273381294964, '#2': 0.004273504273504274, '#3': 0.0, '#4': 0.0010760401721664275, '#5': 0.0014124293785310734, '#6': 0.00036023054755043225, '#7': 0.007874015748031496, '#8': 0.017241379310344827, '#9': 0.009334415584415584}}

To compute more than one metric at once, we can use the compute_quality_metrics function and indicate which metrics we want to compute. This will return a pandas dataframe:

metrics = compute_quality_metrics(we)
print(metrics)

/home/docs/checkouts/readthedocs.org/user_builds/spikeinterface/checkouts/stable/src/spikeinterface/qualitymetrics/misc_metrics.py:848: UserWarning: Units ['#0', '#1', '#2', '#3', '#4', '#5', '#6', '#7', '#8', '#9'] have too few spikes and amplitude_cutoff is set to NaN
  warnings.warn(f"Units {nan_units} have too few spikes and " "amplitude_cutoff is set to NaN")
/home/docs/checkouts/readthedocs.org/user_builds/spikeinterface/checkouts/stable/src/spikeinterface/qualitymetrics/misc_metrics.py:702: UserWarning:
  warnings.warn("")
/home/docs/checkouts/readthedocs.org/user_builds/spikeinterface/checkouts/stable/src/spikeinterface/qualitymetrics/misc_metrics.py:145: UserWarning: Bin duration of 60s is larger than recording duration. Presence ratios are set to NaN.
  warnings.warn(
/home/docs/checkouts/readthedocs.org/user_builds/spikeinterface/checkouts/stable/src/spikeinterface/qualitymetrics/misc_metrics.py:1429: UserWarning: The `sd_ratio` metric require the `spike_amplitudes` waveform extension. Use the `postprocessing.compute_spike_amplitudes()` functions. SD ratio metric will be set to NaN
  warnings.warn(

Computing PCA metrics:   0%|          | 0/10 [00:00<?, ?it/s]
Computing PCA metrics:  10%|█         | 1/10 [00:00<00:01,  7.25it/s]
Computing PCA metrics:  20%|██        | 2/10 [00:00<00:01,  7.68it/s]
Computing PCA metrics:  30%|███       | 3/10 [00:00<00:00,  7.05it/s]
Computing PCA metrics:  40%|████      | 4/10 [00:00<00:01,  5.66it/s]
Computing PCA metrics:  50%|█████     | 5/10 [00:00<00:00,  6.00it/s]
Computing PCA metrics:  60%|██████    | 6/10 [00:00<00:00,  6.58it/s]
Computing PCA metrics:  70%|███████   | 7/10 [00:01<00:00,  5.61it/s]
Computing PCA metrics:  80%|████████  | 8/10 [00:01<00:00,  5.17it/s]
Computing PCA metrics:  90%|█████████ | 9/10 [00:01<00:00,  4.75it/s]
Computing PCA metrics: 100%|██████████| 10/10 [00:01<00:00,  4.64it/s]
Computing PCA metrics: 100%|██████████| 10/10 [00:01<00:00,  5.39it/s]
    amplitude_cutoff  amplitude_cv  ...  nn_hit_rate  nn_miss_rate
#0               NaN           NaN  ...     0.995283      0.003613
#1               NaN           NaN  ...     0.960000      0.002158
#2               NaN           NaN  ...     0.924419      0.004274
#3               NaN           NaN  ...     0.991667      0.000000
#4               NaN           NaN  ...     0.973958      0.001076
#5               NaN           NaN  ...     0.972973      0.001412
#6               NaN           NaN  ...     0.980392      0.000360
#7               NaN           NaN  ...     0.936364      0.007874
#8               NaN           NaN  ...     0.969072      0.017241
#9               NaN           NaN  ...     0.943798      0.009334

[10 rows x 23 columns]