From f9529c78538882879a07cb67e342eade8d2153ab Mon Sep 17 00:00:00 2001 From: Mario Mulansky Date: Mon, 15 Sep 2014 17:01:13 +0200 Subject: isi distance and basic example --- examples/SPIKY_testdata.txt | 40 ++++++++++++++++++++++++ examples/test_data.py | 28 +++++++++++++++++ pyspike/__init__.py | 5 +++ pyspike/distances.py | 75 +++++++++++++++++++++++++++++++++++++++++++++ pyspike/function.py | 60 ++++++++++++++++++++++++++++++++++++ pyspike/spikes.py | 18 +++++++++++ 6 files changed, 226 insertions(+) create mode 100755 examples/SPIKY_testdata.txt create mode 100644 examples/test_data.py create mode 100644 pyspike/__init__.py create mode 100644 pyspike/distances.py create mode 100644 pyspike/function.py create mode 100644 pyspike/spikes.py diff --git a/examples/SPIKY_testdata.txt b/examples/SPIKY_testdata.txt new file mode 100755 index 0000000..8fa3fcf --- /dev/null +++ b/examples/SPIKY_testdata.txt @@ -0,0 +1,40 @@ +64.886 305.81 696 937.77 1059.7 1322.2 1576.1 1808.1 2121.5 2381.1 2728.6 2966.9 3223.7 3473.7 3644.3 3936.3 +65.553 307.49 696.63 948.66 1070.4 1312.2 1712.7 1934.3 2117.6 2356.9 2727.3 2980.6 3226.9 3475.7 3726.4 3944 +69.064 319.1 688.32 947.85 1071.8 1300.8 1697.2 1930.6 2139.4 2354.2 2723.7 2963.6 3221.3 3470.1 +59.955 313.83 692.23 955.95 1070.4 1319.6 1681.9 1963.5 2151.4 2373.8 2729.4 2971.2 3220.2 3475.5 3632.3 3788.9 +59.977 306.84 686.09 935.08 1059.9 1325.9 1543.4 1821.9 2150.2 2390.4 2724.5 2969.6 3222.5 3471.5 3576 3913.9 +66.415 313.41 688.83 931.43 1051.8 1304.6 1555.6 1820.2 2150.5 2383.1 2723.4 2947.7 3196.6 3443.5 3575 3804.9 +66.449 311.02 689.26 947.12 1058.9 1286.6 1708.2 1957.3 2124.8 2375.7 2709.4 2977.6 3191.1 3449.6 3590.4 3831.2 +63.764 318.45 697.48 936.97 1059.3 1325 1687.9 1944.7 2132.5 2377.1 2713.1 2976.6 3196.8 3442.6 3741.6 3998.3 +63.906 314.79 693.26 937.12 1065.9 1315.8 1584.3 1821.5 2126.3 2396.8 2709.1 2967 3197.4 3444 3732.8 3849.5 +69.493 316.62 689.81 943.62 1071.9 1296.3 1654.8 1931.9 2127.5 2390.6 2708.9 2950.4 3194.8 3445.2 3670.1 3903.3 +61.789 317.53 555.82 813.15 1198.7 1448.7 1686.7 1943.5 2060.7 2311.4 2658.2 2900.2 3167.4 3418.2 3617.3 3771 +64.098 309.86 567.27 813.91 1182 1464.3 1576.8 1822.5 2063.1 2311.7 2655.8 2911.7 3168.3 3418.2 3586.4 3999.7 +68.59 315.5 559.52 806.23 1182.5 1441.1 1567.2 1804.8 2074.9 2315.8 2655.1 2913.2 3165.9 3419.5 3648.1 3884.4 +66.507 314.42 556.42 814.83 1182.5 1440.3 1701.3 1911.1 2069.7 2319.3 2662.3 2903.2 3167.4 3418.5 3545 3893.9 +72.744 318.45 554.4 819.64 1186.9 1449.7 1676 1957.4 2051.4 2302.8 2657.8 2916.2 3169.4 3416.7 3570.4 3884.8 +64.779 324.42 560.56 828.99 1174.8 1439.9 1563.7 1790.6 2067.7 2287.6 2657.4 2905.2 3139.2 3389.1 3507.8 3807.5 +64.852 316.63 568.89 815.61 1198.3 1454.1 1710.6 1933.9 2091.5 2309.6 2660.9 2907.5 3137.2 3389.3 3617.2 +63.089 314.52 553.8 827.06 1183.9 1457.6 1558.9 1808.3 2064.5 2337.1 2653.6 2897 3143.7 3385.7 3668.7 3803.8 +62.23 315.16 564.35 812.15 1199.6 1448.9 1562.7 1839.1 2069.7 2308.9 2649.6 2919.7 3141 3389.9 3723.6 3882.2 +69.662 311.93 564.91 805.25 1209.7 1451.4 1691.9 1932.1 2044.2 2329.4 2657.1 2908.5 3142.8 3390.5 3597.3 3991.1 +183.42 431.34 562.41 809.57 1086.3 1308.9 1555.9 1831.3 2057 2326.9 2591.3 2831.4 3113.9 3367.9 3555.3 3956 +188.49 442.39 572.4 810.76 1065 1326.7 1564.3 1803.4 2060.4 2322.4 2607.2 2824.1 3110.2 3363.9 3644.1 3819.6 +177 437.76 569.82 819.66 1064.1 1309.2 1685.7 1957.5 2066.9 2313.8 2593.2 2847 3116.8 3364.5 3727.3 3881.6 +193.9 441.93 586.9 804.98 1062.5 1312.4 1542.4 1793.1 2073.9 2314.7 2587.8 2845.9 3112.4 3359.8 +193.01 440.26 555.64 814.08 1056.3 1315 1689.9 1961.4 2049.1 2305 2593.9 2847.5 3110.6 3361.1 3711.6 3914.7 +194.71 437.57 566.18 806.73 1069.2 1314.6 1682.7 1942.2 2061.8 2304.6 2607.6 2841.7 3082.9 3330.3 3679.7 3848.2 +184.88 441.22 570.92 794.35 1063.7 1309.9 1678.7 1930 2058 2321.3 2606.7 2845 3084.8 3337.3 3640 3952.1 +189.66 443.59 560.67 816.89 1070.4 1303.4 1550.1 1815.5 2057.6 2323.7 2587.1 2843.5 3086.6 3333.6 3618.2 3815.4 +190.41 440.77 568.96 808.56 1073.8 1322.1 1686.5 1952.8 2068.7 2335.7 2595.7 2845.4 3086 3333.5 3635.6 3939.3 +181.16 440.67 577.54 823.52 1052.5 1322.3 1578.4 1822.2 2079.4 2309.1 2596.9 2851.9 3083.5 3335.1 3531.2 3770.6 +181.09 434.97 687.15 943.33 1192.9 1444 1699.4 1942 2194.6 2445.9 2549.4 2785.1 3056.5 3308.2 3620.5 3932.7 +186.7 446.53 688.18 942.86 1186.1 1441.9 1688.1 1922.2 2196.6 2455.3 2534.8 2776.5 3060.3 3309.4 3514.1 3808.6 +196.76 446 681.26 948.27 1195.8 1433.1 1699 1933 2201.2 2461.4 2547.4 2777.8 3055.7 3307.1 3590.6 3952.8 +200.68 427.11 695.67 946.42 1178.6 1440.1 1538.4 1809 2199.8 2432.5 2531.6 2793.2 3056.6 3308.6 3510.6 3928.1 +190.83 429.57 698.73 931.16 1190.6 1428.9 1698.3 1935 2176.8 2424.7 2530.5 2766.9 3062 3309.7 3689.8 +181.47 441.93 682.32 943.01 1190.1 1459.1 1570.6 1819.6 2189.8 2437.9 2543.3 2782.8 3025.9 3280.2 3581 3855.9 +191.38 435.69 702.76 935.62 1188.3 1438.3 1564.2 1823.9 2191.3 2444.9 2531.9 2782.4 3030.7 3275.7 3677.7 3829.2 +191.97 433.85 686.29 932.65 1183.1 1432.7 1563.9 1826.5 2214.1 2436.8 2529.8 2778.3 3028.3 3281.8 3582 3863.4 +189.51 453.21 691.3 940.86 1180.1 1430.1 1567.1 1835 2199 2448.2 2526.7 2773.8 3030.5 3280.1 3576.2 3893.6 +190.88 435.48 692.66 940.51 1189.5 1448.9 1575.1 1824.2 2190.8 2425.9 2530.6 2783.3 3033.3 3279.5 3733 3838.9 diff --git a/examples/test_data.py b/examples/test_data.py new file mode 100644 index 0000000..bddcb15 --- /dev/null +++ b/examples/test_data.py @@ -0,0 +1,28 @@ +# compute the isi distance of some test data +from __future__ import print_function + +import numpy as np +import matplotlib.pyplot as plt + +import pyspike as spk + +# first load the data +spike_trains = [] +spike_file = open("SPIKY_testdata.txt", 'r') +for line in spike_file: + spike_trains.append(spk.spike_train_from_string(line)) + +# plot the spike time +for (i,spikes) in enumerate(spike_trains): + plt.plot(spikes, i*np.ones_like(spikes), 'o') + +f = spk.isi_distance(spike_trains[0], spike_trains[1], 4000) +x, y = f.get_plottable_data() + +plt.figure() +plt.plot(x, y, '-k') + +print("Average: %.8f" % f.avrg()) +print("Absolute average: %.8f" % f.abs_avrg()) + +plt.show() diff --git a/pyspike/__init__.py b/pyspike/__init__.py new file mode 100644 index 0000000..6651eb5 --- /dev/null +++ b/pyspike/__init__.py @@ -0,0 +1,5 @@ +__all__ = ["function", "distances", "spikes"] + +from function import PieceWiseConstFunc +from distances import isi_distance +from spikes import spike_train_from_string diff --git a/pyspike/distances.py b/pyspike/distances.py new file mode 100644 index 0000000..d9790dc --- /dev/null +++ b/pyspike/distances.py @@ -0,0 +1,75 @@ +""" distances.py + +Module containing several function to compute spike distances + +Copyright 2014, Mario Mulansky +""" + +import numpy as np + +from pyspike import PieceWiseConstFunc + +def spike_train_from_string(s, sep=' '): + """ Converts a string of times into a SpikeTrain object. + Params: + - s: the string with (ordered) spike times + - sep: The separator between the time numbers. + Returns: + - array of spike times + """ + return np.fromstring(s, sep=sep) + +def isi_distance(spikes1, spikes2, T_end, T_start=0.0): + """ Computes the instantaneous isi-distance S_isi (t) of the two given spike + trains. + Args: + - spikes1, spikes2: ordered arrays of spike times. + - T_end: end time of the observation interval. + - T_start: begin of the observation interval (default=0.0). + Returns: + - PieceWiseConstFunc describing the isi-distance. + """ + # add spikes at the beginning and end of the interval + s1 = np.empty(len(spikes1)+2) + s1[0] = T_start + s1[-1] = T_end + s1[1:-1] = spikes1 + s2 = np.empty(len(spikes2)+2) + s2[0] = T_start + s2[-1] = T_end + s2[1:-1] = spikes2 + + # compute the interspike interval + nu1 = s1[1:]-s1[:-1] + nu2 = s2[1:]-s2[:-1] + + # compute the isi-distance + spike_events = np.empty(len(nu1)+len(nu2)) + spike_events[0] = T_start + spike_events[-1] = T_end + # the values have one entry less - the number of intervals between events + isi_values = np.empty(len(spike_events)-1) + # add the distance of the first events + # isi_values[0] = nu1[0]/nu2[0] - 1.0 if nu1[0] <= nu2[0] \ + # else 1.0 - nu2[0]/nu1[0] + isi_values[0] = (nu1[0]-nu2[0])/max(nu1[0],nu2[0]) + index1 = 0 + index2 = 0 + index = 1 + while True: + # check which spike is next - from s1 or s2 + if s1[index1+1] <= s2[index2+1]: + index1 += 1 + if index1 >= len(nu1): + break + spike_events[index] = s1[index1] + else: + index2 += 1 + if index2 >= len(nu2): + break + spike_events[index] = s2[index2] + # compute the corresponding isi-distance + isi_values[index] = (nu1[index1]-nu2[index2]) / \ + max(nu1[index1], nu2[index2]) + index += 1 + return PieceWiseConstFunc(spike_events, isi_values) diff --git a/pyspike/function.py b/pyspike/function.py new file mode 100644 index 0000000..c1de9cb --- /dev/null +++ b/pyspike/function.py @@ -0,0 +1,60 @@ +""" function.py + +Module containing classes representing piece-wise constant and piece-wise linear +functions. + +Copyright 2014, Mario Mulansky + +""" +from __future__ import print_function + +import numpy as np + +class PieceWiseConstFunc: + """ A class representing a piece-wise constant function. """ + + def __init__(self, x, y): + """ Constructs the piece-wise const function. + Params: + - x: array of length N+1 defining the edges of the intervals of the pwc + function. + - y: array of length N defining the function values at the intervals. + """ + self.x = x + self.y = y + + def get_plottable_data(self): + """ Returns two arrays containing x- and y-coordinates for immeditate + plotting of the piece-wise function. + """ + + x_plot = np.empty(2*len(self.x)-2) + x_plot[0] = self.x[0] + x_plot[1::2] = self.x[1:] + x_plot[2::2] = self.x[1:-1] + y_plot = np.empty(2*len(self.y)) + y_plot[::2] = self.y + y_plot[1::2] = self.y + + return x_plot, y_plot + + def avrg(self): + """ Computes the average of the piece-wise const function: + a = 1/T int f(x) dx where T is the length of the interval. + Returns: + - the average a. + """ + print(self.x) + print(self.y) + return np.sum((self.x[1:]-self.x[:-1]) * self.y) / \ + (self.x[-1]-self.x[0]) + + def abs_avrg(self): + """ Computes the average of the abs value of the piece-wise const + function: + a = 1/T int |f(x)| dx where T is the length of the interval. + Returns: + - the average a. + """ + return np.sum((self.x[1:]-self.x[:-1]) * np.abs(self.y)) / \ + (self.x[-1]-self.x[0]) diff --git a/pyspike/spikes.py b/pyspike/spikes.py new file mode 100644 index 0000000..42b6501 --- /dev/null +++ b/pyspike/spikes.py @@ -0,0 +1,18 @@ +""" spikes.py + +Module containing several function to load and transform spike trains + +Copyright 2014, Mario Mulansky +""" + +import numpy as np + +def spike_train_from_string(s, sep=' '): + """ Converts a string of times into a SpikeTrain object. + Params: + - s: the string with (ordered) spike times + - sep: The separator between the time numbers. + Returns: + - array of spike times + """ + return np.fromstring(s, sep=sep) -- cgit v1.2.3