Skip to content

03. Automatic Spike Sort

Jump to bottom
Xin Niu edited this page Dec 19, 2024 · 14 revisions

Spike sorting comprises spike detection, spike code calculation, and spike clustering. Before spike detection, the median value of the raw signal for each bundle (8 channels) is subtracted. Spike detection creates *_spikes.mat, spike code calculation creates '*_spikeCodes.mat', and spike clustering creates times_*.mat. This is only for microchannels, as macro ones do not detect spikes due to low sampling frequency (and perhaps the nature of electrodes).

As the input CSC files may have multiple segments (_<suffix>), we will combine the raw CSC voltage signals and do spike sort (spike detection and spike clustering). Spikes across the segments are combined in one file as the output of spike detection.

If not enough spikes are detected, spike clustering will be skipped, and no times_*.mat file will be created. In this case, linear interpolation in LFP extraction still removes the spikes from the raw signals.

Common Average Reference (CAR)

CAR involves calculating the average signal across all recording channels (electrodes) and subtracting this average from each channel's signal. This common average acts as a baseline reference, and subtracting it helps reduce signals common to all electrodes, such as external noise or systemic artifacts.

We apply CAR across channels within each bundle (8 channels) and subtract the median instead of the mean.

Power Line Interference:

For microelectrodes, we may use remove Power Line Interference (removePLI) to increase data quality for spike detection. In this case, signalRemovePLI will be added to the output of unpacked files. This takes extra space and running time. Avoid removePLI by setting arguments in read_CSC.m.

Applying power line interference cancellation tends to cause gaps without spikes in our data. Only use it when normal spike sorting fails due to power line noises.

Storage usage: twice the raw data (int16) as the added signal is in a single type.

Running time: For two hours of microdata (32k Hz), each channel takes 5 minutes to run removePLI on the Hoffman server. For .5 hours of microdata (32k Hz), it takes approximately 200 seconds to run removePLI on Apple Studio with an M1 Max chip.

Use check_removePLI_parameter.m to examine the effects of removePLI.

image

Keshtkaran, M. R., & Yang, Z. (2014). A fast, robust algorithm for power line interference cancellation in neural recording. Journal of neural engineering, 11(2), 026017.

  • M, number of harmonics to remove.
  • B, contains three elements [B0,Binf,Bst]:
    • B0, Initial notch bandwidth of the frequency estimator. Larger values of B0 are preferred (e.g. tens of Hz) to achieve a faster initial convergence.
    • Binf, Asymptotic notch bandwidth of the frequency estimator. Small values of Binf are preferred (e.g. tenths of Hz) to achieve more accurate estimation of the frequency.
    • Bst, Rate of convergence to 95% of the asymptotic bandwidth Binf
  • P, contains three elements [P0,Pinf,Pst]:
    • P0, Initial settling time of the frequency estimator. Smaller values of P0 are preferred (e.g. tenths of seconds) to achieve a faster initial convergence.
    • Pinf, Asymptotic settling time of the frequency estimator. Considering the fact that the power line frequency drifts are slow, larger values of Pinf are preferred (e.g. a few seconds) to obtain a more accurate estimation of the power line frequency.
    • Pst, Rate of convergence to 95% of the asymptotic settling time.
  • W, Settling time of the amplitude and phase estimator. It should be selected reasonably large to obtain an accurate estimation of the interference, thus avoiding the excessive removal of neural signals, while small enough to allow tracking of the interference amplitude fluctuations.
  • f_ac, Optional argument, the nominal AC frequency if known (50 Hz or 60 HZ)

cross-channel spike code:

If spikes are detected as the same time bin across all channels in a bundle, it will potentially be rejected as noise.

SGE job will assign channels in a single bundle to the same task. If there are not enough jobs assigned, channels in 2 or more bundles can be assigned to a task, which will affect how cross-channel spike codes are calculated. If the code is run on a local machine without SGE jobs, cross-channel spike code will be calculated based on all channels instead of channels in each bundle.

Single-Exp Analysis

Input Data

Experiment<exp_id>
  - CSC_micro
    - G[A-D][1-8]-<channel_name>_<suffix1>.mat
    - G[A-D][1-8]-<channel_name>_<suffix2>.mat
    ...

Output Data

Experiment<exp_id>
  - CSC_micro
    - G[A-D][1-8]-bundleMedian_<suffix1>.mat
    - G[A-D][1-8]-bundleMedian_<suffix2>.mat
    ...
  - CSC_micro_spikes
    - G[A-D][1-8]-<channel_name>_spikes.mat
      - spikes
      - spikesTimestamps
    - G[A-D][1-8]-<channel_name>_spikeCodes.mat
      - spikeCodes
      - spikeHist
      - spikeHistPrecise
    - times_G[A-D][1-8]-<channel_name>.mat
      - cluster_class
      - timestampsStart

In spike detection, spikes are extracted from the filtered signals (instead of the raw signals). spikeTimestamps is relative to the start of the experiment. The Unix timestamps can be recovered by adding timestampsStart.

Cluster_class: [class, spikeTimestamps]. Zero-based across all experiments

Multi-Exp Analysis

In multi-exp analysis, we use the CSC data from multiple experiments to calculate the spike detection threshold. The minimum threshold (no less than 18) across multiple experiments is used for spike detection. For details, check getDetectionThresh.

The minimum threshold does not apply to single-exp analysis The detected spikes are then combined before spike clustering. The number of experiments to combine varies across patients from 3 to 8. The procedure for multi-exp analysis is the same as for single-experiment analysis, with multiple segments and the same max gap duration.

Input Data

Experiment<exp_id1>
  - CSC_micro
    - G[A-D][1-8]-<channel_name>_<suffix1>.mat
    - G[A-D][1-8]-<channel_name>_<suffix2>.mat
    ...
Experiment<exp_id2>
  - CSC_micro
    - G[A-D][1-8]-<channel_name>_<suffix1>.mat
    - G[A-D][1-8]-<channel_name>_<suffix2>.mat
    ...
...

Output Data

Experiment<exp_id1>
  - CSC_micro
    - G[A-D][1-8]-bundleMedian_<suffix1>.mat
    - G[A-D][1-8]-bundleMedian_<suffix2>.mat
    ...
Experiment<exp_id2>
  - CSC_micro
    - G[A-D][1-8]-bundleMedian_<suffix1>.mat
    - G[A-D][1-8]-bundleMedian_<suffix2>.mat
    ...
Experiment-<exp_id1>-<exp_id2>-...
  - CSC_micro_spikes
    - G[A-D][1-8]-<channel_name>_spikes.mat
    - G[A-D][1-8]-<channel_name>_spikeCodes.mat
    - times_G[A-D][1-8]-<channel_name>.mat

Clone this wiki locally