[simulations] added random walk neuron figures

plotstyle.py

@@ -52,6 +52,8 @@ lsSpine = {'c': colors['black'], 'linestyle': '-', 'linewidth': 1, 'clip_on': Fa
 lsGrid = {'c': colors['gray'], 'linestyle': '--', 'linewidth': 1}
 lsMarker = {'c': colors['black'], 'linestyle': '-', 'linewidth': 2}
 
+psMarker = dict({'color': colors['black'], 'linestyle': 'none'}, **largemarker)
+
 # line (ls), point (ps), and fill styles (fs).
 
 # Each style is derived from a main color as indicated by the capital letter.

simulations/lecture/randomwalkneuron.py

@@ -0,0 +1,58 @@
+import numpy as np
+import matplotlib.pyplot as plt
+import matplotlib.gridspec as gridspec
+from plotstyle import *
+
+def random_walk(n, p, rng):
+    steps = rng.rand(n)
+    steps[steps>=1.0-p] = 1
+    steps[steps<1.0-p] = -1
+    x = np.hstack(((0.0,), np.cumsum(steps)))
+    return x
+
+def random_walk_neuron(n, p, thresh, rng):
+    x = random_walk(n, p, rng)
+    spikes = []
+    j = -1
+    i = 1
+    while i > 0:
+        if j > 0:
+            spikes.append(j+i)
+        j += i + 1
+        x[j:] -= x[j]
+        i = np.argmax(x[j:] >= thresh - 0.1)
+    return x, spikes
+    
+
+if __name__ == "__main__":
+    fig = plt.figure()
+    spec = gridspec.GridSpec(nrows=1, ncols=2, width_ratios=[2, 1], wspace=0.4,
+                             **adjust_fs(fig, left=6.5, right=1.5))
+    ax = fig.add_subplot(spec[0, 0])
+    rng = np.random.RandomState(52281)
+    p = 0.6
+    thresh = 10.0
+    nmax = 500
+    n = np.arange(0.0, nmax+1, 1.0)
+    x, spikes = random_walk_neuron(nmax, p, thresh, rng)
+    ax.axhline(0.0, **lsGrid)
+    ax.axhline(thresh, **lsAm)
+    ax.plot(n, x, **lsBm)
+    for tspike in spikes:
+        ax.plot([tspike, tspike], [12.0, 16.0], **lsC)
+    ax.set_xlabel('Time')
+    ax.set_ylabel('Potential')
+    ax.set_xlim(0, nmax)
+    ax.set_ylim(-10, 17)
+    ax.set_yticks(np.arange(-10, 11, 10))
+    ax = fig.add_subplot(spec[0, 1])
+    nmax = 100000
+    x, spikes = random_walk_neuron(nmax, p, thresh, rng)
+    isis = np.diff(spikes)
+    ax.hist(isis, np.arange(0.0, 151.0, 10.0), **fsAs)
+    ax.set_xlabel('ISI')
+    ax.set_ylabel('Count')
+    ax.set_xticks(np.arange(0, 151, 50))
+    ax.set_yticks(np.arange(0, 401, 100))
+    fig.savefig("randomwalkneuron.pdf")
+    plt.close()

simulations/lecture/randomwalkone.py

@@ -4,31 +4,39 @@ from plotstyle import *
 
 def random_walk(n, p, rng):
     steps = rng.rand(n)
-    steps[steps>=p] = 1
+    steps[steps>=1.0-p] = 1
-    steps[steps<p] = -1
+    steps[steps<1.0-p] = -1
     x = np.hstack(((0.0,), np.cumsum(steps)))
     return x
 
-if __name__ == "__main__":
+def plot_random_walk(ax, nmax, p, rng, ymin=-20.0, ymax=20.0):
-    fig, ax = plt.subplots()
-    fig.subplots_adjust(**adjust_fs(fig, right=1.0))
-    rng = np.random.RandomState(52281)
-    nmax = 80
     nn = np.linspace(0.0, nmax, 200)
-    ax.fill_between(nn, np.sqrt(nn), -np.sqrt(nn), **fsAa)
+    m = 2*p-1
-    ax.axhline(0.0, **lsAm)
+    v = 4*p*(1-p)
+    ax.fill_between(nn, m*nn+np.sqrt(nn*v), m*nn-np.sqrt(nn*v), **fsAa)
+    ax.plot([0.0, nmax], [0.0, m*nmax], **lsAm)
     lcs = [colors['red'], colors['orange'],
            colors['yellow'], colors['green']]
     n = np.arange(0.0, nmax+1, 1.0)
     for k in range(12):
-        x = random_walk(nmax, 0.5, rng)
+        x = random_walk(nmax, p, rng)
         ls = dict(**lsAm)
         ls['color'] = lcs[k%len(lcs)]
         ax.plot(n, x, **ls)
     ax.set_xlabel('Iteration $n$')
     ax.set_ylabel('Position $x_n$')
     ax.set_xlim(0, nmax)
-    ax.set_ylim(-20, 20)
+    ax.set_ylim(ymin, ymax)
-    ax.set_yticks(np.arange(-20, 21, 10))
+    ax.set_yticks(np.arange(ymin, ymax+1.0, 10.0))
+
+    
+if __name__ == "__main__":
+    fig, (ax1, ax2) = plt.subplots(2, 1, figsize=cm_size(figure_width, 2.0*figure_height))
+    fig.subplots_adjust(**adjust_fs(fig, right=1.0))
+    rng = np.random.RandomState(52281)
+    plot_random_walk(ax1, 80, 0.5, rng)
+    ax1.text(5.0, 20.0, 'symmetric', va='center')
+    plot_random_walk(ax2, 80, 0.6, rng, -10.0, 30.0)
+    ax2.text(5.0, 30.0, 'with drift', va='center')
     fig.savefig("randomwalkone.pdf")
     plt.close()

simulations/lecture/randomwalktwo.py

@@ -27,6 +27,7 @@ if __name__ == "__main__":
         ls = dict(**lsAm)
         ls['color'] = lcs[k%len(lcs)]
         ax.plot(x, y, **ls)
+    ax.plot([0], [0], **psMarker)
     ax.set_xlabel('Position $x_n$')
     ax.set_ylabel('Position $y_n$')
     ax.set_xlim(-30, 30)

simulations/lecture/simulations.tex

@@ -212,31 +212,69 @@ the event does not occur.
 \end{exercise}
 
 
-\subsection{Random walks}
+The gamma distribution (\code{gamrnd()}) phenomenologically describes
-A basic concept for stochastic models is the random walk. A walker
+various types of interspike interval dsitributions
-starts at some initial position $x_0$. At every iteration $n$ it takes
+(chapter~\ref{pointprocesseschapter}). scale and shape. exercise.
-a step $\delta_n$: $x_n = x_{n-1} + \delta_n$. The walker either steps
+
-to the right ($\delta_n = +1$) with some probability $P_r$ or to the
+\subsection{Random integers}
-left ($\delta_n = -1$) with probability $P_l = 1-P_r$.
+\code{randi()}
-
+
-For a symmetric random walk (\figref{randomwalkonefig}), the
+
-probabilities to step to the left or to the right are the same, $P_r =
+\section{Random walks}
-P_l = \frac{1}{2}$. The average position of many walkers is then
+A basic concept for stochastic models is the \enterm{random walk}. A
-independent of the iteration step and equals the initial position
+walker starts at some initial position $x_0$. At every iteration $n$
-$x_0$. The variance of each step is $\sigma_{\delta_n} = 1$. Since
+it takes a step $\delta_n$: $x_n = x_{n-1} + \delta_n$. The walker
-each step is independent of the steps before, the variances of each
+either steps to the right ($\delta_n = +1$) with some probability $p$
-step add up such that the total variance of the walker position at
+or to the left ($\delta_n = -1$) with probability $q = 1-p$.
-step $n$ equals $\sigma_{x_n}^2 = n$. The standard deviation of the
-walker positions grows with the square root of the number of iteration
-steps $n$: $\sigma_{x_n} = \sqrt{n}$.
 
 \begin{figure}[tb]
   \includegraphics{randomwalkone}
   \titlecaption{\label{randomwalkonefig} Random walk in one
-    dimension.}{Twelve symmetric random walks starting at $x_0=0$. The
+    dimension.}{Twelve random walks starting at $x_0=0$. The ensemble
-    ensemble average of the position of the random walkers is zero
+    average of the position of the random walkers marked by the blue
-    (blue line) and the standard deviation increases with the square
+    line and the corresponding standard deviation by the shaded
-    root of the number of steps taken (shaded area).}
+    area. The standard deviation increases with the square root of the
+    number of steps taken. The top plot shows a symmetric random walk
+    where the probability of a step to the right equals the one for a
+    step to the left ($p=q=\frac{1}{2}$). The bottom plot shows a
+    random walk with drift. Here the probability of taking a step to
+    the right ($p=0.6$) is larger than the probability for a step to
+    the left ($q=0.4$).}
+\end{figure}
+
+For a \enterm[random walk!symmetric]{symmetric random walk}
+(\figref{randomwalkonefig} top), the probabilities for a step to the
+left or to the right are the same, $p = q = \frac{1}{2}$. The average
+position of many walkers is then independent of the iteration step and
+equals the initial position $x_0$. The variance of each step is
+$\sigma_{\delta_n}^2 = p \cdot 1^2 + q \cdot(-1)^2 = 1$. Since each
+step is independent of the steps before, the variances of all the
+steps add up such that the total variance of the walker position at
+step $n$ equals $\sigma_{x_n}^2 = \sum_{i=1}^n \sigma_{\delta_n}^2 =
+n$. The standard deviation of the walker positions grows with the
+square root of the number of iteration steps $n$: $\sigma_{x_n} =
+\sqrt{n}$.
+
+\enterm[random walk!with drift]{Random walks with drift} are biased to
+one direction (\figref{randomwalkonefig} bottom). The probability to
+step to the left does not equal the probability to step to the right,
+$p \ne q$. Consequently, the average step is no longer zero. Instead
+it equals $\mu_{\delta_n} = p \cdot 1 + q \cdot (-1) = 2p - 1$. The
+average position grows proportional to the number of iterations:
+$\mu_{x_n} = n(2p - 1)$. The variance of each step turns out to be
+$\sigma_{\delta_n}^2 = 4pq$ (this indeed equals one for
+$p=q=\frac{1}{2}$). Again, because of the independence of the steps,
+the variance of the position at the $n$-th iteration, $\sigma_{x_n}^2
+= 4pqn$, is proportional to $n$ and the standard deviation,
+$\sigma_{x_n} = 2\sqrt{pqn}$, grows with the square root of $n$.
+
+\begin{exercise}{}{}
+  Random walk exercise.
+\end{exercise}
+
+\begin{figure}[tb]
+  \includegraphics{randomwalkneuron}
+  \titlecaption{\label{randomwalkneuronfig} Random walk spiking neuron.}{.}
 \end{figure}
 
 The sub-threshold membrane potential of a neuron can be modeled as a
@@ -253,23 +291,23 @@ times are again exponentially distributed. ??? Is that so ???
 
 % Gerstein and Mandelbrot 1964
 
-Higher dimensions (\figref{randomwalktwofig}), Brownian motion. Crystal growth.
+\begin{exercise}{}{}
+  Random walk neuron exercise.
+\end{exercise}
 
 \begin{figure}[tb]
   \includegraphics{randomwalktwo}
   \titlecaption{\label{randomwalktwofig} Random walk in two
-    dimensions.}{Fivehundred steps of five random walkers that in each
+    dimensions.}{Fivehundred steps of five random walkers starting at
-    iteration take a step either up, down, to the left, or to the
+    the origin (black dot) and taking steps of distance one either up,
-    right.}
+    down, to the left, or to the right with equal probability.}
 \end{figure}
 
+Higher dimensions (\figref{randomwalktwofig}), Brownian motion. Crystal growth.
 
-The gamma distribution (\code{gamrnd()}) phenomenologically describes
+\begin{exercise}{}{}
-various types of interspike interval dsitributions
+  Random walk in 2D exercise.
-(chapter~\ref{pointprocesseschapter}). scale and shape. exercise.
+\end{exercise}
-
-\subsection{Random integers}
-\code{randi()}
 
 
 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%