braintools.analysis.spike.spike_distance

Contents

braintools.analysis.spike.spike_distance#

braintools.analysis.spike.spike_distance(a, b, t_start=None, t_stop=None)#

Compute the SPIKE-distance between two populations, neuron by neuron.

The SPIKE-distance of Kreuz et al. [1] compares spike timing. For a train \(x\) at time \(t\), let \(t_P^x \le t < t_F^x\) be the previous and following spikes, \(x_{ISI} = t_F^x - t_P^x\), and \(\Delta t_P^x\), \(\Delta t_F^x\) the distances of those spikes to the nearest spike of the other train. The local spike-time difference is interpolated between them,

\[S_x(t) = \frac{\Delta t_P^x (t_F^x - t) + \Delta t_F^x (t - t_P^x)}{x_{ISI}},\]

and the profile weights both trains by their local ISIs,

\[S(t) = \frac{S_a(t)\, b_{ISI} + S_b(t)\, a_{ISI}} {2 \langle x_{ISI} \rangle^2}, \qquad \langle x_{ISI} \rangle = \frac{a_{ISI} + b_{ISI}}{2}.\]

The distance is the time average of \(S(t)\) over the recording. It is 0 for identical trains and at most 1.

Parameters:

  • a (SpikeTrains) – The two populations, with the same number of neurons N (their capacities may differ). Neuron i of a is compared with neuron i of b. Both must carry compatible time units, or both none.

  • b (SpikeTrains) – The two populations, with the same number of neurons N (their capacities may differ). Neuron i of a is compared with neuron i of b. Both must carry compatible time units, or both none.

  • t_start (Any | None) – Averaging interval. A Quantity is converted to the unit of the spike times; a plain number is in milliseconds (the brainstate convention), so it is used as is for unitless trains and converted from milliseconds for trains that carry units. Defaults to min(a.t_start, b.t_start) and max(a.t_stop, b.t_stop), as PySpike does. Spikes outside the interval are ignored.

  • t_stop (Any | None) – Averaging interval. A Quantity is converted to the unit of the spike times; a plain number is in milliseconds (the brainstate convention), so it is used as is for unitless trains and converted from milliseconds for trains that carry units. Defaults to min(a.t_start, b.t_start) and max(a.t_stop, b.t_stop), as PySpike does. Spikes outside the interval are ignored.

Returns:

SPIKE-distances, shape (N,), dimensionless. nan when t_stop <= t_start.

Return type:

Array

Raises:

  • TypeError – If a or b is not a SpikeTrains, if only one carries units, or if t_start/t_stop is a Quantity that is not a time.

  • ValueError – If the numbers of neurons differ, or a SpikeTrains is batched.

See also

isi_distance

The ISI-distance, sensitive to local firing rates.

van_rossum_distance

Kernel-based distance between two spike trains.

braintools.spike.SpikeTrains

The padded spike-time container.

Notes

\(S(t)\) is linear on every interval between the merged spike times of the two trains, so the time average is integrated exactly with the trapezoid rule over those intervals.

The edge conventions are those of PySpike 0.9 [3], which follow Kreuz et al. [2], so the values agree with pyspike.spike_distance to round-off:

  • auxiliary spikes are placed at min(t_start, s_0 - (s_1 - s_0)) and max(t_stop, s_{n-1} + (s_{n-1} - s_{n-2})) (at t_start and t_stop for a single spike); they serve as the previous/following spikes at the edges and as candidates for the nearest spike of the other train;

  • before the first and after the last spike, \(S_x\) is held at the distance of that first or last spike, and the edge ISIs are those of isi_distance();

  • a train whose only spike is at t_start uses t_stop as its following spike;

  • an empty train is replaced by spikes at t_start and t_stop; repeated spike times within a train count once.

One deviation: PySpike returns nan when both trains have a single spike exactly at t_stop; here the distance is 0.

The function works under jax.jit and jax.vmap (over trials, or over neurons to build a pair matrix, see the examples). Gradients with respect to the spike times are not supported (the nearest-neighbour search is piecewise). Memory and time are O(N (K_a + K_b) log(K_a + K_b)).

References

Examples

>>> import jax
>>> import jax.numpy as jnp
>>> import numpy as np
>>> from braintools.spike import SpikeTrains
>>> from braintools.analysis.spike import spike_distance
>>> a = SpikeTrains.from_list([[1., 5., 9.], [3.]], t_start=0., t_stop=10.)
>>> b = SpikeTrains.from_list([[1., 5., 9.], [6.]], t_start=0., t_stop=10.)
>>> np.round(np.asarray(spike_distance(a, b)), 4)
array([0.    , 0.5566], dtype=float32)
>>>
>>> # From a dense raster (time on axis 0), with units.
>>> import brainunit as u
>>> raster = jnp.zeros((100, 2)).at[::10, 0].set(1.).at[5::10, 1].set(1.)
>>> st = SpikeTrains.from_dense(raster, dt=1. * u.ms)
>>> d = spike_distance(st.select(0), st.select(1))
>>> d.shape
(1,)
>>>
>>> # All pairs of neurons of one population: vmap over the first neuron.
>>> n = st.n_neurons
>>> pair = jax.vmap(lambda i: spike_distance(st.select(jnp.full((n,), i)), st))
>>> pair(jnp.arange(n)).shape
(2, 2)
Contents