Add Spike Order example and Tutorial section

Adds an example computing spike order profile and the optimal
spike train order. Also adds a section on spike train order to the
tutorial.

3 files changed, 119 insertions, 1 deletions

diff --git a/doc/tutorial.rst b/doc/tutorial.rst
index aff03a8..377c0a2 100644
--- a/doc/tutorial.rst
+++ b/doc/tutorial.rst
@@ -231,3 +231,69 @@ The following example computes and plots the ISI- and SPIKE-distance matrix as w
     plt.title("SPIKE-Sync")
 
     plt.show()
+
+
+Quantifying Leaders and Followers: Spike Train Order
+---------------------------------------
+
+PySpike provides functionality to quantify how much a set of spike trains
+resembles a synfire pattern (ie perfect leader-follower pattern). For details
+on the algorithms please see
+`our article in NJP <http://iopscience.iop.org/article/10.1088/1367-2630/aa68c3>`_.
+
+The following example computes the Spike Order profile and Synfire Indicator
+of two Poissonian spike trains.
+
+.. code:: python
+    import numpy as np
+    from matplotlib import pyplot as plt
+    import pyspike as spk
+
+
+    st1 = spk.generate_poisson_spikes(1.0, [0, 20])
+    st2 = spk.generate_poisson_spikes(1.0, [0, 20])
+
+    d = spk.spike_directionality(st1, st2)
+
+    print "Spike Directionality of two Poissonian spike trains:", d
+
+    E = spk.spike_train_order_profile(st1, st2)
+
+    plt.figure()
+    x, y = E.get_plottable_data()
+    plt.plot(x, y, '-ob')
+    plt.ylim(-1.1, 1.1)
+    plt.xlabel("t")
+    plt.ylabel("E")
+    plt.title("Spike Train Order Profile")
+
+    plt.show()
+
+Additionally, PySpike can also compute the optimal ordering of the spike trains,
+ie the ordering that most resembles a synfire pattern. The following example
+computes the optimal order of a set of 20 Poissonian spike trains:
+
+.. code:: python
+
+    M = 20
+    spike_trains = [spk.generate_poisson_spikes(1.0, [0, 100]) for m in xrange(M)]
+
+    F_init = spk.spike_train_order(spike_trains)
+    print "Initial Synfire Indicator for 20 Poissonian spike trains:", F_init
+
+    D_init = spk.spike_directionality_matrix(spike_trains)
+    phi, _ = spk.optimal_spike_train_sorting(spike_trains)
+    F_opt = spk.spike_train_order(spike_trains, indices=phi)
+    print "Synfire Indicator of optimized spike train sorting:", F_opt
+
+    D_opt = spk.permutate_matrix(D_init, phi)
+
+    plt.figure()
+    plt.imshow(D_init)
+    plt.title("Initial Directionality Matrix")
+
+    plt.figure()
+    plt.imshow(D_opt)
+    plt.title("Optimized Directionality Matrix")
+
+    plt.show()
diff --git a/examples/spike_train_order.py b/examples/spike_train_order.py
new file mode 100644
index 0000000..3a42472
--- /dev/null
+++ b/examples/spike_train_order.py
@@ -0,0 +1,52 @@
+import numpy as np
+from matplotlib import pyplot as plt
+import pyspike as spk
+
+
+st1 = spk.generate_poisson_spikes(1.0, [0, 20])
+st2 = spk.generate_poisson_spikes(1.0, [0, 20])
+
+d = spk.spike_directionality(st1, st2)
+
+print "Spike Directionality of two Poissonian spike trains:", d
+
+E = spk.spike_train_order_profile(st1, st2)
+
+plt.figure()
+x, y = E.get_plottable_data()
+plt.plot(x, y, '-ob')
+plt.ylim(-1.1, 1.1)
+plt.xlabel("t")
+plt.ylabel("E")
+plt.title("Spike Train Order Profile")
+
+
+###### Optimize spike train order of 20 Random spike trains #######
+
+M = 20
+
+spike_trains = [spk.generate_poisson_spikes(1.0, [0, 100]) for m in xrange(M)]
+
+F_init = spk.spike_train_order(spike_trains)
+
+print "Initial Synfire Indicator for 20 Poissonian spike trains:", F_init
+
+D_init = spk.spike_directionality_matrix(spike_trains)
+
+phi, _ = spk.optimal_spike_train_sorting(spike_trains)
+
+F_opt = spk.spike_train_order(spike_trains, indices=phi)
+
+print "Synfire Indicator of optimized spike train sorting:", F_opt
+
+D_opt = spk.permutate_matrix(D_init, phi)
+
+plt.figure()
+plt.imshow(D_init)
+plt.title("Initial Directionality Matrix")
+
+plt.figure()
+plt.imshow(D_opt)
+plt.title("Optimized Directionality Matrix")
+
+plt.show()
diff --git a/pyspike/spike_directionality.py b/pyspike/spike_directionality.py
index 2cb2bf1..248862c 100644
--- a/pyspike/spike_directionality.py
+++ b/pyspike/spike_directionality.py
@@ -501,7 +501,7 @@ def optimal_spike_train_sorting(spike_trains,  indices=None, interval=None,
     D = spike_directionality_matrix(spike_trains, normalize=False,
                                     indices=indices, interval=interval,
                                     max_tau=max_tau)
-    return _optimal_spike_train_order_from_matrix(D, full_output)
+    return _optimal_spike_train_sorting_from_matrix(D, full_output)
 
 ############################################################
 # permutate_matrix