Source Code - SPIKE- and ISI-Distance


 

SPIKE- and ISI-Distance between two or more spike trains


For a detailed description of the methods please refer to:

 

Kreuz T, Chicharro D, Houghton C, Andrzejak RG, Mormann F:
Monitoring spike train synchrony.
J Neurophysiol 109, 1457 (2013) [PDF].


and

Kreuz T:
SPIKE-distance.
Scholarpedia 7(12), 30652 (2012).

For an overview please refer to:

Mulansky M, Bozanic N, Sburlea A, Kreuz T:
A guide to time-resolved and parameter-free measures of spike train synchrony.
IEEE Proceeding on Event-based Control, Communication, and Signal Processing (EBCCSP), 1-8 and arXiv [PDF] (2015).



Powerpoint presentation

Kreuz T:
Monitoring spike train synchrony.
[PDF].


Further information on the SPIKE-distance can also be found on the  SPIKY Facebook page. New!


The best demonstration of the SPIKE-distance can be seen in the following movie. This and many more movies can now also be seen on the SPIKY Youtube Channel. New!


Movie: How the SPIKE-distance traces spike train clustering

Here we explain how this time-resolved and parameter-free measure traces spike train clustering with both its normal and its realtime version.

Top: Rasterplot with artificially generated spike trains whose clustering changes every 500 ms.
Middle: Dissimilarity matrices for the regular (1st column) and the real-time (2nd column) SPIKE-distance. The spike train group matrices (3rd and 4th column) are obtained by averaging over the respective submatrices of the original dissimilarity matrices.
Bottom: Corresponding dendrograms.

In the first part of the movie the instantaneous dissimilarities and the clustering dendrograms of both measures are updated as the green time bar moves from left to right. The second part comprises several examples of selective temporal averaging of individual or combined intervals as well as triggered averaging.


 

A second example:

Movie: Demonstration of externally triggered averaging.

Here we explain how the instantaneous nature of the SPIKE-distance can be very useful when this method is combined with triggered averaging.

In this simulated setup a population of 50 neurons is recorded upon presentation of a time-varying non-periodic stimulus, in this case a chirp function. It is assumed that half of the neurons are sensitive to negative amplitudes and accordingly exhibit higher (lower) reliability for local minima (maxima) of the chirp function.

As the amplitude of the chirp function varies, so does the spike train synchrony of half of the neurons, but due to the non-periodicity of the stimulus this is quite difficult to detect. Detection is facilitated by externally triggered averaging where the triggering is performed on common stimulus features, here the amplitude of the chirp function. As can be seen in the dissimilarity matrices of the SPIKE-distance and even better in the dendrogram, a decrease of this stimulus amplitude leads to the emergence of a spike train cluster consisting of the modulated spike trains which indicates their increase in reliability.




The definiton of the SPIKE-distance builds on these earlier papers:

Kreuz T, Chicharro D, Greschner M, Andrzejak RG:
Time-resolved and time-scale adaptive measures of spike train synchrony.
J Neurosci Methods 195, 92 (2011) [PDF].

Kreuz T, Chicharro D, Andrzejak RG, Haas JS, Abarbanel HDI:
Measuring multiple spike train synchrony.
J Neurosci Methods 183, 287 (2009) [PDF].

Kreuz T, Haas JS, Morelli A, Abarbanel HDI, Politi A:
Measuring spike train synchrony.
J Neurosci Methods 165, 151 (2007) [PDF].


The Matlab code is now complemented by PySpike, an open source Python library written by Mario Mulansky and available on github. Its core functionality is the implementation of the bivariate ISI- and SPIKE-distance. Additionally, it provides functions to compute multivariate ISI- and SPIKE-distances, as well as averaging and general spike train processing. All computation intensive parts are implemented in C (via cython) to reach a competitive performance (factor 100-200 over plain Python). There is also documentation, currently consisting of a brief introduction and the API reference.


See also:

Python-Implementation of the pairwise SPIKE-distance (written by Jeremy Fix)

Python-Implementation of the pairwise ISI-distance (maintained by Michael Chary)

 

===================================


Matlab code to calculate and visualize the SPIKE- and the ISI-Distance between two or more spike trains:

 

Kreuz_et_al_SPIKE-distance.zip

(Last updated version 1.2: April 29, 2013;  Copyright:  Thomas Kreuz, Nebojsa Bozanic)

 Important note: This code will not be updated anymore. In principle it is replaced by the program 'SPIKY_loop' which is included in the zip-file of the new graphical user interface SPIKY which makes use of the same MEX-files. The Beta-version of these source codes can now be found here.


This version uses MEX- instead of m-files for the most time-consuming calculations. This reduces the computational cost considerably (on my notebook the code running all the examples below was faster by a factor of 85). It also uses a new output-structure 'results' which allows easy access not only to the overall dissimilarity value but also to the measure profiles and the pairwise distance matrices for all the selected measures.

 


 

Note of warning (see also 'Notice' at the end):

This program calculates the SPIKE- and the ISI-Distance between two or more spike trains. It also allows the time-resolved visualization via measure profiles and the spatial visualization via dissimilarity matrices (plots and figures). This, however, makes it a rather complex program so it might happen initially that there are some bugs or that the program is not as error-tolerant as it should be. Whenever you encounter any problem please provide feedback and we will try to resolve the problem and improve the source codes accordingly.

The name of the zip-file will always include the date of release which will allow users to ensure that they have the most recent version. As the development of the measure evolved, more and more applications came to mind and were integrated, and this might not be the end of it. Thus it might be worthwhile to check for updates now and then, in particular in these early alpha development stages. If you want to be informed about possible future releases, please send an email to Thomas Kreuz at "thomas.kreuz (at) cnr.it". You will then be included in a mailing list.

 

Description:

The zip-file contains many m- and eight MEX-files which all should be stored in the same directory. Before the first run of the program the MEX-files should be compiled. To do so please run from this directory the program 'Compile_MEX.m'.

 

Distances_Main_Demo: Demo program, used to create all the figures below

f_get_data_demo: Main function (sets spike-matrix as well as the parameters describing these data), called by 'Distances_Main_Demo'

Distances_Main: Main program; loads or creates data, sets parameters, and calls main function 'f_distances'

f_get_data: Main function (sets spike-matrix as well as the parameters describing these data), called by 'Distances_Main'

f_distances_MEX: Main function (calculates distances and does the plotting), called by 'Distances_Main' or 'Distances_Main_Demo', uses MEX-functions for computationally expensive loops

This function calls the following functions (included):

f_measure_profiles: Function that creates individual measure profiles, called by 'f_distances_MEX'

f_moving_average (& variants): Functions to generate the (weighted) moving average of a measure profile, depending on the measure profile the non-causal or the causal variant is used

f_compute_gauss_smooth: Function that creates a smoothed version of the PSTH

(Copyright: Jude Mitchell, http://www.snl.salk.edu/~jude/sfn2008/index.html)

as well as a few more functions (some of them from Matlab toolboxes)

 

To get started, try either Distances_Main_Demo or Distances_Main. However, before the first start of the program a few variables (see also below) defining the work environment should be set:

f_para.imagespath: Path where images (ps) will be stored
f_para.moviespath: Path where movies (avi) will be stored
f_para.matpath:      Path where Matlab files (mat) will be stored
 

 

Main function f_distances_MEX:

Input:

Matrix spikes with two or more spike trains (if trains have different numbers of spikes, fill with zeros)
Structure d_para that describe the data (see below)
Structure f_para that determines the appearance of the figures and the movie

Output:

Structure 'results' with these fields:

all of these variables have one entry for each measure selected via the variable 'f_para.subplot_posi' (see below):

measures: Name of selected measures (helps to identify the order within all other variables)

dissimilarity_profiles: cell with overall dissimilarity profiles obtained by averaging over spike train pairs (either piecewise constant (size  'num_isi') or sampled (length 'len'), see additional figure 1 below for an example)

distance_matrices: pairwise distance matrices, obtained by averaging over time

overall_dissimilarities: overall dissimilarities, obtained by averaging over both spike trains and time

In case a sampled measure was selected:

time: sampling values of length 'len' (defines x-axis for sampled measures)

In case a piecewise constant measure was selected:

isi: length of successive isi within the pooled spike train, length 'num_isi' (defines x-axis for piecewise constant measures)

 

For the piecewise constant profiles (here for the first) the function f_pico can be used to obtain the average value as well as x- and y-vectors for plotting:

[overall_dissimilarity,plot_x_values,plot_y_values] = f_pico(results.isi,results.dissimilarity_profiles{1},d_para.dts,d_para.tmin);

 


Example call:

After all other parameters have been set to standard values, i.e.

d_para=d_para_default;
f_para=f_para_default;


d_para.tmin=0; d_para.tmax=101; d_para.dts=0.001;
num_trains=4; num_spikes=100;
spikes=zeros(num_trains,num_spikes);
spikes(1,1:num_spikes)=1:num_spikes;
for rc=2:num_trains
      spikes(rc,1:num_spikes-rc)=spikes(1,1:num_spikes-rc)+(rand(1,num_spikes-rc)-0.5);
end
results=
f_distances_MEX(spikes,d_para,f_para);

 

Parameters:

The program uses three sets of parameters : Each of these are stored in a Matlab-struct.

1. d_para:   Parameters that describe the data (set within the function f_get_data and then used in the call of the main function f_distances_MEX from the main program)

2. f_para:    Parameters that determine the appearance of the figures and the movie (set in the main program and then used in the call of the main function f_distances_MEX from the main program)

3. s_para:   Parameters that determine the appearance of the individual subplots in the figure (set in the function f_distances_MEX and then used in the call of the function f_measure_profiles)

 

Here is a description of the individual structures:

1. Structure 'd_para': parameters that describe the data

tmin: Beginning of recording
tmax: End of recording
dts: Sampling interval, precision of spike times [!!! Please take care that this value is not larger than the actual sampling size, otherwise two spikes can occur at the same time instant and this can lead to problems in the algorithm !!!]
all_train_group_names: Names of spike train groups
all_train_group_sizes: Sizes of respective spike train groups
select_train_mode: Selection of spike trains (1-all,2-selected groups,3-selected trains)
select_train_groups: Selected spike train groups (if 'select_train_mode==2')
select_trains: Selected spike trains (if 'select_train_mode==3')
select_averages: One or more continuous intervals for selective temporal averaging
trigger_averages: Non-continuous time-instants for triggered temporal averaging, external (e.g. certain stimulus properties) but also internal (e.g. certain event times)
markers: Relevant time instants
markers2: Even more relevant time instants
separators: Relevant seperations between groups of spike trains
separators2: Even more relevant seperations between groups of spike trains
interval_separators: Edges of subsections
interval_strings: Captions for subsections

comment_string: Comment on the data, will be used in file and figure names

 

2. Structure 'f_para': parameters that determine the appearance of the figures (and the movie)

imagespath: Path where images (postscript) will be stored
moviespath: Path where movies (avi) will be stored
matpath: Path where Matlab file (mat) will be stored
filename: Name under which images, movies and Matlab files will be stored
title_string: Appears in the figure titles
saving: Saving to postscript file? (0-no,1-yes)
print_mode: Saving to postscript file? (0-no,1-yes)
publication: Omits otherwise helpful information, such as additional comments (0-no,1-yes)
comment_string: Additional comment on the example, will be used in file and figure names
num_fig: Number of figure
pos_fig: Position of figure
font_size: Font size of labels (and headlines)
multi_figure: Open many figures (0-no,1-yes)
timeunit_string: string with time unit, used in labels
xfact: Conversion of time unit
ma_mode: Moving average mode: (0-no,1-only,2-both)
spike_mao: Order of moving average (pooled ISI)
time_mao: Order of moving average (time)
dtm: Sampling of measure profile, setting this to larger values allows downsampling in order to facilitate memory management
mov_step: Step size for movie frames
mov_frames_per_second: Well, frames per second
mov_num_average_frames: Number of mov_frames the averages are shown at the end of the movie (if this is too small they are hardly visable)
mov_frames: Selection of individual frames to be shown in the movie (replaces use of mov_step as step size if non-empty)
plot_mode: +1:vs time,+2-different measures and cuts,+4-different measures,+8-different cuts,+16:different cuts-Movie (binary addition allows all combinations)
norm_mode: normalization of averaged bivariate measure profiles (1-Absolute maximum value 'one',2-Overall maximum value,3-Individual maximum value)
color_norm_mode: normalization of pairwise color matrices (1-Absolute maximum value,2-Overall maximum value,3-Each one individually)
block_matrices: Allows tracing the overall synchronization among groups of spike trains (0-no,1-yes)
dendrograms: Cluster trees from distance matrices (0-no,1-yes)
dendro_color_mode: coloring of dendrograms (0-no,1-groups,2-trains)
subplot_size: Vector with relative size of subplots
subplot_posi: Vector with order of subplots
 

3. Structure 's_para': parameters that describe the appearance of the individual subplots in the figures (measure time profiles)

window_mode: time interval selection (1-all (recording limits),2:outer spikes,3-inner spikes,4-smaller analysis window)
colors: Colors for 1:Stim/STs/ISI,2:Stim/STs/ISI-outside,3:Mean values,4:Mean values-outside,5:Moving averages,6:Moving averages-outside
font_size: Font size of labels (and headlines)
line_width: Line width
line_style: Line style (1-,2:,3-.,4--,5:none)
nma: Selected moving averages (depends on f_para.ma_mode)
spike_mao: Order of moving average (pooled ISI)
time_mao: Order of moving average (time)
causal: determines kind of moving average, set automatically for each measure (0-no,1-yes)
itmin: Beginning of recording
itmax: End of recording
wmin: Beginning of selected window (if window_mode==4)
wmax: End of selected window (if window_mode==4)
num_subplots: depends on f_para.subplot_posi, set automatically
xl: x-limits for plotting, set automatically
yl: y-limits for plotting, set automatically
dtm: equals f_para.dtm
plot_mode: equals f_para.plot_mode
 

 

Some of the more important parameters are here explained in more detail:

f_para.plot_mode  - This variable selects the types of plots (and movie)

+1: time-resolved measure profile (see Figs. 2 and 3 of the 2013 paper and below)

+2: different measures and cuts (see Figs. 5 and 6 of the 2013 paper and below)

+4: different measures (creates one figures of different cuts for each measure)

+8: different cuts (basically screenshots of the movie, see Fig. 9 in the 2013 paper and below)

+16: different cuts-Movie (basically a succession of plots, see Fig. 9 in the 2013 paper and below as well as the supplementary movie)

Binary addition allows all combinations of these plots/movies.

 

For the measure profiles (created by adding +1 to f_para.plot_mode):

f_para.subplot_posi and f_para.subplot_size - These variables allow you to choose the selection and order the subplots as well as their relative size.

Use f_para.subplot_posi to select the order of the subplots. Use 0 if a subplot is not needed. Make sure that the range from 1 to the desired number of subplots is covered, no (positive) number should appear twice.

Use f_para.subplot_size to select the relative size of the selected subplots.

The first variable f_para.subplot_posi refers to the standard order of the subplots which is defined below.

For each of the following variables there will be an entry in the output variable (in the order of the list below.

The second variable f_para.subplot_size then refers to the subplots selected by "f_para.subplot_posi" and gives their relative sizes from top to bottom (its length should be equal to the maximum number in "f_para.subplot_posi").

The subplot "Stimulus" allows you to add any representation of the stimulus that triggered the spike trains (put 0 if not needed).

 

In f_distances_MEX the function f_measure_profiles is called once for each selected measure profile. This is a list of the measures provided.

Stimulus, Spike trains, Classic, Rate, ISI:

1             2        3       4         5     6       7      8            9    10
Stimulus  Spikes  PSTH  GPSTH  ISI  <ISI>  Rate  <Rate>   Ia   I^a

SPIKE:

11  12       13       14          15    16           17      18         19         20         21      22
Sa  S^a     Sa_r   S_r^a     Sa_f   S_f^a     Sta     St^a     Sta_r  St_r^a      Sta_f   St_f^a

From 9 to 22 these are six pairs where the first entry (Xa) refers to all pairwise measure profiles and the second entry (X^a) to the averaged bivariate measure profiles.

Piecewise constant and sampled variants:

The dissimilarity profile S (t) of the SPIKE-distance is piecewise linear rather than piecewise constant (as is the case for the ISI-distance). Therefore, when the localized visualization is desired, a new value has to be calculated for each sampling point and not just once per each interval in the pooled spike train.

In cases where the distance value itself is sufficient, the short computation time can be even further decreased by representing each interval by the value of its center and weighting it by its length. This is not only faster, but it actually gives the exact result, whereas the time-resolved calculation is a very good approximation only for sufficiently small sampling intervals dt (imagine the example of a rectangular function, at some point any sampled representation has to cut the right angle). The dissimilarity profile S_r (t) of the real-time SPIKE-distance is hyperbolic and not linear but here also the exact result can be obtained by piecewise integration over all intervals of the pooled spike train.

 

Measures 9-16 are piecewise constant (one value for each interval in the pooled spike train), while measures 17-22 are sampled with the sampling interval given by the variable f_para.dtm.

9-10: ISI-distance (per definition piecewise constant)

11-12: SPIKE-distance (piecewise constant)

13-14: Real-time SPIKE-distance (piecewise constant)

15-16: Future SPIKE-distance (piecewise constant)

17-18: SPIKE-distance (sampled with variable 'time')

19-20: Real-time SPIKE-distance (sampled with variable 'time')

21-22: Future SPIKE-distance (sampled with variable 'time')

For each selected variable a value in the variable 'mean_values' is returned (the order is defined by f_para.subplot_posi).

 

There exist a number of different possibilities how the first and the last inter-spike interval are defined. For longer spike trains with sufficient statistics these differences become rather irrelevant.

s_para.window_mode (1-all, 2-outer spikes, 3-inner spikes, 4-chosen analysis window)

1. All - time profiles are calculated within this interval (which should cover the spike train such that in each spike train the first interval is from d_para.tmin to the first spike, the last interval is from the last spike to d_para.tmax). This corresponds to the case where spike trains were recorded for a certain time. The variable d_para.min then denotes the onset of the recording, d_para.tmax its end.

2. Outer spikes - time profiles are averaged from the very first spike till the very last spike.

3. Inner spikes - time profiles are calculated on the common support, i.e., from the latest first spike till the earliest last spike. This option does not make much sense if the spike trains do not have overlapping support (i.e., if one spike train starts after the other has already ended).

4. Window - restricts the averaging to a smaller window [s_para.wmin s_para.wmax] to be defined. These values should be within the range [tmin tmax].

s_para.wmin and s_para.wmax - beginning and end of the analysis window (for window_mode=4).

 

s_para.colors - colors used to display the stimulus/spike trains/ISI, the mean values and the moving averages both for the whole time profile and the selected analysis window (for window_mode>1)

f_para.publication (0-no, 1-yes) - If yes, the figure will have a title with all the selected distance values and the postscript-file will be a large landscape print, otherwise no title and a smaller portrait file.

 

For the movies (created by adding +16 to f_para.plot_mode):

The movies will always append individual figures in this order:

1. Individual frames (either individually selected or sampled with a certain step size)

2. Selected averages

3. Triggered averages

 

There are two possibilities to select individual frames in the first part of the movie.

1. f_para.mov_frames - Selection of individual frames to be shown in the movie (replaces use of mov_step as step size if non-empty)

2. f_para.mov_step - Step size for movie frames (will only be used if  f_para.mov_frames is empty).

 

f_para.mov_num_average_frames - Number of frames the selected and triggerad averages are shown at the end of the movie (in case the first part of the movie contains a sequence of snapshots with high temporal resolution and f_para.mov_frames_per_second is rather large this number has to be rather large (e.g., of the order of f_para.mov_frames_per_second) in order to make the averages visible)

f_para.block_matrices: Allows tracing the overall synchronization among groups of spike trains (0-no,1-yes) (see Fig. 9 of the 2013 paper and below)

f_para.dendrograms: Cluster trees from distance matrices (0-no,1-yes), (see Figs. 8 and 9 of the 2013 paper and below)

d_para.select_averages: Continuous intervals for selective temporal averaging (see Fig. 7 of the 2013 paper and below)

d_para.trigger_averages: Non-continuous time-instants for triggered temporal averaging, external (e.g. certain stimulus properties) but also internal (e.g. certain event times) (see Fig. 8 of the 2013 paper and below)

 

To control the position and size of the matrix or dendrogram subplots the variable supos is set in f_distances_MEX. Standard settings for up to 8 different subplots are provided, however, further adjustments might be necessary.

 

In case you want to create your own plots the pairwise results are stored in the four-dimensional variable 'mov_mat' with dimensions [#(selected measures)] * [#(selected frames) + #(selected averages) + #(selected triggers)] *[#(spike trains)[ * [#(spike trains)]. The respective variables for the second dimensions are num_frame_select, num_select_averages, and num_trigger_averages.

 

Memory management:

For large numbers of spike trains with many spikes the interval is cut into smaller segments, and the averaging over all pairs of spike trains is performed for each segment separately. During calculations this is indicated on the screen by 'run_info = #run  #(overall runs)'.

In the beginning of the function the variable 'max_memo_init' is set. This variable should be big enough to hold the basic matrices for one run of the segment loop but small enough to not run out of memory.

 


===================================================

In the following you will find an explanation of the most relevant parameter settings needed in order to reproduce some of the plots from the 2013 paper:

To reproduce these figures run 'Distance_Main_Demo' and select the figures by setting the variable 'examples' accordingly

 

 

Fig. 2a:  SPIKE-distance applied to a bivariate example with varying spike match

(example 21)

f_para.subplot_posi=[0 1   0 0 0 0 0 0   0 0     0 0 0 0 0 0   0 2 0 0 0 0]; % Vector with order of subplots
f_para.plot_mode=1; % +1:vs time,+2-different measures and cuts,+4-different measures,+8-different cuts,+16:different cuts-Movie
f_para.comment_string=['Fig2a']; % Additional comment on the example, will be used in file and figure names

 

 

 

 

Fig. 2b:  SPIKE-distance applied to a multivariate example

(example 22)

f_para.subplot_posi=[0 1   0 0 0 0 0 0   0 0     0 0 0 0 0 0   0 2 0 0 0 0]; % Vector with order of subplots
f_para.plot_mode=1; % +1:vs time,+2-different measures and cuts,+4-different measures,+8-different cuts,+16:different cuts-Movie
f_para.comment_string=['Fig2b']; % Additional comment on the example, will be used in file and figure names

 

 

 

 

Fig. 3a:  Real-time SPIKE-distance applied to a bivariate example with varying spike match

(example 31)

f_para.subplot_posi=[0 1   0 0 0 0 0 0   0 0     0 0 0 0 0 0   0 0 0 2 0 0]; % Vector with order of subplots
f_para.plot_mode=1; % +1:vs time,+2-different measures and cuts,+4-different measures,+8-different cuts,+16:different cuts-Movie
f_para.ma_mode=2; % Moving average mode: (0-no,1-only,2-both)
f_para.time_mao=120; % Order of moving average (time)
f_para.comment_string=['Fig3a_Realtime']; % Additional comment on the example, will be used in file and figure names

 

 

 

 

Fig. 3b:  SPIKE-distance applied to a multivariate example

(example 32)

f_para.subplot_posi=[0 1   0 0 0 0 0 0   0 0     0 0 0 0 0 0   0 0 0 2 0 0]; % Vector with order of subplots
f_para.plot_mode=1; % +1:vs time,+2-different measures and cuts,+4-different measures,+8-different cuts,+16:different cuts-Movie
f_para.ma_mode=2; % Moving average mode: (0-no,1-only,2-both)
f_para.time_mao=190; % Order of moving average (time)
f_para.comment_string=['Fig3b_Realtime']; % Additional comment on the example, will be used in file and figure names

 

 

 

 

Fig. 4:  Real-time SPIKE-distance: Peaks during reliable spiking events
are not spurious

(example 41)

f_para.subplot_posi=[0 1   0 0 0 0 0 0   0 0     0 0 0 0 0 0   0 0 0 2 0 0]; % Vector with order of subplots
f_para.plot_mode=1; % +1:vs time,+2-different measures and cuts,+4-different measures,+8-different cuts,+16:different cuts-Movie
f_para.comment_string=['Fig4']; % Additional comment on the example, will be used in file and figure names

 

 

 

 

Fig. 5A:  Instantaneous clustering for artificially generated spike trains

(example 51)

d_para.interval_separators=500:500:d_para.tmax-500; % Edges of subsections
f_para.subplot_posi=[0 1   0 0 0 0 0 0   0 0     0 0 0 0 0 0   0 2 0 3 0 0]; % Vector with order of subplots
f_para.plot_mode=2; % +1:vs time,+2-different measures and cuts,+4-different measures,+8-different cuts,+16:different cuts-Movie
f_para.mov_frames=[250 750 1190 1510]; % Selection of individual frames
f_para.xfact=1000; % Conversion of time unit
f_para.timeunit_string='[s]'; % Time unit, used in labels
f_para.comment_string=['Fig5a']; % Additional comment on the example, will be used in file and figure names

 

 

 

 

Fig. 5B:  Further examples of instantaneous clustering for artificially generated spike trains

(example 61)

d_para.interval_separators=500:500:d_para.tmax-500; % Edges of subsections
f_para.subplot_posi=[0 1   0 0 0 0 0 0   0 0     0 0 0 0 0 0   0 2 0 3 0 0]; % Vector with order of subplots
f_para.plot_mode=2; % +1:vs time,+2-different measures and cuts,+4-different measures,+8-different cuts,+16:different cuts-Movie
f_para.mov_frames=[250 750 1250 1750]; % Selection of individual frames f_para.xfact=1000; % Conversion of time unit
f_para.timeunit_string='[s]'; % Time unit, used in labels
f_para.comment_string=['Fig5b']; % Additional comment on the example, will be used in file and figure names

 

 

 

 

Fig. 6:  Selective temporal averaging

(example 71)

d_para.interval_separators=500:500:d_para.tmax-500; % Edges of subsections
d_para.select_averages={[0 500];[1000 2000];[2000 2500 3000 3500];[0 4000]}; % Selected average over different intervals
f_para.subplot_posi=[0 1   0 0 0 0 0 0   0 0     0 0 0 0 0 0   0 2 0 3 0 0]; % Vector with order of subplots
f_para.plot_mode=2; % +1:vs time,+2-different measures and cuts,+4-different measures,+8-different cuts,+16:different cuts-Movie
f_para.xfact=1000; % Conversion of time unit
f_para.timeunit_string='[s]'; % Time unit, used in labels
f_para.comment_string=['Fig6']; % Additional comment on the example, will be used in file and figure names

 

 

 

 

Fig. 7a:  Triggered temporal averaging

(example 81)

d_para.select_averages={[d_para.tmin d_para.tmax]};  % One continuous intervals for selective temporal averaging
f_para.subplot_posi=[0 1   0 0 0 0 0 0   0 0     0 0 0 0 0 0   0 2 0 0 0 0]; % Vector with order of subplots
f_para.plot_mode=8; % +1:vs time,+2-different measures and cuts,+4-different measures,+8-different cuts,+16:different cuts-Movie
f_para.dendrograms=1; % Cluster trees from distance matrices (0-no,1-yes)
f_para.comment_string=['Fig7a']; % Additional comment on the example, will be used in file and figure names

 

 

 

 

Fig. 7b:  SPIKE-distance applied to a bivariate example with varying spike match

(example 81)

d_para.trigger_averages{1}=spikes(1,1:num_spikes(1)); % Triggered averaging over all time instants when a certain neuron fires
f_para.subplot_posi=[0 1   0 0 0 0 0 0   0 0     0 0 0 0 0 0   0 2 0 0 0 0]; % Vector with order of subplots

f_para.plot_mode=8; % +1:vs time,+2-different measures and cuts,+4-different measures,+8-different cuts,+16:different cuts-Movie
f_para.dendrograms=1; % Cluster trees from distance matrices (0-no,1-yes)
f_para.comment_string=['Fig7b']; % Additional comment on the example, will be used in file and figure names

 

 

 

 

Fig. 8:  Screenshot from the supplementary movie

(example 91)

d_para.markers=[500:500:d_para.tmax-500]; % Relevant time instants
d_para.all_train_group_names={'G1';'G2';'G3';'G4'}; % Names of spike train groups
d_para.all_train_group_sizes=[10 10 10 10]; % Sizes of respective spike train groups
d_para.select_averages={500+[0 500 1000 1500]}; % Two continuous intervals for selective temporal averaging
f_para.subplot_posi=[0 1   0 0 0 0 0 0   0 0     0 0 0 0 0 0   0 2 0 3 0 0]; % Vector with order of subplots
f_para.plot_mode=8; % +1:vs time,+2-different measures and cuts,+4-different measures,+8-different cuts,+16:different cuts-Movie
f_para.block_matrices=1; % Allows tracing the overall synchronization among groups of spike trains (0-no,1-yes)
f_para.dendrograms=1; % Cluster trees from distance matrices (0-no,1-yes)
f_para.dendro_color_mode=1; % Coloring of dendrograms (0-no,1-groups,2-trains)
f_para.comment_string=['Fig8']; % Additional comment on the example, will be used in file and figure names

 

 

 

 

Additional Figure 1:  Comparison of piecewise constant and sampled dissimilarity profile

(example 21)

f_para.subplot_posi=[0 1   0 0 0 0 0 0   0 0     0 2 0 0 0 0   0 2 0 0 0 0]; % Vector with order of subplots

f_para.subplot_size=[1 1.5]; % Vector with relative size of subplots
f_para.plot_mode=1; % +1:vs time,+2-different measures and cuts,+4-different measures,+8-different cuts,+16:different cuts-Movie
f_para.publication=0; % Omits otherwise helpful information, such as additional comments (0-no,1-yes)

[The dissimilarity profile with the normal line width is sampled (with one value each sample point), whereas the bold dissimilarity profile is piecewise constant (with one value for each interval in the pooled spike train). This profile is obtained by representing each interval by the value of its center and weighting it by its length. This is not only faster, but it actually gives the exact result, whereas the sampled calculation is a very good approximation only for sufficiently small sampling intervals (imagine the example of a rectangular function, at some point any sampled representation has to cut the right angle). On the other hand, the sampled profile has a higher temporal resolution in the visualization. If this is not essential or if the averaged values are sufficient, the calculation of the piecewise constant distances is preferable.

 

Reminder:

9-10: ISI-distance (per definition piecewise constant)

11-12: SPIKE-distance (piecewise constant)

13-14: Real-time SPIKE-distance (piecewise constant)

15-16:Future SPIKE-distance (piecewise constant)

17-18: SPIKE-distance (sampled with variable 'time')

19-20: Real-time SPIKE-distance (sampled with variable 'time')

21-22: Future SPIKE-distance (sampled with variable 'time')

In the source codes the sampled variables can be recognized by the ending '_t'.

Regarding the plotting parameters: The sampled profile was generated because f_para.subplot_posi(18)>0 and the piecewise constant profile was generated because f_para.subplot_posi(12)>0. Since both values were 2 they were plotted in the same subplot. And since f_para.subplot_size was set to [1 1.5], the second subplot is 1.5 times larger than the first.]

 

 

 

 

Additional Figure 2:  Poisson with different firing rate - Distance value only depends on distance from diagonal (distances are scale-free)

(example 121)

d_para.select_averages={[d_para.tmin d_para.tmax]}; % Continuous intervals for selective temporal averaging
f_para.subplot_posi=[0 1   0 0 0 0 0 0   0 0     0 2 0 3 0 0   0 0 0 0 0 0]; % Vector with order of subplots
f_para.plot_mode=8; % +1:vs time,+2-different measures and cuts,+4-different measures,+8-different cuts,+16:different cuts-Movie
f_para.color_norm_mode=3; % normalization of pairwise color matrices (1-Absolute maximum value,2-Overall maximum value,3-Each one individually)

 

 

 

Additional Figure 3:  Poisson with different firing rate - Distance value only depends on distance from diagonal (distances are scale-free)

(example 121)

d_para.select_averages={[d_para.tmin d_para.tmax]}; % Continuous intervals for selective temporal averaging
f_para.subplot_posi=[0 1   0 0 0 0 0 0   0 0     0 2 0 3 0 0   0 0 0 0 0 0]; % Vector with order of subplots
f_para.plot_mode=8; % +1:vs time,+2-different measures and cuts,+4-different measures,+8-different cuts,+16:different cuts-Movie
f_para.color_norm_mode=3; % normalization of pairwise color matrices (1-Absolute maximum value,2-Overall maximum value,3-Each one individually)

 

 

 

Additional Figure 4:  Expectation values for Poisson spike trains with identical rate (but these values do not depend on the rate)

(example 131)

f_para.subplot_posi=[0 1   0 0 0 0 0 0   0 0     0 2 0 3 0 0   0 0 0 0 0 0]; % Vector with order of subplots
f_para.plot_mode=1; % +1:vs time,+2-different measures and cuts,+4-different measures,+8-different cuts,+16:different cuts-Movie
f_para.publication=0; % Omits otherwise helpful information, such as additional comments (0-no,1-yes)

[The firing rate of the Poisson spike trains is 1000 Hz so these are almost 1.000.000 spikes. It takes quite some time (~1h on my notebook) to run this example. The interval is cut into smaller segments, and the averaging over all pairs of spike trains is performed for each segment separately.]

 

See also:

Measuring spike train synchrony I:   SPIKY (graphical user interface)

Graphical user interface (Matlab) which can be used to calculate and visualize both the SPIKE- and the ISI-distance between two (or more) spike trains. This is an extension and update of the code in II (there you also find the links to the relevant articles).

Measuring spike train synchrony II:   cSPIKE

This Matlab command line software utilizes MEX files with C++ back ends for improved performance and abstraction. It offers the same basic functionality as SPIKY and PySpike and can be easily incorporated as part of your own Matlab program. You can also use cSPIKE as a stand-alone program to calculate the ISI-distance, SPIKE-distance and SPIKE synchronization measures and many more for your spike train sets.

Measuring spike train synchrony IV:   Event Synchronization

Matlab-Code to measure the event synchronization and the event delay between two given spike trains

For a detailed description of the method please refer to:

Quian Quiroga R, Kreuz T, and Grassberger P:
Event Synchronization: A simple and fast method to measure synchronicity and time delay patterns.
Phys.Rev. E 66, 041904 (2002) [PDF].

Measuring spike train synchrony V:   Directionality

Matlab code to calculate the directionality measure L between two given spike trains (or between two continuous datasets or between a spike train and a continuous dataset)

For a detailed description of the method please refer to:

Andrzejak RG, Kreuz T:
Characterizing unidirectional couplings between point processes and flows.
European Physics Letters 96, 50012 (2011) [PDF].

Measuring spike train synchrony VI:   van Rossum distance and multi-neuron extension

Matlab codes to calculate the spike train metric by van Rossum and the multi-neuron extension by Houghton and Sen.

For a detailed description of the method please refer to:

Houghton C, Kreuz T:
On the efficient calculation of van Rossum distances.
Network: Computation in Neural Systems 23, 48 (2012) [PDF].

Measuring spike train synchrony VII:  Victor-Purpura distance and multi-neuron extension

Matlab codes to calculate the spike train metric by Victor-Purpura and the multi-neuron extension by Victor-Purpura-Aronov.

(Homepage of Prof. Jonathan D. Victor, Cornell, NY, USA)

See also:

Matlab code for the Victor-Purpura distance which in addition calculates the percentage of spikes that have been matched by a time shift as well as the average time shift

For a detailed description of the algorithm please refer to:

Chicharro D, Kreuz T, Andrzejak RG:
What can spike train distances tell us about the neural code?
J Neurosci Methods 199, 146 (2011) [PDF].


For questions and comments please contact me at "thomas.kreuz (at) cnr.it".

 


These are some keywords to help people to find this page.

spike train analysis, measuring spike train synchrony, monitoring spike train synchrony, ISI-distance, SPIKE-distance, open Matlab source codes, free download


Notices

FOR SCIENTIFIC USE ONLY

These codes are free of charge for research and education purposes only. Any commercial or military use of this software is prohibited.

The software on this site is provided "as-is," without any expressed or implied warranty. In no event am I or my host institution liable for any damages arising from the use of the software. Since it is distributed for free, I do also not take responsibility for any eventual error in it.


activate the plugin aktuelle Homepage