more directionality functions

+ spike delay asymmetry profile
+ spike delay dual profile

3 files changed, 102 insertions, 4 deletions

diff --git a/pyspike/directionality/__init__.py b/pyspike/directionality/__init__.py
index e6de1de..9a45b54 100644
--- a/pyspike/directionality/__init__.py
+++ b/pyspike/directionality/__init__.py
@@ -8,5 +8,6 @@ __all__ = ["spike_delay_asymmetry"]
 
 from spike_delay_asymmetry import spike_delay_asymmetry_profile, \
     spike_delay_asymmetry, spike_delay_asymmetry_profile_multi, \
-    spike_delay_asymmetry_matrix, optimal_asymmetry_order, \
-    optimal_asymmetry_order_from_D, reorder_asymmetry_matrix
+    spike_delay_asymmetry_matrix, spike_delay_asymmetry_full, \
+    optimal_asymmetry_order, optimal_asymmetry_order_from_D, \
+    permutate_asymmetry_matrix
diff --git a/pyspike/directionality/cython/cython_directionality.pyx b/pyspike/directionality/cython/cython_directionality.pyx
index f5ea752..00e3b86 100644
--- a/pyspike/directionality/cython/cython_directionality.pyx
+++ b/pyspike/directionality/cython/cython_directionality.pyx
@@ -130,6 +130,52 @@ def spike_delay_asymmetry_profile_cython(double[:] spikes1, double[:] spikes2,
 
 
 ############################################################
+# spike_delay_dual_profile_cython
+############################################################
+def spike_delay_dual_profile_cython(double[:] spikes1,
+                                    double[:] spikes2,
+                                    double t_start, double t_end,
+                                    double max_tau):
+
+    cdef int N1 = len(spikes1)
+    cdef int N2 = len(spikes2)
+    cdef int i = -1
+    cdef int j = -1
+    cdef double[:] a1 = np.zeros(N1)  # asymmetry values
+    cdef double[:] a2 = np.zeros(N2)  # asymmetry values
+    cdef double interval = t_end - t_start
+    cdef double tau
+    while i + j < N1 + N2 - 2:
+        if (i < N1-1) and (j == N2-1 or spikes1[i+1] < spikes2[j+1]):
+            i += 1
+            tau = get_tau(spikes1, spikes2, i, j, interval, max_tau)
+            if j > -1 and spikes1[i]-spikes2[j] < tau:
+                # coincidence between the current spike and the previous spike
+                # spike from spike train 1 after spike train 2
+                # leading spike gets +1, following spike -1
+                a1[i] = -1
+                a2[j] = +1
+        elif (j < N2-1) and (i == N1-1 or spikes1[i+1] > spikes2[j+1]):
+            j += 1
+            tau = get_tau(spikes1, spikes2, i, j, interval, max_tau)
+            if i > -1 and spikes2[j]-spikes1[i] < tau:
+                # coincidence between the current spike and the previous spike
+                # spike from spike train 1 before spike train 2
+                # leading spike gets +1, following spike -1
+                a1[i] = +1
+                a2[j] = -1
+        else:   # spikes1[i+1] = spikes2[j+1]
+            # advance in both spike trains
+            j += 1
+            i += 1
+            # equal spike times: zero asymmetry value
+            a1[i] = 0
+            a2[j] = 0
+
+    return a1, a2
+
+
+############################################################
 # spike_delay_asymmetry_cython
 ############################################################
 def spike_delay_asymmetry_cython(double[:] spikes1, double[:] spikes2,
diff --git a/pyspike/directionality/spike_delay_asymmetry.py b/pyspike/directionality/spike_delay_asymmetry.py
index 7d59601..0676f6d 100644
--- a/pyspike/directionality/spike_delay_asymmetry.py
+++ b/pyspike/directionality/spike_delay_asymmetry.py
@@ -147,6 +147,57 @@ def spike_delay_asymmetry_matrix(spike_trains, indices=None,
 
 
 ############################################################
+# spike_delay_asymmetry_full
+############################################################
+def spike_delay_asymmetry_full(spike_trains, indices=None,
+                               interval=None, max_tau=None):
+    """ Computes the spike train symmetry value for each spike in each spike
+    train.
+    """
+    if indices is None:
+        indices = np.arange(len(spike_trains))
+    indices = np.array(indices)
+    # check validity of indices
+    assert (indices < len(spike_trains)).all() and (indices >= 0).all(), \
+        "Invalid index list."
+    # list of arrays for reulting asymmetry values
+    asymmetry_list = [np.zeros_like(st.spikes) for st in spike_trains]
+    # generate a list of possible index pairs
+    pairs = [(indices[i], j) for i in range(len(indices))
+             for j in indices[i+1:]]
+
+    # cython implementation
+    try:
+        from cython.cython_directionality import \
+            spike_delay_dual_profile_cython as \
+            sda_dual_profile_impl
+    except ImportError:
+        raise NotImplementedError()
+#         if not(pyspike.disable_backend_warning):
+#             print("Warning: spike_distance_cython not found. Make sure that \
+# PySpike is installed by running\n 'python setup.py build_ext --inplace'!\n \
+# Falling back to slow python backend.")
+#         # use python backend
+#         from cython.python_backend import coincidence_python \
+#             as coincidence_profile_impl
+
+    if max_tau is None:
+        max_tau = 0.0
+
+    for i, j in pairs:
+        a1, a2 = sda_dual_profile_impl(spike_trains[i].spikes,
+                                       spike_trains[j].spikes,
+                                       spike_trains[i].t_start,
+                                       spike_trains[i].t_end,
+                                       max_tau)
+        asymmetry_list[i] += a1
+        asymmetry_list[j] += a2
+    for a in asymmetry_list:
+        a /= len(spike_trains)-1
+    return asymmetry_list
+
+
+############################################################
 # optimal_asymmetry_order_from_D
 ############################################################
 def optimal_asymmetry_order_from_D(D, full_output=False):
@@ -188,7 +239,7 @@ def optimal_asymmetry_order_from_D(D, full_output=False):
 
 
 ############################################################
-# _optimal_asymmetry_order
+# optimal_asymmetry_order
 ############################################################
 def optimal_asymmetry_order(spike_trains,  indices=None, interval=None,
                             max_tau=None, full_output=False):
@@ -203,7 +254,7 @@ def optimal_asymmetry_order(spike_trains,  indices=None, interval=None,
 ############################################################
 # reorder_asymmetry_matrix
 ############################################################
-def reorder_asymmetry_matrix(D, p):
+def permutate_asymmetry_matrix(D, p):
     N = len(D)
     D_p = np.empty_like(D)
     for n in xrange(N):