ASSET is a statistical method [1] for the detection of repeating sequences of synchronous spiking events in parallel spike trains. Given a list sts of spike trains, the analysis comprises the following steps:
Build the intersection matrix imat (optional) and the associated probability matrix pmat with the desired bin size:
>>> binsize = 5 * pq.ms
>>> dt = 1 * pq.s
>>> imat, xedges, yedges = intersection_matrix(sts, binsize, dt, norm=2)
>>> pmat, xedges, yedges = probability_matrix_analytical(
sts, binsize, dt)
Compute the joint probability matrix jmat, using a suitable filter:
>>> filter_shape = (5,2) # filter shape
>>> nr_neigh = 5 # nr of largest neighbors
>>> jmat = joint_probability_matrix(pmat, filter_shape, nr_neigh)
Create from pmat and jmat a masked version of the intersection matrix:
>>> alpha1 = 0.99
>>> alpha2 = 0.99999
>>> mask = mask_matrices([pmat, jmat], [alpha1, alpha2])
Cluster significant elements of imat into diagonal structures (“DSs”):
>>> epsilon = 10
>>> minsize = 2
>>> stretch = 5
>>> cmat = asset.cluster_matrix_entries(mask, epsilon, minsize, stretch)
Extract sequences of synchronous events associated to each worm
>>> extract_sse(sts, x_edges, y_edges, cmat)
References:
[1] Torre, Canova, Denker, Gerstein, Helias, Gruen (submitted)
elephant.asset.
cluster_matrix_entries
(mat, eps=10, min=2, stretch=5)[source]¶Given a matrix mat, replaces its positive elements with integers representing different cluster ids. Each cluster comprises close-by elements.
In ASSET analysis, mat is a thresholded (“masked”) version of an intersection matrix imat, whose values are those of imat only if considered statistically significant, and zero otherwise.
A cluster is built by pooling elements according to their distance, via the DBSCAN algorithm (see sklearn.cluster.dbscan()). Elements form a neighbourhood if at least one of them has a distance not larger than eps from the others, and if they are at least min. Overlapping neighborhoods form a cluster.
- Clusters are assigned integers from 1 to the total number k of clusters
- Unclustered (“isolated”) positive elements of mat are assigned value -1
- Non-positive elements are assigned the value 0.
The distance between the positions of two positive elements in mat is given by an Euclidean metric which is stretched if the two positions are not aligned along the 45 degree direction (the main diagonal direction), as more, with maximal stretching along the anti-diagonal. Specifically, the Euclidean distance between positions (i1, j1) and (i2, j2) is stretched by a factor
where is the angle between the pixels and the 45deg
direction. The stretching factor thus varies between 1 and stretch.
Parameters: | mat : numpy.ndarray
eps : float >=0, optional
min : int, optional
stretch : float > 1, optional
|
---|---|
Returns: | cmat : numpy.ndarray of integers
|
elephant.asset.
extract_sse
(spiketrains, x_edges, y_edges, cmat, ids=None)[source]¶Given a list of spike trains, two arrays of bin edges and a clustered intersection matrix obtained from those spike trains via worms analysis using the specified edges, extracts the sequences of synchronous events (SSEs) corresponding to clustered elements in the cluster matrix.
Parameters: | spiketrains : list of neo.SpikeTrain
x_edges : quantities.Quantity
y_edges : quantities.Quantity
cmat: numpy.ndarray
ids : list or None, optional
|
---|---|
Returns: | sse : dict
|
elephant.asset.
intersection_matrix
(spiketrains, binsize, dt, t_start_x=None, t_start_y=None, norm=None)[source]¶Generates the intersection matrix from a list of spike trains.
Given a list of SpikeTrains, consider two binned versions of them differing for the starting time of the binning (t_start_x and t_start_y respectively; the two times can be identical). Then calculate the intersection matrix M of the two binned data, where M[i,j] is the overlap of bin i in the first binned data and bin j in the second binned data (i.e. the number of spike trains spiking both at bin i and at bin j). The matrix entries can be normalized to values between 0 and 1 via different normalizations (see below).
Parameters: | spiketrains : list of neo.SpikeTrains
binsize : quantities.Quantity
dt : quantities.Quantity
t_start_x : quantities.Quantity, optional
t_start_y : quantities.Quantity, optional
norm : int, optional
|
---|---|
Returns: | imat : numpy.ndarray of floats
x_edges : numpy.ndarray
y_edges : numpy.ndarray
|
elephant.asset.
joint_probability_matrix
(pmat, filter_shape, nr_largest=None, alpha=0, pvmin=1e-05)[source]¶Map a probability matrix pmat to a joint probability matrix jmat, where jmat[i, j] is the joint p-value of the largest neighbors of pmat[i, j].
The values of pmat are assumed to be uniformly distributed in the range [alpha, 1] (alpha=0 by default). Centered a rectangular kernel of shape filter_shape=(l, w) around each entry pmat[i, j], aligned along the diagonal where pmat[i, j] lies into, extracts the nr_largest highest values falling within the kernel and computes their joint p-value jmat[i, j] (see [1]).
Parameters: | pmat : ndarray
filter_shape : tuple
nr_largest : int, optional
alpha : float in [0, 1), optional
pvmin : flaot in [0, 1), optional
|
---|---|
Returns: | jmat : numpy.ndarray
|
References
[1] Torre et al (in prep) ...
elephant.asset.
mask_matrices
(matrices, thresholds)[source]¶Given a list of matrices and a list of thresholds, return a boolean matrix B (“mask”) such that B[i,j] is True if each input matrix in the list strictly exceeds the corresponding threshold at that position.
Parameters: | matrices : list of numpy.ndarrays
thresholds : list of floats
|
---|---|
Returns: | mask : numpy.ndarray of bools
|
elephant.asset.
probability_matrix_analytical
(spiketrains, binsize, dt, t_start_x=None, t_start_y=None, fir_rates='estimate', kernel_width=array(100.) * ms, verbose=False)[source]¶Given a list of spike trains, approximates the cumulative probability of each entry in their intersection matrix (see: intersection_matrix()).
The approximation is analytical and works under the assumptions that the input spike trains are independent and Poisson. It works as follows:
- Bin each spike train at the specified binsize: this yields a binary array of 1s (spike in bin) and 0s (no spike in bin) (clipping used)
- If required, estimate the rate profile of each spike train by convolving the binned array with a boxcar kernel of user-defined length
- For each neuron k and each pair of bins i and j, compute the probability p_ijk that neuron k fired in both bins i and j.
- Approximate the probability distribution of the intersection value at (i, j) by a Poisson distribution with mean parameter l = sum_k (p_ijk), justified by Le Cam’s approximation of a sum of independent Bernouilli random variables with a Poisson distribution.
Parameters: | spiketrains : list of neo.SpikeTrains
binsize : quantities.Quantity
dt : quantities.Quantity
t_start_x, t_start_y : quantities.Quantity, optional
fir_rates: list of neo.AnalogSignals or ‘estimate’, optional
kernel_width : quantities.Quantity, optional
verbose : bool, optional
|
---|---|
Returns: | pmat : numpy.ndarray
x_edges : numpy.ndarray
y_edges : numpy.ndarray
|
elephant.asset.
probability_matrix_montecarlo
(spiketrains, binsize, dt, t_start_x=None, t_start_y=None, surr_method='dither_spike_train', j=None, n_surr=100, verbose=False)[source]¶by a Monte Carlo approach using surrogate data. Contrarily to the analytical version (see: probability_matrix_analytical()) the Monte Carlo one does not incorporate the assumptions of Poissonianity in the null hypothesis.
The method produces surrogate spike trains (using one of several methods at disposal, see below) and calculates their intersection matrix M. For each entry (i, j), the intersection cdf P[i, j] is then given by:
P[i, j] = #(spike_train_surrogates such that M[i, j] < I[i, j]) / #(spike_train_surrogates)
If P[i, j] is large (close to 1), I[i, j] is statistically significant: the probability to observe an overlap equal to or larger then I[i, j] under the null hypothesis is 1-P[i, j], very small.
Parameters: | sts : list of neo.SpikeTrains
binsize : quantities.Quantity
dt : quantities.Quantity
t_start_x, t_start_y : quantities.Quantity, optional
surr_method : str, optional
j : quantities.Quantity, optional
n_surr : int, optional
|
---|---|
Returns: | pmat : ndarray
|
See also
elephant.asset.
sse_difference
(sse1, sse2, difference='linkwise')[source]¶Given two sequences of synchronous events (SSEs) sse1 and sse2, each consisting of a pool of pixel positions and associated synchronous events (see below), computes the difference between sse1 and sse2. The difference can be performed ‘pixelwise’ or ‘linkwise’:
- if ‘pixelwise’, it yields a new SSE which contains all (and only) the events in sse1 whose pixel position doesn’t match any pixel in sse2.
- if ‘linkwise’, for each pixel (i, j) in sse1 and corresponding synchronous event S1, if (i, j) is a pixel in sse2 corresponding to the event S2, it retains the set difference S1 - S2. If (i, j) is not a pixel in sse2, it retains the full set S1.
Note that in either case the difference is a non-symmetric operation: intersection(sse1, sse2) != intersection(sse2, sse1).
Both sse1 and sse2 must be provided as dictionaries of the type
{(i1, j1): S1, (i2, j2): S2, ..., (iK, jK): SK},
where each i, j is an integer and each S is a set of neuron ids. (See also: extract_sse() that extracts SSEs from given spiketrains).
Parameters: | sse1, sse2 : each a dict
difference : str, optional
|
---|---|
Returns: | sse : dict
|
elephant.asset.
sse_intersection
(sse1, sse2, intersection='linkwise')[source]¶Given two sequences of synchronous events (SSEs) sse1 and sse2, each consisting of a pool of positions (iK, jK) of matrix entries and associated synchronous events SK, finds the intersection among them. The intersection can be performed ‘pixelwise’ or ‘linkwise’.
- if ‘pixelwise’, it yields a new SSE which retains only events in sse1 whose pixel position matches a pixel position in sse2. This operation is not symmetric: intersection(sse1, sse2) != intersection(sse2, sse1).
- if ‘linkwise’, an additional step is performed where each retained synchronous event SK in sse1 is intersected with the corresponding event in sse2. This yields a symmetric operation: intersection(sse1, sse2) = intersection(sse2, sse1).
Both sse1 and sse2 must be provided as dictionaries of the type
{(i1, j1): S1, (i2, j2): S2, ..., (iK, jK): SK},
where each i, j is an integer and each S is a set of neuron ids. (See also: extract_sse() that extracts SSEs from given spiketrains).
Parameters: | sse1, sse2 : each a dict
intersection : str, optional
|
---|---|
Returns: | sse : dict
|
elephant.asset.
sse_isdisjoint
(sse1, sse2)[source]¶Given two sequences of synchronous events (SSEs) sse1 and sse2, each consisting of a pool of pixel positions and associated synchronous events (see below), determines whether sse1 and sse2 are disjoint. Two SSEs are disjoint if they don’t share pixels, or if the events associated to common pixels are disjoint.
Both sse1 and sse2 must be provided as dictionaries of the type
{(i1, j1): S1, (i2, j2): S2, ..., (iK, jK): SK},
where each i, j is an integer and each S is a set of neuron ids. (See also: extract_sse() that extracts SSEs from given spiketrains).
Parameters: | sse1, sse2 : each a dictionary
|
---|---|
Returns: | is_disjoint : bool
|
elephant.asset.
sse_isequal
(sse1, sse2)[source]¶Given two sequences of synchronous events (SSEs) sse1 and sse2, each consisting of a pool of pixel positions and associated synchronous events (see below), determines whether sse1 is strictly contained in sse2. sse1 is strictly contained in sse2 if all its pixels are pixels of sse2, if its associated events are subsets of the corresponding events in sse2, and if sse2 contains events, or neuron ids in some event, which do not belong to sse1 (i.e. sse1 and sse2 are not identical)
Both sse1 and sse2 must be provided as dictionaries of the type
{(i1, j1): S1, (i2, j2): S2, ..., (iK, jK): SK},
where each i, j is an integer and each S is a set of neuron ids. (See also: extract_sse() that extracts SSEs from given spiketrains).
Parameters: | sse1, sse2 : each a dict
|
---|---|
Returns: | is_equal : bool
|
elephant.asset.
sse_issub
(sse1, sse2)[source]¶Given two sequences of synchronous events (SSEs) sse1 and sse2, each consisting of a pool of pixel positions and associated synchronous events (see below), determines whether sse1 is strictly contained in sse2. sse1 is strictly contained in sse2 if all its pixels are pixels of sse2, if its associated events are subsets of the corresponding events in sse2, and if sse2 contains non-empty events, or neuron ids in some event, which do not belong to sse1 (i.e. sse1 and sse2 are not identical)
where each i, j is an integer and each S is a set of neuron ids. (See also: extract_sse() that extracts SSEs from given spiketrains).
Parameters: | sse1, sse2 : each a dict
|
---|---|
Returns: | is_sub : bool
|
elephant.asset.
sse_issuper
(sse1, sse2)[source]¶Given two sequences of synchronous events (SSEs) sse1 and sse2, each consisting of a pool of pixel positions and associated synchronous events (see below), determines whether sse1 strictly contains sse2. sse1 strictly contains sse2 if it contains all pixels of sse2, if all associated events in sse1 contain those in sse2, and if sse1 additionally contains other pixels / events not contained in sse2.
where each i, j is an integer and each S is a set of neuron ids. (See also: extract_sse() that extracts SSEs from given spiketrains).
Note: sse_issuper(sse1, sse2) is identical to sse_issub(sse2, sse1).
Parameters: | sse1, sse2 : each a dict
|
---|---|
Returns: | is_super : bool
|
elephant.asset.
sse_overlap
(sse1, sse2)[source]¶Given two sequences of synchronous events (SSEs) sse1 and sse2, each consisting of a pool of pixel positions and associated synchronous events (see below), determines whether sse1 strictly contains sse2. sse1 strictly contains sse2 if it contains all pixels of sse2, if all associated events in sse1 contain those in sse2, and if sse1 additionally contains other pixels / events not contained in sse2.
where each i, j is an integer and each S is a set of neuron ids. (See also: extract_sse() that extracts SSEs from given spiketrains).
Note: sse_issuper(sse1, sse2) is identical to sse_issub(sse2, sse1).
Parameters: | sse1, sse2 : each a dict
|
---|---|
Returns: | is_super : bool
|