From ac8e9ca85a889ee2d9fe984d19976917ffdfea46 Mon Sep 17 00:00:00 2001 From: Mario Mulansky Date: Thu, 23 Oct 2014 12:48:19 +0200 Subject: +averaging example, removed empty averaging sec from docs --- examples/averages.py | 42 ++++++++++++++++++++++++++++++++++++++++++ 1 file changed, 42 insertions(+) create mode 100644 examples/averages.py (limited to 'examples/averages.py') diff --git a/examples/averages.py b/examples/averages.py new file mode 100644 index 0000000..a06871a --- /dev/null +++ b/examples/averages.py @@ -0,0 +1,42 @@ +""" averages.py + +Simple example showing how to compute averages of distance profiles + +Copyright 2014, Mario Mulansky + +Distributed under the BSD License +""" + +from __future__ import print_function + +import pyspike as spk + +spike_trains = spk.load_spike_trains_from_txt("PySpike_testdata.txt", + time_interval=(0, 4000)) + +f = spk.isi_profile(spike_trains[0], spike_trains[1]) + +print("ISI-distance: %.8f" % f.avrg()) + +isi1 = f.avrg(interval=(0, 1000)) +isi2 = f.avrg(interval=(1000, 2000)) +isi3 = f.avrg(interval=(2000, 3000)) +isi4 = f.avrg(interval=(3000, 4000)) + +print("ISI-distance (0-1000): %.8f" % isi1) +print("ISI-distance (1000-2000): %.8f" % isi2) +print("ISI-distance (2000-3000): %.8f" % isi3) +print("ISI-distance (3000-4000): %.8f" % isi4) +print() + +print("SPIKE-distance: %.8f" % f.avrg()) + +spike1 = f.avrg(interval=(0, 1000)) +spike2 = f.avrg(interval=(1000, 2000)) +spike3 = f.avrg(interval=(2000, 3000)) +spike4 = f.avrg(interval=(3000, 4000)) + +print("SPIKE-distance (0-1000): %.8f" % spike1) +print("SPIKE-distance (1000-2000): %.8f" % spike2) +print("SPIKE-distance (2000-3000): %.8f" % spike3) +print("SPIKE-distance (3000-4000): %.8f" % spike4) -- cgit v1.2.3 From 3a36f81d52137435910a4b7c656a478ae5b38ece Mon Sep 17 00:00:00 2001 From: Mario Mulansky Date: Thu, 23 Oct 2014 15:33:06 +0200 Subject: support for multiple averaging intervals --- examples/averages.py | 26 ++++----- pyspike/function.py | 144 ++++++++++++++++++++++++++++++++++++-------------- test/test_function.py | 10 +++- 3 files changed, 126 insertions(+), 54 deletions(-) (limited to 'examples/averages.py') diff --git a/examples/averages.py b/examples/averages.py index a06871a..6413ab2 100644 --- a/examples/averages.py +++ b/examples/averages.py @@ -20,23 +20,25 @@ print("ISI-distance: %.8f" % f.avrg()) isi1 = f.avrg(interval=(0, 1000)) isi2 = f.avrg(interval=(1000, 2000)) -isi3 = f.avrg(interval=(2000, 3000)) -isi4 = f.avrg(interval=(3000, 4000)) +isi3 = f.avrg(interval=[(0, 1000), (2000, 3000)]) +isi4 = f.avrg(interval=[(1000, 2000), (3000, 4000)]) -print("ISI-distance (0-1000): %.8f" % isi1) -print("ISI-distance (1000-2000): %.8f" % isi2) -print("ISI-distance (2000-3000): %.8f" % isi3) -print("ISI-distance (3000-4000): %.8f" % isi4) +print("ISI-distance (0-1000): %.8f" % isi1) +print("ISI-distance (1000-2000): %.8f" % isi2) +print("ISI-distance (0-1000) and (2000-3000): %.8f" % isi3) +print("ISI-distance (1000-2000) and (3000-4000): %.8f" % isi4) print() +f = spk.spike_profile(spike_trains[0], spike_trains[1]) + print("SPIKE-distance: %.8f" % f.avrg()) spike1 = f.avrg(interval=(0, 1000)) spike2 = f.avrg(interval=(1000, 2000)) -spike3 = f.avrg(interval=(2000, 3000)) -spike4 = f.avrg(interval=(3000, 4000)) +spike3 = f.avrg(interval=[(0, 1000), (2000, 3000)]) +spike4 = f.avrg(interval=[(1000, 2000), (3000, 4000)]) -print("SPIKE-distance (0-1000): %.8f" % spike1) -print("SPIKE-distance (1000-2000): %.8f" % spike2) -print("SPIKE-distance (2000-3000): %.8f" % spike3) -print("SPIKE-distance (3000-4000): %.8f" % spike4) +print("SPIKE-distance (0-1000): %.8f" % spike1) +print("SPIKE-distance (1000-2000): %.8f" % spike2) +print("SPIKE-distance (0-1000) and (2000-3000): %.8f" % spike3) +print("SPIKE-distance (1000-2000) and (3000-4000): %.8f" % spike4) diff --git a/pyspike/function.py b/pyspike/function.py index d53042f..0c0a391 100644 --- a/pyspike/function.py +++ b/pyspike/function.py @@ -11,6 +11,7 @@ Distributed under the BSD License from __future__ import print_function import numpy as np +import collections ############################################################## @@ -74,21 +75,18 @@ class PieceWiseConstFunc(object): return x_plot, y_plot - def avrg(self, interval=None): - """ Computes the average of the piece-wise const function: - :math:`a = 1/T int_0^T f(x) dx` where T is the length of the interval. + def integral(self, interval=None): + """ Returns the integral over the given interval. - :param interval: averaging interval given as a pair of floats, if None - the average over the whole function is computed. + :param interval: integration interval given as a pair of floats, if + None the integral over the whole function is computed. :type interval: Pair of floats or None. - :returns: the average a. - :rtype: double - + :returns: the integral + :rtype: float """ if interval is None: - # no interval given, average over the whole spike train - a = np.sum((self.x[1:]-self.x[:-1]) * self.y) / \ - (self.x[-1]-self.x[0]) + # no interval given, integrate over the whole spike train + a = np.sum((self.x[1:]-self.x[:-1]) * self.y) else: # find the indices corresponding to the interval start_ind = np.searchsorted(self.x, interval[0], side='right') @@ -103,7 +101,39 @@ class PieceWiseConstFunc(object): a += (self.x[start_ind]-interval[0]) * self.y[start_ind-1] # correction from last index to end a += (interval[1]-self.x[end_ind]) * self.y[end_ind] - a /= (interval[1]-interval[0]) + return a + + def avrg(self, interval=None): + """ Computes the average of the piece-wise const function: + :math:`a = 1/T int_0^T f(x) dx` where T is the length of the interval. + + :param interval: averaging interval given as a pair of floats, a + sequence of pairs for averaging multiple intervals, or + None, if None the average over the whole function is + computed. + :type interval: Pair, sequence of pairs, or None. + :returns: the average a. + :rtype: float + """ + if interval is None: + # no interval given, average over the whole spike train + return self.integral() / (self.x[-1]-self.x[0]) + + # check if interval is as sequence + assert isinstance(interval, collections.Sequence), \ + "Invalid value for `interval`. None, Sequence or Tuple expected." + # check if interval is a sequence of intervals + if not isinstance(interval[0], collections.Sequence): + # just one interval + a = self.integral(interval) / (interval[1]-interval[0]) + else: + # several intervals + a = 0.0 + int_length = 0.0 + for ival in interval: + a += self.integral(ival) + int_length += ival[1] - ival[0] + a /= int_length return a def add(self, f): @@ -203,16 +233,14 @@ class PieceWiseLinFunc: y_plot[1::2] = self.y2 return x_plot, y_plot - def avrg(self, interval=None): - """ Computes the average of the piece-wise linear function: - :math:`a = 1/T int_0^T f(x) dx` where T is the length of the interval. + def integral(self, interval=None): + """ Returns the integral over the given interval. - :param interval: averaging interval given as a pair of floats, if None - the average over the whole function is computed. + :param interval: integration interval given as a pair of floats, if + None the integral over the whole function is computed. :type interval: Pair of floats or None. - :returns: the average a. - :rtype: double - + :returns: the integral + :rtype: float """ def intermediate_value(x0, x1, y0, y1, x): @@ -220,9 +248,8 @@ class PieceWiseLinFunc: return y0 + (y1-y0)*(x-x0)/(x1-x0) if interval is None: - # no interval given, average over the whole spike train - a = np.sum((self.x[1:]-self.x[:-1]) * 0.5*(self.y1+self.y2)) / \ - (self.x[-1]-self.x[0]) + # no interval given, integrate over the whole spike train + integral = np.sum((self.x[1:]-self.x[:-1]) * 0.5*(self.y1+self.y2)) else: # find the indices corresponding to the interval start_ind = np.searchsorted(self.x, interval[0], side='right') @@ -230,26 +257,61 @@ class PieceWiseLinFunc: assert start_ind > 0 and end_ind < len(self.x), \ "Invalid averaging interval" # first the contribution from between the indices - a = np.sum((self.x[start_ind+1:end_ind+1] - - self.x[start_ind:end_ind]) * - 0.5*(self.y1[start_ind:end_ind] + - self.y2[start_ind:end_ind])) + integral = np.sum((self.x[start_ind+1:end_ind+1] - + self.x[start_ind:end_ind]) * + 0.5*(self.y1[start_ind:end_ind] + + self.y2[start_ind:end_ind])) # correction from start to first index - a += (self.x[start_ind]-interval[0]) * 0.5 * \ - (self.y2[start_ind-1] + - intermediate_value(self.x[start_ind-1], self.x[start_ind], - self.y1[start_ind-1], - self.y2[start_ind-1], - interval[0] - )) + integral += (self.x[start_ind]-interval[0]) * 0.5 * \ + (self.y2[start_ind-1] + + intermediate_value(self.x[start_ind-1], + self.x[start_ind], + self.y1[start_ind-1], + self.y2[start_ind-1], + interval[0] + )) # correction from last index to end - a += (interval[1]-self.x[end_ind]) * 0.5 * \ - (self.y1[end_ind] + - intermediate_value(self.x[end_ind], self.x[end_ind+1], - self.y1[end_ind], self.y2[end_ind], - interval[1] - )) - a /= (interval[1]-interval[0]) + integral += (interval[1]-self.x[end_ind]) * 0.5 * \ + (self.y1[end_ind] + + intermediate_value(self.x[end_ind], self.x[end_ind+1], + self.y1[end_ind], self.y2[end_ind], + interval[1] + )) + return integral + + def avrg(self, interval=None): + """ Computes the average of the piece-wise linear function: + :math:`a = 1/T int_0^T f(x) dx` where T is the length of the interval. + + :param interval: averaging interval given as a pair of floats, a + sequence of pairs for averaging multiple intervals, or + None, if None the average over the whole function is + computed. + :type interval: Pair, sequence of pairs, or None. + :returns: the average a. + :rtype: float + + """ + + if interval is None: + # no interval given, average over the whole spike train + return self.integral() / (self.x[-1]-self.x[0]) + + # check if interval is as sequence + assert isinstance(interval, collections.Sequence), \ + "Invalid value for `interval`. None, Sequence or Tuple expected." + # check if interval is a sequence of intervals + if not isinstance(interval[0], collections.Sequence): + # just one interval + a = self.integral(interval) / (interval[1]-interval[0]) + else: + # several intervals + a = 0.0 + int_length = 0.0 + for ival in interval: + a += self.integral(ival) + int_length += ival[1] - ival[0] + a /= int_length return a def add(self, f): diff --git a/test/test_function.py b/test/test_function.py index ed0b2ed..ba87db8 100644 --- a/test/test_function.py +++ b/test/test_function.py @@ -39,6 +39,10 @@ def test_pwc(): a = f.avrg([1.0, 4.0]) assert_almost_equal(a, (-0.5*1.0+0.5*1.5+1.5*0.75)/3.0, decimal=16) + # averaging over multiple intervals + a = f.avrg([(0.5, 1.5), (1.5, 3.5)]) + assert_almost_equal(a, (0.5-0.5+0.5*1.5+1.0*0.75)/3.0, decimal=16) + def test_pwc_add(): # some random data @@ -116,12 +120,16 @@ def test_pwl(): a = f.avrg([1.5, 3.5]) assert_almost_equal(a, (-0.425*0.5 + 0.75 + (0.75+0.75-0.5/1.5)/2) / 2.0, decimal=16) - a = f.avrg([1.0, 3.5]) + a = f.avrg((1.0, 3.5)) assert_almost_equal(a, (-0.45 + 0.75 + (0.75+0.75-0.5/1.5)/2) / 2.5, decimal=16) a = f.avrg([1.0, 4.0]) assert_almost_equal(a, (-0.45 + 0.75 + 1.5*0.5) / 3.0, decimal=16) + # averaging over multiple intervals + a = f.avrg([(0.5, 1.5), (1.5, 2.5)]) + assert_almost_equal(a, (1.375*0.5 - 0.45 + 0.75)/2.0, decimal=16) + def test_pwl_add(): x = [0.0, 1.0, 2.0, 2.5, 4.0] -- cgit v1.2.3 From 7bec7f0fe1c40146e8b45757d4c156a4f9b49001 Mon Sep 17 00:00:00 2001 From: Mario Mulansky Date: Tue, 20 Jan 2015 14:39:22 +0100 Subject: add spike sync to some examples --- examples/averages.py | 15 +++++++++++++++ examples/spike_sync.py | 3 +-- 2 files changed, 16 insertions(+), 2 deletions(-) (limited to 'examples/averages.py') diff --git a/examples/averages.py b/examples/averages.py index 6413ab2..c3e81e2 100644 --- a/examples/averages.py +++ b/examples/averages.py @@ -42,3 +42,18 @@ print("SPIKE-distance (0-1000): %.8f" % spike1) print("SPIKE-distance (1000-2000): %.8f" % spike2) print("SPIKE-distance (0-1000) and (2000-3000): %.8f" % spike3) print("SPIKE-distance (1000-2000) and (3000-4000): %.8f" % spike4) +print() + +f = spk.spike_sync_profile(spike_trains[0], spike_trains[1]) + +print("SPIKE-Synchronization: %.8f" % f.avrg()) + +spike_sync1 = f.avrg(interval=(0, 1000)) +spike_sync2 = f.avrg(interval=(1000, 2000)) +spike_sync3 = f.avrg(interval=[(0, 1000), (2000, 3000)]) +spike_sync4 = f.avrg(interval=[(1000, 2000), (3000, 4000)]) + +print("SPIKE-Sync (0-1000): %.8f" % spike_sync1) +print("SPIKE-Sync (1000-2000): %.8f" % spike_sync2) +print("SPIKE-Sync (0-1000) and (2000-3000): %.8f" % spike_sync3) +print("SPIKE-Sync (1000-2000) and (3000-4000): %.8f" % spike_sync4) diff --git a/examples/spike_sync.py b/examples/spike_sync.py index 535f19f..a72dbbf 100644 --- a/examples/spike_sync.py +++ b/examples/spike_sync.py @@ -11,7 +11,6 @@ spike_trains = spk.load_spike_trains_from_txt("../test/SPIKE_Sync_Test.txt", plt.figure() f = spk.spike_sync_profile(spike_trains[0], spike_trains[1]) -# f = spk.spike_sync_profile(spikes1, spikes2) x, y = f.get_plottable_data() plt.plot(x, y, '--ok', label="SPIKE-SYNC profile") print(f.x) @@ -33,7 +32,7 @@ plt.figure() f = spk.spike_sync_profile_multi(spike_trains) x, y = f.get_plottable_data() -plt.plot(x, y, '-k', label="SPIKE-SYNC profile") +plt.plot(x, y, '-k', label="SPIKE-Sync profile") print("Average:", f.avrg()) -- cgit v1.2.3 From ee0e980b72c299eed12b7a3afc542fc470dd6d98 Mon Sep 17 00:00:00 2001 From: Mario Mulansky Date: Wed, 9 Mar 2016 12:30:35 +0100 Subject: updated examples to use new unified interface removed all occasions of *multi functions from examples as they are considered deprecated now. Uses unified interface everywhere. --- examples/averages.py | 2 +- examples/merge.py | 6 +++--- examples/multivariate.py | 4 ++-- examples/performance.py | 27 +++++++++++++++------------ examples/plot.py | 5 +++-- examples/profiles.py | 4 ++-- examples/spike_sync.py | 2 +- 7 files changed, 27 insertions(+), 23 deletions(-) (limited to 'examples/averages.py') diff --git a/examples/averages.py b/examples/averages.py index c3e81e2..8b405d0 100644 --- a/examples/averages.py +++ b/examples/averages.py @@ -12,7 +12,7 @@ from __future__ import print_function import pyspike as spk spike_trains = spk.load_spike_trains_from_txt("PySpike_testdata.txt", - time_interval=(0, 4000)) + edges=(0, 4000)) f = spk.isi_profile(spike_trains[0], spike_trains[1]) diff --git a/examples/merge.py b/examples/merge.py index 2ea96ea..b4437a3 100644 --- a/examples/merge.py +++ b/examples/merge.py @@ -21,9 +21,9 @@ merged_spike_train = spk.merge_spike_trains([spike_trains[0], spike_trains[1]]) print(merged_spike_train.spikes) -plt.plot(spike_trains[0].spikes, np.ones_like(spike_trains[0].spikes), 'o') -plt.plot(spike_trains[1].spikes, np.ones_like(spike_trains[1].spikes), 'x') +plt.plot(spike_trains[0], np.ones_like(spike_trains[0]), 'o') +plt.plot(spike_trains[1], np.ones_like(spike_trains[1]), 'x') plt.plot(merged_spike_train.spikes, - 2*np.ones_like(merged_spike_train.spikes), 'o') + 2*np.ones_like(merged_spike_train), 'o') plt.show() diff --git a/examples/multivariate.py b/examples/multivariate.py index 93f8516..e9579a5 100644 --- a/examples/multivariate.py +++ b/examples/multivariate.py @@ -28,7 +28,7 @@ num_of_spikes = sum([len(spike_trains[i]) print("Number of spikes: %d" % num_of_spikes) # calculate the multivariate spike distance -f = spk.spike_profile_multi(spike_trains) +f = spk.spike_profile(spike_trains) t_spike = time.clock() @@ -39,7 +39,7 @@ print("Spike distance from average: %.8f" % avrg) t_avrg = time.clock() # compute average distance directly, should give the same result as above -spike_dist = spk.spike_distance_multi(spike_trains) +spike_dist = spk.spike_distance(spike_trains) print("Spike distance directly: %.8f" % spike_dist) t_dist = time.clock() diff --git a/examples/performance.py b/examples/performance.py index ec6c830..30691f8 100644 --- a/examples/performance.py +++ b/examples/performance.py @@ -31,38 +31,41 @@ for i in range(M): t_end = datetime.now() runtime = (t_end-t_start).total_seconds() +sort_by = 'tottime' +# sort_by = 'cumtime' + print("Spike generation runtime: %.3fs" % runtime) print() print("================ ISI COMPUTATIONS ================") print(" MULTIVARIATE DISTANCE") -cProfile.run('spk.isi_distance_multi(spike_trains)', 'performance.stat') +cProfile.run('spk.isi_distance(spike_trains)', 'performance.stat') p = pstats.Stats('performance.stat') -p.strip_dirs().sort_stats('tottime').print_stats(5) +p.strip_dirs().sort_stats(sort_by).print_stats(5) print(" MULTIVARIATE PROFILE") -cProfile.run('spk.isi_profile_multi(spike_trains)', 'performance.stat') +cProfile.run('spk.isi_profile(spike_trains)', 'performance.stat') p = pstats.Stats('performance.stat') -p.strip_dirs().sort_stats('tottime').print_stats(5) +p.strip_dirs().sort_stats(sort_by).print_stats(5) print("================ SPIKE COMPUTATIONS ================") print(" MULTIVARIATE DISTANCE") -cProfile.run('spk.spike_distance_multi(spike_trains)', 'performance.stat') +cProfile.run('spk.spike_distance(spike_trains)', 'performance.stat') p = pstats.Stats('performance.stat') -p.strip_dirs().sort_stats('tottime').print_stats(5) +p.strip_dirs().sort_stats(sort_by).print_stats(5) print(" MULTIVARIATE PROFILE") -cProfile.run('spk.spike_profile_multi(spike_trains)', 'performance.stat') +cProfile.run('spk.spike_profile(spike_trains)', 'performance.stat') p = pstats.Stats('performance.stat') -p.strip_dirs().sort_stats('tottime').print_stats(5) +p.strip_dirs().sort_stats(sort_by).print_stats(5) print("================ SPIKE-SYNC COMPUTATIONS ================") print(" MULTIVARIATE DISTANCE") -cProfile.run('spk.spike_sync_multi(spike_trains)', 'performance.stat') +cProfile.run('spk.spike_sync(spike_trains)', 'performance.stat') p = pstats.Stats('performance.stat') -p.strip_dirs().sort_stats('tottime').print_stats(5) +p.strip_dirs().sort_stats(sort_by).print_stats(5) print(" MULTIVARIATE PROFILE") -cProfile.run('spk.spike_sync_profile_multi(spike_trains)', 'performance.stat') +cProfile.run('spk.spike_sync_profile(spike_trains)', 'performance.stat') p = pstats.Stats('performance.stat') -p.strip_dirs().sort_stats('tottime').print_stats(5) +p.strip_dirs().sort_stats(sort_by).print_stats(5) diff --git a/examples/plot.py b/examples/plot.py index 1922939..a0e04da 100644 --- a/examples/plot.py +++ b/examples/plot.py @@ -24,7 +24,8 @@ spike_trains = spk.load_spike_trains_from_txt("PySpike_testdata.txt", for (i, spike_train) in enumerate(spike_trains): plt.scatter(spike_train, i*np.ones_like(spike_train), marker='|') -f = spk.isi_profile(spike_trains[0], spike_trains[1]) +# profile of the first two spike trains +f = spk.isi_profile(spike_trains, indices=[0, 1]) x, y = f.get_plottable_data() plt.figure() @@ -32,7 +33,7 @@ plt.plot(x, np.abs(y), '--k', label="ISI-profile") print("ISI-distance: %.8f" % f.avrg()) -f = spk.spike_profile(spike_trains[0], spike_trains[1]) +f = spk.spike_profile(spike_trains, indices=[0, 1]) x, y = f.get_plottable_data() plt.plot(x, y, '-b', label="SPIKE-profile") diff --git a/examples/profiles.py b/examples/profiles.py index 05494bd..8412ffb 100644 --- a/examples/profiles.py +++ b/examples/profiles.py @@ -29,7 +29,7 @@ print("Average ISI distance:", f.avrg()) print() # compute the multivariate ISI profile -f = spk.isi_profile_multi(spike_trains) +f = spk.isi_profile(spike_trains) t = 1200 print("Multivariate ISI value at t =", t, ":", f(t)) @@ -56,7 +56,7 @@ print("Average SPIKE distance:", f.avrg()) print() # compute the multivariate SPIKE profile -f = spk.spike_profile_multi(spike_trains) +f = spk.spike_profile(spike_trains) # SPIKE values at certain points t = 1200 diff --git a/examples/spike_sync.py b/examples/spike_sync.py index 37dbff4..13ca0ce 100644 --- a/examples/spike_sync.py +++ b/examples/spike_sync.py @@ -31,7 +31,7 @@ plt.figure() plt.subplot(211) -f = spk.spike_sync_profile_multi(spike_trains) +f = spk.spike_sync_profile(spike_trains) x, y = f.get_plottable_data() plt.plot(x, y, '-b', alpha=0.7, label="SPIKE-Sync profile") -- cgit v1.2.3