Source Code - cSPIKE


cSPIKE-Logo

cSPIKE is an easy to use spike train analysis software. It is run on Matlab command line and it uses MEX files with C++ backends for speed. cSPIKE implements functions such as ISI-distance, SPIKE-distance, SPIKe synchronization and their adaptive variants as well as basic plot functions for plotting spike trains and profiles.

The latest version 1.3 also contains the complementary directional method Spike Train Order and SPIKE-order (for details see Ref. [1]) and two algorithms for population analysis (for details see Ref. [2]). See below for examples.

[1] Kreuz T, Satuvuori E, Pofahl M, Mulansky M:
Leaders and followers: Quantifying consistency in spatio-temporal propagation patterns
New J. Phys., 19, 043028 [PDF] and arXiv [PDF] (2017).

[2] Satuvuori E, Mulansky M, Daffertshofer A, Kreuz T:
Using spike train distances to identify the most discriminative neuronal subpopulation JNeurosci Methods, 308, 354 [PDF] and arXiv [PDF] (2018).


Download:

New     cSPIKE.zip        New      
(Version 1.3, January 31, 2019)


Installation:


cSPIKE allows easy incorporation as a part of a Matlab program. All you need to do is, in this order, extract the files into your project folder, compile the MEX files and set the work path (see below).




1. Compiling MEX files:

In order to build the MEX files and the C++ backends you need a C++ compiler in your Matlab. Here are few links that could help:

What you need to build MEX files
About compilers for Matlab R2016a

Once you have your C++ compiler go to the folder cSPIKEmex and run the script MEX_compile. This will compile your MEX files and C++ backends in the folder. The process may take few minutes.

As a summary you need to:

1) Extract the cSPIKE to folder of your choise
2) Install a C++ compiler for your Matlab
3) Compile your MEX files and C++ backends
4) Add the folder cSPIKEmex to your work path
5) Get to cSPIKE main folder and create your first SpikeTrainSet

Known problems:

1. In case you are using Matlab 2018b or newer and you get some error messages with the MEX-compilation please check this link for assistance:

https://se.mathworks.com/help/matlab/matlab_external/what-if-i-do-not-upgrade.html


2. If you get an error message while compiling that resembles the following:

MEX_compile

Building with 'g++'.
Error using mex
In file included from /usr/include/c++/4.9/bits/stl_algobase.h:61:0,
                 from /usr/include/c++/4.9/vector:60,
                 from /home/jr3830/Insync/proj/lesica/cspike/cSPIKEmex/mexAdaptiveISIDistance.cpp:12:
/usr/include/c++/4.9/bits/cpp_type_traits.h:206:12: error: redefinition of ‘struct std::__is_integer<short unsigned int>’
     struct __is_integer<unsigned short>
            ^
/usr/include/c++/4.9/bits/cpp_type_traits.h:184:12: error: previous definition of ‘struct std::__is_integer<short unsigned int>’
     struct __is_integer<char16_t>
            ^
In file included from /usr/include/c++/4.9/bits/move.h:57:0,
                 from /usr/include/c++/4.9/bits/stl_pair.h:59,
                 from /usr/include/c++/4.9/bits/stl_algobase.h:64,
                 from /usr/include/c++/4.9/vector:60,
                 from /home/jr3830/Insync/proj/lesica/cspike/cSPIKEmex/mexAdaptiveISIDistance.cpp:12:
/usr/include/c++/4.9/type_traits:218:12: error: redefinition of ‘struct std::__is_integral_helper<short unsigned int>’
     struct __is_integral_helper<unsigned short>
            ^
/usr/include/c++/4.9/type_traits:206:12: error: previous definition of ‘struct std::__is_integral_helper<short unsigned int>’
     struct __is_integral_helper<char16_t>
            ^
/usr/include/c++/4.9/type_traits:1711:12: error: redefinition of ‘struct std::__make_signed<short unsigned int>’
     struct __make_signed<char16_t> : __make_signed<uint_least16_t>
            ^
/usr/include/c++/4.9/type_traits:1688:12: error: previous definition of ‘struct std::__make_signed<short unsigned int>’
     struct __make_signed<unsigned short>
            ^
In file included from /usr/include/c++/4.9/bits/stl_bvector.h:1181:0,
                 from /usr/include/c++/4.9/vector:65,
                 from /home/jr3830/Insync/proj/lesica/cspike/cSPIKEmex/mexAdaptiveISIDistance.cpp:12:
/usr/include/c++/4.9/bits/functional_hash.h:113:3: error: redefinition of ‘struct std::hash<short unsigned int>’
   _Cxx_hashtable_define_trivial_hash(unsigned short)
   ^
/usr/include/c++/4.9/bits/functional_hash.h:95:3: error: previous definition of ‘struct std::hash<short unsigned int>’
   _Cxx_hashtable_define_trivial_hash(char16_t)
   ^

Error in MEX_compile (line 5)
mex -Dchar16_t=uint16_T mexAdaptiveISIDistance.cpp Spiketrain.cpp DataReader.cpp Spiketrains.cpp Pair.cpp ISIProfile.cpp
SPIKEProfile.cpp

You have a compiler incompatibility with the code. The key components in the error message are char16_t  and unsigned short. If these two appear in the error message there are a few datatypes you need to change manually.
This is the likely fix:


1.    In file “MEX_compile” remove from all lines the following parameter: “-Dchar16_t=uint16_T”
2.    In file "mexAdaptiveSPIKESynchroProfile.cpp" change in line 66  notion "int NOST = profile1.size();" with " const long unsigned int NOST = profile1.size();"

This should fix the issue without changing the functionality of the code.


2. Setting work path:

The MEX files are located in cSPIKEmex folder and in order to run cSPIKE you need to include that folder into your work path. The easiest way to do so is to call this script at the beginning of your program:

>> initializecSPIKE

This will not only add the MEX-folder but also two other folders containing the SPIKE-order (SpikeTrainOrder) and the population analysis (PopulationAnalysis) to the work path.

Now cSPIKE is ready to be used.



Using cSPIKE:

cSPIKE is built around the main class of SpikeTrainSet. The class requires your spike time data in the form of a cell array in which each cell is a series of spike times. Each array of spike times needs to be a double precision row vector. This cell array together with the recording start and end times form a SpikeTrainSet. You can create your SpikeTrainSet object as:

STS = SpikeTrainSet(MyCellArrayOfSpikes,recordingStartTime,recordingEndTime)

Once the spike train set object has been created you can start calling the methods.
Not all of cSPIKE functionalities are decribed here. Please refer to the documentation (in particular this link here) and papers for more information. Also you can see a short description of all the methods by:

>> help SpikeTrainSet



But below are few examples. The example data used here can be found in the main folder of the cSPIKE.

We have three spike trains with 15 evenly distributed random spikes from 0 to 20 in cell matrix ExampleSet. To use cSPIKE we create a SpikeTrainSet with this information and initialize cSPIKE:

>> STS = SpikeTrainSet(ExampleSet,0,20);
>> InitializecSPIKE;


We can calculate the ISI-distance for the whole data set:

>> STS.ISIdistance()

ans =

    0.5264


ISI-distance from time 13 to time 15:

>> STS.ISIdistance(13,15)

ans =

    0.5417


Adaptive ISI-Distance from time 13 to time 15 with automated threshold:

>> STS.AdaptiveISIdistance(13,15)

ans =

    0.5285


And with threshold of 0.1:

>> STS.AdaptiveISIdistance(13,15,0.1)

ans =

    0.5417


cSPIKE also includes fast MEX files for calculating pairwise distance matrices where each index is a spike train number. Because the distances are symmetric, the matrix is symmetric as well:

>> STS.SPIKEdistanceMatrix()

ans =

         0    0.2626    0.3191
    0.2626         0    0.3101
    0.3191    0.3101         0


You can also call the methods for the SpikeTrainSet by giving a method the SpikeTrainSet as first parameter. The following calls are identical:

>> STS.SPIKEdistance();
>> SPIKEdistance(STS);


>> STS.ISIdistance(13,15)
>> ISIdistance(STS,13,15)

In addition to the distance measures, the SpikeTrainSet provides the profiles for the ISI- and SPIKE-distances. Once called, the method returns a Profile-object with three main methods. You can plot the profile to current axis with the optional colour and opacity(alpha) given:

>> Profile = STS.SPIKEdistanceProfile(13,15)

Profile.Plot(colour, alpha);

>> Profile.Plot();



ISI-profile

You can also get separately the X-values of the profile and the Y-values of the profile as double arrays:

>> Profile.PlotProfileX

ans =

  Columns 1 through 4

   13.0000   13.1020   13.1020   13.5941

  Columns 5 through 8

   13.5941   13.9815   13.9815   14.1873

  Columns 9 through 10

   14.1873   15.0000


>> Profile.PlotProfileY

ans =

  Columns 1 through 4

    0.4532    0.4535    0.4061    0.4012

  Columns 5 through 8

    0.3869    0.3183    0.4433    0.3685

  Columns 9 through 10

    0.3330    0.0566


These are vectors of equal length with each index corresponding a point in the graph.

The SpikeTrainSet can also plot the spike trains to current axis with optional colour and spike width
of your choise:

STS.plotSpikeTrainSet(colour, width)

>> STS.plotSpikeTrainSet()

SpikeTrains

Spike Train Order

Recently we have added the
newly introduced method to sort spike trains from leader to follower to cSPIKE. The function call is very simple. In order to perform the basic analysis, the SpikeTrainSet is built as before. Then the analysis is called as

[initialIteration,optimalIteration,synf] = STS.SpikeTrainOrderWithSurrogates(numberOfSurrogates)

The optional parameter numberOfSurrogates lets you do the test based on a number of surrogates you want. The default value is 19 for 5% confidence.

The return values are two stuctures, where the output is stored and an array containing the F values of the surrogates. Please see the paper "Leaders and Followers: Quantifying consistency in spatio-temporal propagation patterns" for more details. Also the cSPIKE package contains an example for how to create an analysis using the structure in file SpikeTrainOrderExample.m. that produces the following figure.

Spike Train Order example


For more information on Spike Train Order and SPIKE-order please refer to:


Kreuz T, Satuvuori E, Pofahl M, Mulansky M:
Leaders and followers: Quantifying consistency in spatio-temporal propagation patterns
New J. Phys., 19, 043028 [PDF] and arXiv [PDF] (2017).


Population analysis

Our latest addition to cSPIKE is these recently introduced algorithms to find in a population of neurons among all possible subpopulations the one that distinguishes a given set of stimuli best. The analysis can be performed in two ways, either assuming a summed population coding or following the labelled line hypothesis. To get started have a look at the file 'population_main.m' in the folder 'PopulationAnalysis'. Please note that this part of cSPIKE is still in some kind of 'alpha mode'. The functions should work correctly but we did not yet have the time to perform software optimization and proper beta-testing. If you encounter any problems please just contact us.

For more information on the population analysis please refer to:


Satuvuori E, Mulansky M, Daffertshofer A, Kreuz T:
Using spike train distances to identify the most discriminative neuronal subpopulation JNeurosci Methods, 308, 354 [PDF] and arXiv [PDF] (2018).


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 III:   SPIKE- and ISI-Distance

Matlab codes to calculate both the SPIKE- and the 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].

Also submitted to the arXiv [PDF].

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

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 "eero.satuvuori (at) gmail.com".

 


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, SPIKE-synchronization, Adaptive ISI-distance,Adaptive SPIKE-distance, Adaptive SPIKE-synchronization, 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.

BSD license:

Copyright (c) 2016 Eero Satuvuori
All rights reserved.

Redistribution and use in source and binary forms, with or without modification, are permitted provided that the following conditions are met:
* Redistributions of source code must retain the above copyright notice, this list of conditions and the following disclaimer.
* Redistributions in binary form must reproduce the above copyright notice, this list of conditions and the following disclaimer in the documentation and/or other materials provided with the distribution.
* Neither the name of the author nor the names of its contributors may be used to endorse or promote products derived from this software without specific prior written permission.

THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL <COPYRIGHT HOLDER> BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
 


activate the plugin aktuelle Homepage