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 neuronsN(their capacities may differ). Neuroniofais compared with neuroniofb. Both must carry compatible time units, or both none.b (
SpikeTrains) – The two populations, with the same number of neuronsN(their capacities may differ). Neuroniofais compared with neuroniofb. 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 (thebrainstateconvention), so it is used as is for unitless trains and converted from milliseconds for trains that carry units. Defaults tomin(a.t_start, b.t_start)andmax(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 (thebrainstateconvention), so it is used as is for unitless trains and converted from milliseconds for trains that carry units. Defaults tomin(a.t_start, b.t_start)andmax(a.t_stop, b.t_stop), as PySpike does. Spikes outside the interval are ignored.
- Returns:
ISI-distances, shape
(N,), dimensionless.nanwhent_stop <= t_start.- Return type:
Array- Raises:
TypeError – If
aorbis not aSpikeTrains, if only one carries units, or ift_start/t_stopis a Quantity that is not a time.ValueError – If the numbers of neurons differ, or a
SpikeTrainsis batched.
See also
spike_distanceThe SPIKE-distance, sensitive to spike timing.
victor_purpura_distanceEdit distance between two spike trains.
braintools.spike.SpikeTrainsThe 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_distanceto round-off:before the first spike, the ISI is
max(s_0 - t_start, s_1 - s_0), and after the last spikemax(t_stop - s_{n-1}, s_{n-1} - s_{n-2}); with a single spike these ares_0 - t_startandt_stop - s_0;an empty train is replaced by spikes at
t_startandt_stop, so two empty trains are at distance 0 and an empty train against a train with one spike att_start + T/2is at distance 0.5;repeated spike times within a train count once.
One deviation: PySpike returns
nanwhen both trains have a single spike exactly att_stop; here the zero-length interval after it is ignored and the distance is 0.The function works under
jax.jitandjax.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 areO(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)