BinnedSpikeTrain (conversion)

This module allows to convert standard data representations (e.g., a spike train stored as Neo SpikeTrain object) into other representations useful to perform calculations on the data. An example is the representation of a spike train as a sequence of 0-1 values (binned spike train).

class elephant.conversion.BinnedSpikeTrain(spiketrains, binsize=None, num_bins=None, t_start=None, t_stop=None)[source]

Class which calculates a binned spike train and provides methods to transform the binned spike train to a boolean matrix or a matrix with counted time points.

A binned spike train represents the occurrence of spikes in a certain time frame. I.e., a time series like [0.5, 0.7, 1.2, 3.1, 4.3, 5.5, 6.7] is represented as [0, 0, 1, 3, 4, 5, 6]. The outcome is dependent on given parameter such as size of bins, number of bins, start and stop points.

A boolean matrix represents the binned spike train in a binary (True/False) manner. Its rows represent the number of spike trains and the columns represent the binned index position of a spike in a spike train. The calculated matrix entry containing True indicates a spike.

A matrix with counted time points is calculated the same way, but its entries contain the number of spikes that occurred in the given bin of the given spike train.

Parameters:
spiketrainsneo.SpikeTrain or list of neo.SpikeTrain or np.ndarray

Spike train(s) to be binned.

binsizepq.Quantity, optional

Width of a time bin. Default: None

num_binsint, optional

Number of bins of the binned spike train. Default: None

t_startpq.Quantity, optional

Time of the left edge of the first bin (left extreme; included). Default: None

t_stoppq.Quantity, optional

Time of the right edge of the last bin (right extreme; excluded). Default: None

Raises:
AttributeError

If less than 3 optional parameters are None.

TypeError

If spiketrains is an np.ndarray with dimensionality different than NxM or if type of num_bins is not an int or num_bins < 0.

ValueError

When number of bins calculated from t_start, t_stop and binsize differs from provided num_bins or if t_stop of any spike train is smaller than any t_start or if any spike train does not cover the full [t_start, t_stop`] range.

Warns:
UserWarning

If some spikes fall outside of [t_start, t_stop] range

See also

_convert_to_binned
spike_indices
to_bool_array
to_array

Notes

There are four minimal configurations of the optional parameters which have to be provided, otherwise a ValueError will be raised: * t_start, num_bins, binsize * t_start, num_bins, t_stop * t_start, bin_size, t_stop * t_stop, num_bins, binsize

If spiketrains is a neo.SpikeTrain or a list thereof, it is enough to explicitly provide only one parameter: num_bins or binsize. The t_start and t_stop will be calculated from given spiketrains (max t_start and min t_stop of `neo.SpikeTrain`s). Missing parameter will be calculated automatically. All parameters will be checked for consistency. A corresponding error will be raised, if one of the four parameters does not match the consistency requirements.

property bin_centers

Returns each center time point of all bins between t_start and t_stop points.

The center of each bin of all time steps between start and stop.

Returns:
bin_edgespq.Quantity

All center edges in interval (start, stop).

property bin_edges

Returns all time edges as a quantity array with num_bins bins.

The borders of all time steps between t_start and t_stop with a step binsize. It is crucial for many analyses that all bins have the same size, so if t_stop - t_start is not divisible by binsize, there will be some leftover time at the end (see https://github.com/NeuralEnsemble/elephant/issues/255). The length of the returned array should match num_bins.

Returns:
bin_edgespq.Quantity

All edges in interval [t_start, t_stop] with num_bins bins are returned as a quantity array.

binarize(self, copy=True)[source]

Clip the internal array (no. of spikes in a bin) to 0 (no spikes) or 1 (at least one spike) values only.

Parameters:
copybool

Perform the clipping in-place (False) or on a copy (True). Default: True.

Returns:
bstBinnedSpikeTrain

BinnedSpikeTrain with both sparse and dense (if present) array representation clipped to 0 (no spike) or 1 (at least one spike) entries.

get_num_of_spikes(self, axis=None)[source]

Compute the number of binned spikes.

Parameters:
axisint, optional

If None, compute the total num. of spikes. Otherwise, compute num. of spikes along axis. If axis is 1, compute num. of spikes per spike train (row). Default is None.

Returns:
n_spikes_per_rowint or np.ndarray

The number of binned spikes.

property is_binary

Checks and returns True if given input is a binary input. Beware, that the function does not know if the input is binary because e.g to_bool_array() was used before or if the input is just sparse (i.e. only one spike per bin at maximum).

Returns:
bool

True for binary input, False otherwise.

remove_stored_array(self)[source]

Unlinks the matrix with counted time points from memory.

property sparsity
Returns:
float

Matrix sparsity defined as no. of nonzero elements divided by the matrix size

property spike_indices

A list of lists for each spike train (i.e., rows of the binned matrix), that in turn contains for each spike the index into the binned matrix where this spike enters.

In contrast to to_sparse_array().nonzero(), this function will report two spikes falling in the same bin as two entries.

Examples

>>> import elephant.conversion as conv
>>> import neo as n
>>> import quantities as pq
>>> st = n.SpikeTrain([0.5, 0.7, 1.2, 3.1, 4.3, 5.5, 6.7] * pq.s,
...                   t_stop=10.0 * pq.s)
>>> x = conv.BinnedSpikeTrain(st, num_bins=10, binsize=1 * pq.s,
...                           t_start=0 * pq.s)
>>> print(x.spike_indices)
[[0, 0, 1, 3, 4, 5, 6]]
>>> print(x.to_sparse_array().nonzero()[1])
[0 1 3 4 5 6]
>>> print(x.to_array())
[[2, 1, 0, 1, 1, 1, 1, 0, 0, 0]]
to_array(self, store_array=False)[source]

Returns a dense matrix, calculated from the sparse matrix, with counted time points of spikes. The rows correspond to spike trains and the columns correspond to bins in a BinnedSpikeTrain. Entries contain the count of spikes that occurred in the given bin of the given spike train. If the boolean store_array is set to True, the matrix will be stored in memory.

Returns:
matrixnp.ndarray

Matrix with spike counts. Columns represent the index positions of the binned spikes and rows represent the spike trains.

See also

scipy.sparse.csr_matrix
scipy.sparse.csr_matrix.toarray

Examples

>>> import elephant.conversion as conv
>>> import neo as n
>>> a = n.SpikeTrain([0.5, 0.7, 1.2, 3.1, 4.3, 5.5, 6.7] * pq.s,
...                  t_stop=10.0 * pq.s)
>>> x = conv.BinnedSpikeTrain(a, num_bins=10, binsize=1 * pq.s,
...                           t_start=0 * pq.s)
>>> print(x.to_array())
[[2 1 0 1 1 1 1 0 0 0]]
to_bool_array(self)[source]

Returns a matrix, in which the rows correspond to the spike trains and the columns correspond to the bins in the BinnedSpikeTrain. True indicates a spike in given bin of given spike train and False indicates lack of spikes.

Returns:
numpy.ndarray

Returns a dense matrix representation of the sparse matrix, with True indicating a spike and False indicating a no-spike. The columns represent the index position of the bins and rows represent the number of spike trains.

See also

scipy.sparse.csr_matrix
scipy.sparse.csr_matrix.toarray

Examples

>>> import elephant.conversion as conv
>>> import neo as n
>>> import quantities as pq
>>> a = n.SpikeTrain([0.5, 0.7, 1.2, 3.1, 4.3, 5.5, 6.7] * pq.s,
...                  t_stop=10.0 * pq.s)
>>> x = conv.BinnedSpikeTrain(a, num_bins=10, binsize=1 * pq.s,
...                           t_start=0 * pq.s)
>>> print(x.to_bool_array())
[[ True  True False  True  True  True  True False False False]]
to_sparse_array(self)[source]

Getter for sparse matrix with time points.

Returns:
scipy.sparse.csr_matrix

Sparse matrix, version with spike counts.

See also

scipy.sparse.csr_matrix
to_array
to_sparse_bool_array(self)[source]

Getter for boolean version of the sparse matrix, calculated from sparse matrix with counted time points.

Returns:
scipy.sparse.csr_matrix

Sparse matrix, binary, boolean version.

See also

scipy.sparse.csr_matrix
to_bool_array
elephant.conversion.binarize(spiketrain, sampling_rate=None, t_start=None, t_stop=None, return_times=False)[source]

Return an array indicating if spikes occurred at individual time points.

The array contains boolean values identifying whether at least one spike occurred in the corresponding time bin. Time bins start at t_start and end at t_stop, spaced in 1/sampling_rate intervals.

Accepts either a neo.SpikeTrain, a pq.Quantity array, or a plain np.ndarray. Returns a boolean array with each element indicating the presence or absence of a spike in that time bin.

Optionally also returns an array of time points corresponding to the elements of the boolean array. The units of this array will be the same as the units of the neo.SpikeTrain, if any.

Parameters:
spiketrainneo.SpikeTrain or pq.Quantity or np.ndarray

The spike times. Does not have to be sorted.

sampling_ratefloat or pq.Quantity, optional

The sampling rate to use for the time points. If not specified, retrieved from the sampling_rate attribute of spiketrain. Default: None.

t_startfloat or pq.Quantity, optional

The start time to use for the time points. If not specified, retrieved from the t_start attribute of spiketrain. If this is not present, defaults to 0. Any element of spiketrain lower than t_start is ignored. Default: None.

t_stopfloat or pq.Quantity, optional

The stop time to use for the time points. If not specified, retrieved from the t_stop attribute of spiketrain. If this is not present, defaults to the maximum value of spiketrain. Any element of spiketrain higher than t_stop is ignored. Default: None.

return_timesbool, optional

If True, also return the corresponding time points. Default: False.

Returns:
valuesnp.ndarray of bool

A True value at a particular index indicates the presence of one or more spikes at the corresponding time point.

timesnp.ndarray or pq.Quantity, optional

The time points. This will have the same units as spiketrain. If spiketrain has no units, this will be an np.ndarray array.

Raises:
TypeError

If spiketrain is an np.ndarray and t_start, t_stop, or sampling_rate is a pq.Quantity.

ValueError

If sampling_rate is not explicitly defined and not present as an attribute of spiketrain.

Notes

Spike times are placed in the bin of the closest time point, going to the higher bin if exactly between two bins.

So in the case where the bins are 5.5 and 6.5, with the spike time being 6.0, the spike will be placed in the 6.5 bin.

The upper edge of the last bin, equal to t_stop, is inclusive. That is, a spike time exactly equal to t_stop will be included.

If spiketrain is a pq.Quantity or neo.SpikeTrain and t_start, t_stop or sampling_rate is not, then the arguments that are not pq.Quantity will be assumed to have the same units as spiketrain.