1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
|
""" directionality_python_backend.py
Collection of python functions that can be used instead of the cython
implementation.
Copyright 2015, Mario Mulansky <mario.mulansky@gmx.net>
Distributed under the BSD License
"""
import numpy as np
############################################################
# spike_train_order_python
############################################################
def spike_train_order_python(spikes1, spikes2, t_start, t_end, max_tau):
def get_tau(spikes1, spikes2, i, j, max_tau):
m = t_end - t_start # use interval as initial tau
if i < len(spikes1)-1 and i > -1:
m = min(m, spikes1[i+1]-spikes1[i])
if j < len(spikes2)-1 and j > -1:
m = min(m, spikes2[j+1]-spikes2[j])
if i > 0:
m = min(m, spikes1[i]-spikes1[i-1])
if j > 0:
m = min(m, spikes2[j]-spikes2[j-1])
m *= 0.5
if max_tau > 0.0:
m = min(m, max_tau)
return m
N1 = len(spikes1)
N2 = len(spikes2)
i = -1
j = -1
n = 0
st = np.zeros(N1 + N2 + 2) # spike times
a = np.zeros(N1 + N2 + 2) # coincidences
mp = np.ones(N1 + N2 + 2) # multiplicity
while i + j < N1 + N2 - 2:
if (i < N1-1) and (j == N2-1 or spikes1[i+1] < spikes2[j+1]):
i += 1
n += 1
tau = get_tau(spikes1, spikes2, i, j, max_tau)
st[n] = spikes1[i]
if j > -1 and spikes1[i]-spikes2[j] < tau:
# coincidence between the current spike and the previous spike
# both get marked with 1
a[n] = -1
a[n-1] = -1
elif (j < N2-1) and (i == N1-1 or spikes1[i+1] > spikes2[j+1]):
j += 1
n += 1
tau = get_tau(spikes1, spikes2, i, j, max_tau)
st[n] = spikes2[j]
if i > -1 and spikes2[j]-spikes1[i] < tau:
# coincidence between the current spike and the previous spike
# both get marked with 1
a[n] = 1
a[n-1] = 1
else: # spikes1[i+1] = spikes2[j+1]
# advance in both spike trains
j += 1
i += 1
n += 1
# add only one event with zero asymmetry value and multiplicity 2
st[n] = spikes1[i]
a[n] = 0
mp[n] = 2
st = st[:n+2]
a = a[:n+2]
mp = mp[:n+2]
st[0] = t_start
st[len(st)-1] = t_end
if N1 + N2 > 0:
a[0] = a[1]
a[len(a)-1] = a[len(a)-2]
mp[0] = mp[1]
mp[len(mp)-1] = mp[len(mp)-2]
else:
a[0] = 1
a[1] = 1
return st, a, mp
|
