From 1b0421a207f5a1eb43b12bb18d5e783e753b739e Mon Sep 17 00:00:00 2001 From: Mario Mulansky Date: Tue, 21 Oct 2014 19:05:39 +0200 Subject: doc: distance matrix, improved example --- examples/distance_matrix.py | 33 +++++++++++++++++++++++++++++++++ 1 file changed, 33 insertions(+) create mode 100644 examples/distance_matrix.py (limited to 'examples/distance_matrix.py') diff --git a/examples/distance_matrix.py b/examples/distance_matrix.py new file mode 100644 index 0000000..38bd9c8 --- /dev/null +++ b/examples/distance_matrix.py @@ -0,0 +1,33 @@ +""" distance_matrix.py + +Simple example showing how to compute the isi distance matrix of a set of spike +trains. + +Copyright 2014, Mario Mulansky + +Distributed under the BSD License +""" + + +from __future__ import print_function + +import matplotlib.pyplot as plt + +import pyspike as spk + +# first load the data, interval ending time = 4000, start=0 (default) +spike_trains = spk.load_spike_trains_from_txt("PySpike_testdata.txt", 4000) + +print(len(spike_trains)) + +plt.figure() +isi_distance = spk.isi_distance_matrix(spike_trains) +plt.imshow(isi_distance, interpolation='none') +plt.title("ISI-distance") + +plt.figure() +spike_distance = spk.spike_distance_matrix(spike_trains) +plt.imshow(spike_distance, interpolation='none') +plt.title("SPIKE-distance") + +plt.show() -- cgit v1.2.3 From 4a295e6045abc7564a2e72d1a2173bf2b04c5950 Mon Sep 17 00:00:00 2001 From: Mario Mulansky Date: Thu, 23 Oct 2014 12:46:55 +0200 Subject: docs: added interval averaging explanations --- Readme.rst | 38 ++++++++++++++++++++++++++++---------- examples/distance_matrix.py | 4 ++-- examples/plot.py | 3 ++- pyspike/distances.py | 14 +++++++------- 4 files changed, 39 insertions(+), 20 deletions(-) (limited to 'examples/distance_matrix.py') diff --git a/Readme.rst b/Readme.rst index 4482d01..c6ded74 100644 --- a/Readme.rst +++ b/Readme.rst @@ -87,8 +87,8 @@ Both the ISI and the SPIKE distance computation require the presence of auxiliar # if you provide only a single value, it is interpreted as T_end, while T_start=0 spike_train = spk.add_auxiliary_spikes(spike_train, T_end) -Computing bi-variate distances ------------------------------- +Computing bi-variate distances profiles +--------------------------------------- **Important note:** @@ -124,12 +124,25 @@ The following code loads some exemplary spike trains, computes the dissimilarity plt.show() The ISI-profile is a piece-wise constant function, there the function :code:`isi_profile` returns an instance of the :code:`PieceWiseConstFunc` class. -As shown above, this class allows you to obtain arrays that can be used to plot the function with :code":`plt.plt`, but also to compute the average, which amounts to the final scalar ISI-distance. +As shown above, this class allows you to obtain arrays that can be used to plot the function with :code:`plt.plt`, but also to compute the time average, which amounts to the final scalar ISI-distance. +By default, the time average is computed for the whole :code:`PieceWiseConstFunc` function. +However, it is also possible to obtain the average of some interval by providing a pair of floats defining the start and end of the interval. +In the above example, the following code computes the ISI-distances obtained from averaging the ISI-profile over four different intervals: + +.. code:: python + + isi1 = isi_profil.avrg(interval=(0,1000)) + isi2 = isi_profil.avrg(interval=(1000,2000)) + isi3 = isi_profil.avrg(interval=(2000,3000)) + isi4 = isi_profil.avrg(interval=(3000,4000)) + If you are only interested in the scalar ISI-distance and not the profile, you can simly use: .. code:: python - isi_dist = spk.isi_distance(spike_trains[0], spike_trains[1]) + isi_dist = spk.isi_distance(spike_trains[0], spike_trains[1], interval) + +where :code:`interval` is optional, as above, and if omitted the ISI-distance is computed for the complete spike trains. Furthermore, PySpike provides the :code:`average_profile` function that can be used to compute the average profile of a list of given :code:`PieceWiseConstFunc` instances. @@ -176,14 +189,16 @@ This short example computes and plots the SPIKE-profile of the first two spike t In contrast to the ISI-profile, a SPIKE-profile is a piece-wise *linear* function and thusly represented by a :code:`PieceWiseLinFunc` object. Just like the :code:`PieceWiseconstFunc` for the ISI-profile, the :code:`PieceWiseLinFunc` provides a :code:`get_plottable_data` member function that returns array that can be used directly to plot the function. Furthermore, the :code:`avrg` member function returns the average of the profile defined as the overall SPIKE distance. +As above, you can provide an interval as a pair of floats to :code:`avrg` to specify the averaging interval if required. Again, you can use .. code:: python - spike_dist = spk.spike_distance(spike_trains[0], spike_trains[1]) + spike_dist = spk.spike_distance(spike_trains[0], spike_trains[1], interval) to compute the SPIKE distance directly, if you are not interested in the profile at all. +:code:`interval` is optional and defines the averaging interval, if neglected the whole spike train is used. Furthmore, you can use the :code:`average_profile` function to compute an average profile of a list of SPIKE-profiles: .. code:: python @@ -195,7 +210,7 @@ Furthmore, you can use the :code:`average_profile` function to compute an averag Computing multi-variate profiles and distances ---------------------------------- +---------------------------------------------- To compute the multi-variate ISI- or SPIKE-profile of a set of spike trains, you can compute all bi-variate profiles separately and then use the :code:`average_profile` function above. However, PySpike provides convenience functions for that purpose. @@ -210,11 +225,12 @@ The following example computes the multivariate ISI- and SPIKE-profile for a lis Both functions take an optional parameter :code:`indices`, a list of indices that allows to define the spike trains that should be used for the multi-variate profile. As before, if you are only interested in the distance values, and not in the profile, PySpike offers the functions: :code:`isi_distance_multi` and :code:`spike_distance_multi`, that return the scalar multi-variate ISI- and SPIKE-distance. +Both distance functions also accept an :code:`interval` parameter that can be used to specify the averaging interval as a pair of floats, if neglected the complete interval is used. Another option to address large sets of spike trains are distance matrices. Each entry in the distance matrix represents a bi-variate distance of the spike trains. Hence, the distance matrix is symmetric and has zero values at the diagonal. -The following example computes and plots the ISI- and SPIKE-distance matrix. +The following example computes and plots the ISI- and SPIKE-distance matrix, where for the latter one only the time interval T=0..1000 is used for the averaging. .. code:: python @@ -226,15 +242,17 @@ The following example computes and plots the ISI- and SPIKE-distance matrix. plt.title("ISI-distance") plt.figure() - spike_distance = spk.spike_distance_matrix(spike_trains) + spike_distance = spk.spike_distance_matrix(spike_trains, interval=(0,1000)) plt.imshow(spike_distance, interpolation='none') plt.title("SPIKE-distance") plt.show() -Averaging ---------- +Time Averages +------------- + + .. _ISI: http://www.scholarpedia.org/article/Measures_of_spike_train_synchrony#ISI-distance diff --git a/examples/distance_matrix.py b/examples/distance_matrix.py index 38bd9c8..142db2c 100644 --- a/examples/distance_matrix.py +++ b/examples/distance_matrix.py @@ -26,8 +26,8 @@ plt.imshow(isi_distance, interpolation='none') plt.title("ISI-distance") plt.figure() -spike_distance = spk.spike_distance_matrix(spike_trains) +spike_distance = spk.spike_distance_matrix(spike_trains, interval=(0, 1000)) plt.imshow(spike_distance, interpolation='none') -plt.title("SPIKE-distance") +plt.title("SPIKE-distance, T=0-1000") plt.show() diff --git a/examples/plot.py b/examples/plot.py index 59334c9..d32c464 100644 --- a/examples/plot.py +++ b/examples/plot.py @@ -1,6 +1,7 @@ """ plot.py -Simple example showing how to load and plot spike trains and their distances. +Simple example showing how to load and plot spike trains and their distance +profiles. Copyright 2014, Mario Mulansky diff --git a/pyspike/distances.py b/pyspike/distances.py index 5135b9b..34f7d78 100644 --- a/pyspike/distances.py +++ b/pyspike/distances.py @@ -57,14 +57,13 @@ def isi_distance(spikes1, spikes2, interval=None): isi-distance is the integral over the isi distance profile :math:`S_{isi}(t)`: - .. math:: I = \int_0^T S_{isi}(t) dt. + .. math:: I = \int_{T_0}^{T_1} S_{isi}(t) dt. :param spikes1: ordered array of spike times with auxiliary spikes. :param spikes2: ordered array of spike times with auxiliary spikes. - :param interval: averaging interval given as a pair of floats, if None - the average over the whole function is computed. + :param interval: averaging interval given as a pair of floats (T0, T1), + if None the average over the whole function is computed. :type interval: Pair of floats or None. - :returns: The isi-distance I. :rtype: double """ @@ -114,12 +113,13 @@ Falling back to slow python backend.") def spike_distance(spikes1, spikes2, interval=None): """ Computes the spike-distance S of the given spike trains. The spike-distance is the integral over the isi distance profile S_spike(t): - :math:`S = \int_^T S_spike(t) dt`. + + .. math:: S = \int_{T_0}^{T_1} S_{spike}(t) dt. :param spikes1: ordered array of spike times with auxiliary spikes. :param spikes2: ordered array of spike times with auxiliary spikes. - :param interval: averaging interval given as a pair of floats, if None - the average over the whole function is computed. + :param interval: averaging interval given as a pair of floats (T0, T1), + if None the average over the whole function is computed. :type interval: Pair of floats or None. :returns: The spike-distance. :rtype: double -- cgit v1.2.3 From f4266628dc89b8747923bbddcf1b7f6b13b701d8 Mon Sep 17 00:00:00 2001 From: Mario Mulansky Date: Tue, 20 Jan 2015 14:30:34 +0100 Subject: added spike sync desc to docs --- Readme.rst | 94 +++++++++++++++++++++++++++------------------ examples/distance_matrix.py | 5 +++ pyspike/distances.py | 41 ++++++++++++++++---- 3 files changed, 95 insertions(+), 45 deletions(-) (limited to 'examples/distance_matrix.py') diff --git a/Readme.rst b/Readme.rst index 4053995..bcf6fa9 100644 --- a/Readme.rst +++ b/Readme.rst @@ -5,8 +5,8 @@ PySpike :target: https://travis-ci.org/mariomulansky/PySpike PySpike is a Python library for the numerical analysis of spike train similarity. -Its core functionality is the implementation of the bivariate ISI_ and SPIKE_ distance [#]_ [#]_. -Additionally, it provides functions to compute multivariate SPIKE and ISI distances, as well as averaging and general spike train processing. +Its core functionality is the implementation of the bivariate ISI_ and SPIKE_ distance [#]_ [#]_ as well as SPIKE-Synchronization_ [#]_. +Additionally, it provides functions to compute multivariate profiles, distance matrices, 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). PySpike provides the same fundamental functionality as the SPIKY_ framework for Matlab, which additionally contains spike-train generators, more spike train distance measures and many visualization routines. @@ -17,6 +17,8 @@ All source codes are published under the BSD_License_. .. [#] Kreuz T, Chicharro D, Houghton C, Andrzejak RG, Mormann F, *Monitoring spike train synchrony.* J Neurophysiol 109, 1457 (2013) `[pdf] `_ +.. [#] Kreuz T, Mulansky M and Bozanic N, *SPIKY: A graphical user interface for monitoring spike train synchrony*, tbp (2015) + Requirements and Installation ----------------------------- @@ -162,32 +164,11 @@ If you are only interested in the scalar ISI-distance and not the profile, you c where :code:`interval` is optional, as above, and if omitted the ISI-distance is computed for the complete spike trains. -Furthermore, PySpike provides the :code:`average_profile` function that can be used to compute the average profile of a list of given :code:`PieceWiseConstFunc` instances. - -.. code:: python - - isi_profile1 = spk.isi_profile(spike_trains[0], spike_trains[1]) - isi_profile2 = spk.isi_profile(spike_trains[0], spike_trains[2]) - isi_profile3 = spk.isi_profile(spike_trains[1], spike_trains[2]) - - avrg_profile = spk.average_profile([isi_profile1, isi_profile2, isi_profile3]) - x, y = avrg_profile.get_plottable_data() - plt.plot(x, y, label="Average ISI profile") - -Note the difference between the :code:`average_profile` function, which returns a :code:`PieceWiseConstFunc` (or :code:`PieceWiseLinFunc`, see below), and the :code:`avrg` member function above, that computes the integral over the time profile resulting in a single value. -So to obtain overall average ISI-distance of a list of ISI profiles you can first compute the average profile using :code:`average_profile` and the use - -.. code:: python - - avrg_isi = avrg_profile.avrg() - -to obtain the final, scalar average ISI distance of the whole set (see also "Computing multivariate distance" below). - SPIKE-distance .............. -To compute for the spike distance you use the function :code:`spike_profile` instead of :code:`isi_profile` above. +To compute for the spike distance profile you use the function :code:`spike_profile` instead of :code:`isi_profile` above. But the general approach is very similar: .. code:: python @@ -217,22 +198,52 @@ Again, you can use to compute the SPIKE distance directly, if you are not interested in the profile at all. The parameter :code:`interval` is optional and if neglected the whole spike train is used. -Furthmore, you can use the :code:`average_profile` function to compute an average profile of a list of SPIKE-profiles: + + +SPIKE synchronization +.............. + +**Important note:** + +------------------------------ + + SPIKE-Synchronization measures *similarity*. + That means, a value of zero indicates absence of synchrony, while a value of one denotes the presence of synchrony. + This is exactly opposite to the other two measures: ISI- and SPIKE-distance. + +---------------------- + + +SPIKE synchronization is another approach to measure spike synchrony. +In contrast to the SPIKE- and ISI-distance, it measures similarity instead of dissimilarity, i.e. higher values represent larger synchrony. +Another difference is that the SPIKE synchronization profile is only defined exactly at the spike times, not for the whole interval of the spike trains. +Therefore, it is represented by a :code:`DiscreteFunction`. + +To compute for the spike synchronization profile, PySpike provides the function :code:`spike_sync_profile`. +The general handling of the profile, however, is similar to the other profiles above: .. code:: python + + import matplotlib.pyplot as plt + import pyspike as spk - avrg_profile = spk.average_profile([spike_profile1, spike_profile2, - spike_profile3]) - x, y = avrg_profile.get_plottable_data() - plt.plot(x, y, label="Average SPIKE profile") + spike_trains = spk.load_spike_trains_from_txt("PySpike_testdata.txt", + time_interval=(0, 4000)) + spike_profile = spk.spike_sync_profile(spike_trains[0], spike_trains[1]) + x, y = spike_profile.get_plottable_data() + +For the direct computation of the overall spike synchronization value within some interval, the :code:`spike_sync` function can be used: + +.. code:: python + + spike_sync = spk.spike_sync(spike_trains[0], spike_trains[1], interval) Computing multivariate profiles and distances ---------------------------------------------- -To compute the multivariate ISI- or SPIKE-profile of a set of spike trains, you can compute all bivariate profiles separately and then use the :code:`average_profile` function above. -However, PySpike provides convenience functions for that purpose. -The following example computes the multivariate ISI- and SPIKE-profile for a list of spike trains: +To compute the multivariate ISI-profile, SPIKE-profile or SPIKE-Synchronization profile f a set of spike trains, PySpike provides multi-variate version of the profile function. +The following example computes the multivariate ISI-, SPIKE- and SPIKE-Sync-profile for a list of spike trains: .. code:: python @@ -240,15 +251,16 @@ The following example computes the multivariate ISI- and SPIKE-profile for a lis time_interval=(0, 4000)) avrg_isi_profile = spk.isi_profile_multi(spike_trains) avrg_spike_profile = spk.spike_profile_multi(spike_trains) + avrg_spike_sync_profile = spk.spike_sync_profile_multi(spike_trains) -Both functions take an optional parameter :code:`indices`, a list of indices that allows to define the spike trains that should be used for the multivariate profile. -As before, if you are only interested in the distance values, and not in the profile, PySpike offers the functions: :code:`isi_distance_multi` and :code:`spike_distance_multi`, that return the scalar multivariate ISI- and SPIKE-distance. -Both distance functions also accept an :code:`interval` parameter that can be used to specify the begin and end of the averaging interval as a pair of floats, if neglected the complete interval is used. +All functions take an optional parameter :code:`indices`, a list of indices that allows to define the spike trains that should be used for the multivariate profile. +As before, if you are only interested in the distance values, and not in the profile, PySpike offers the functions: :code:`isi_distance_multi`, :code:`spike_distance_multi` and :code:`spike_sync_multi`, that return the scalar overall multivariate ISI- and SPIKE-distance as well as the SPIKE-Synchronization value. +Those functions also accept an :code:`interval` parameter that can be used to specify the begin and end of the averaging interval as a pair of floats, if neglected the complete interval is used. Another option to characterize large sets of spike trains are distance matrices. -Each entry in the distance matrix represents a bivariate distance of two spike trains. -The distance matrix is symmetric and has zero values at the diagonal. -The following example computes and plots the ISI- and SPIKE-distance matrix, where for the latter one the averaging is performed only the time interval T=0..1000. +Each entry in the distance matrix represents a bivariate distance (similarity for SPIKE-Synchronization) of two spike trains. +The distance matrix is symmetric and has zero values (ones) at the diagonal. +The following example computes and plots the ISI- and SPIKE-distance matrix as well as the SPIKE-Synchronization-matrix, with different intervals. .. code:: python @@ -264,6 +276,11 @@ The following example computes and plots the ISI- and SPIKE-distance matrix, whe plt.imshow(spike_distance, interpolation='none') plt.title("SPIKE-distance") + plt.figure() + spike_sync = spk.spike_sync_matrix(spike_trains, interval=(2000,4000)) + plt.imshow(spike_sync, interpolation='none') + plt.title("SPIKE-Sync") + plt.show() @@ -284,6 +301,7 @@ Curie Initial Training Network* `Neural Engineering Transformative Technologies .. _ISI: http://www.scholarpedia.org/article/Measures_of_spike_train_synchrony#ISI-distance .. _SPIKE: http://www.scholarpedia.org/article/SPIKE-distance +.. _SPIKE-Synchronization: http://www.scholarpedia.org/article/Measures_of_spike_train_synchrony#SPIKE_synchronization .. _cython: http://www.cython.org .. _SPIKY: http://wwwold.fi.isc.cnr.it/users/thomas.kreuz/Source-Code/SPIKY.html .. _BSD_License: http://opensource.org/licenses/BSD-2-Clause diff --git a/examples/distance_matrix.py b/examples/distance_matrix.py index 142db2c..0e4d825 100644 --- a/examples/distance_matrix.py +++ b/examples/distance_matrix.py @@ -30,4 +30,9 @@ spike_distance = spk.spike_distance_matrix(spike_trains, interval=(0, 1000)) plt.imshow(spike_distance, interpolation='none') plt.title("SPIKE-distance, T=0-1000") +plt.figure() +spike_sync = spk.spike_sync_matrix(spike_trains, interval=(2000, 4000)) +plt.imshow(spike_sync, interpolation='none') +plt.title("SPIKE-Sync, T=2000-4000") + plt.show() diff --git a/pyspike/distances.py b/pyspike/distances.py index 8a14a8d..2cac4bc 100644 --- a/pyspike/distances.py +++ b/pyspike/distances.py @@ -133,6 +133,19 @@ def spike_distance(spikes1, spikes2, interval=None): # spike_sync_profile ############################################################ def spike_sync_profile(spikes1, spikes2): + """ Computes the spike-synchronization profile S_sync(t) of the two given + spike trains. Returns the profile as a DiscreteFunction object. The S_sync + values are either 1 or 0, indicating the presence or absence of a + coincidence. The spike trains are expected to have auxiliary spikes at the + beginning and end of the interval. Use the function add_auxiliary_spikes to + add those spikes to the spike train. + + :param spikes1: ordered array of spike times with auxiliary spikes. + :param spikes2: ordered array of spike times with auxiliary spikes. + :returns: The spike-distance profile :math:`S_{sync}(t)`. + :rtype: :class:`pyspike.function.DiscreteFunction` + + """ # cython implementation try: @@ -155,6 +168,20 @@ Falling back to slow python backend.") # spike_sync ############################################################ def spike_sync(spikes1, spikes2, interval=None): + """ Computes the spike synchronization value SYNC of the given spike + trains. The spike synchronization value is the computed as the total number + of coincidences divided by the total number of spikes: + + .. math:: SYNC = \sum_n C_n / N. + + :param spikes1: ordered array of spike times with auxiliary spikes. + :param spikes2: ordered array of spike times with auxiliary spikes. + :param interval: averaging interval given as a pair of floats (T0, T1), + if None the average over the whole function is computed. + :type interval: Pair of floats or None. + :returns: The spike synchronization value. + :rtype: double + """ return spike_sync_profile(spikes1, spikes2).avrg(interval) @@ -340,16 +367,16 @@ def spike_distance_multi(spike_trains, indices=None, interval=None): ############################################################ def spike_sync_profile_multi(spike_trains, indices=None): """ Computes the multi-variate spike synchronization profile for a set of - spike trains. That is the average spike-distance of all pairs of spike - trains: - :math:`S_ss(t) = 2/((N(N-1)) sum_{} S_{ss}^{i, j}`, - where the sum goes over all pairs + spike trains. For each spike in the set of spike trains, the multi-variate + profile is defined as the number of coincidences divided by the number of + spike trains pairs involving the spike train of containing this spike, + which is the number of spike trains minus one (N-1). :param spike_trains: list of spike trains :param indices: list of indices defining which spike trains to use, if None all given spike trains are used (default=None) :type indices: list or None - :returns: The averaged spike profile :math:`(t)` + :returns: The multi-variate spike sync profile :math:`(t)` :rtype: :class:`pyspike.function.DiscreteFunction` """ @@ -374,7 +401,7 @@ def spike_sync_multi(spike_trains, indices=None, interval=None): :param interval: averaging interval given as a pair of floats, if None the average over the whole function is computed. :type interval: Pair of floats or None. - :returns: The averaged spike synchronization value SYNC. + :returns: The multi-variate spike synchronization value SYNC. :rtype: double """ return spike_sync_profile_multi(spike_trains, indices).avrg(interval) @@ -459,7 +486,7 @@ def spike_distance_matrix(spike_trains, indices=None, interval=None): # spike_sync_matrix ############################################################ def spike_sync_matrix(spike_trains, indices=None, interval=None): - """ Computes the time averaged spike-synchronization value of all pairs of + """ Computes the overall spike-synchronization value of all pairs of spike-trains. :param spike_trains: list of spike trains -- cgit v1.2.3