braintools.analysis.spike.isi_distance

Contents

braintools.analysis.spike.isi_distance#

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

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

The ISI-distance of Kreuz et al. [1] compares the local firing rates of two spike trains. At every instant \(t\), let \(\nu_a(t)\) and \(\nu_b(t)\) be the inter-spike intervals that contain \(t\) in each train. The ISI-profile is

\[I(t) = \frac{|\nu_a(t) - \nu_b(t)|}{\max(\nu_a(t), \nu_b(t))},\]

and the distance is its time average over the recording,

\[D_I = \frac{1}{t_{stop} - t_{start}} \int_{t_{start}}^{t_{stop}} I(t)\, dt .\]

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:

ISI-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

spike_distance

The SPIKE-distance, sensitive to spike timing.

victor_purpura_distance

Edit distance between two spike trains.

braintools.spike.SpikeTrains

The padded spike-time container.

Notes

The profile is piecewise constant between the merged spike times of the two trains, so the integral is computed exactly as a sum over those intervals. The edge conventions are those of PySpike 0.9 [3] (Kreuz et al. [2]), so the values agree with pyspike.isi_distance to round-off:

  • before the first spike, the ISI is max(s_0 - t_start, s_1 - s_0), and after the last spike max(t_stop - s_{n-1}, s_{n-1} - s_{n-2}); with a single spike these are s_0 - t_start and t_stop - s_0;

  • an empty train is replaced by spikes at t_start and t_stop, so two empty trains are at distance 0 and an empty train against a train with one spike at t_start + T/2 is at distance 0.5;

  • 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 zero-length interval after it is ignored and 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). It is not differentiable in a useful sense: the profile is piecewise constant in the spike times. 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 isi_distance
>>> a = SpikeTrains.from_list([[1., 5., 9.], [2., 4., 6., 8.]], t_start=0., t_stop=10.)
>>> b = SpikeTrains.from_list([[1., 5., 9.], []], t_start=0., t_stop=10.)
>>> isi_distance(a, b)
Array([0. , 0.8], 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[::20, 1].set(1.)
>>> st = SpikeTrains.from_dense(raster, dt=1. * u.ms)
>>> d = isi_distance(st.select(0), st.select(1))
>>> np.round(np.asarray(d), 4)
array([0.5], dtype=float32)
>>>
>>> # All pairs of neurons of one population: vmap over the first neuron.
>>> n = st.n_neurons
>>> pair = jax.vmap(lambda i: isi_distance(st.select(jnp.full((n,), i)), st))
>>> np.round(np.asarray(pair(jnp.arange(n))), 4)
array([[0. , 0.5],
       [0.5, 0. ]], dtype=float32)