diff --git a/plotstyle.py b/plotstyle.py index 694bdb0..5cbf30e 100644 --- a/plotstyle.py +++ b/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. diff --git a/simulations/lecture/randomwalkneuron.py b/simulations/lecture/randomwalkneuron.py new file mode 100644 index 0000000..da248ce --- /dev/null +++ b/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() diff --git a/simulations/lecture/randomwalkone.py b/simulations/lecture/randomwalkone.py index 6fec50f..9b4724a 100644 --- a/simulations/lecture/randomwalkone.py +++ b/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<1.0-p] = -1 x = np.hstack(((0.0,), np.cumsum(steps))) return x -if __name__ == "__main__": - fig, ax = plt.subplots() - fig.subplots_adjust(**adjust_fs(fig, right=1.0)) - rng = np.random.RandomState(52281) - nmax = 80 +def plot_random_walk(ax, nmax, p, rng, ymin=-20.0, ymax=20.0): nn = np.linspace(0.0, nmax, 200) - ax.fill_between(nn, np.sqrt(nn), -np.sqrt(nn), **fsAa) - ax.axhline(0.0, **lsAm) + m = 2*p-1 + 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_yticks(np.arange(-20, 21, 10)) + ax.set_ylim(ymin, ymax) + 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() diff --git a/simulations/lecture/randomwalktwo.py b/simulations/lecture/randomwalktwo.py index 39a58b9..79221d1 100644 --- a/simulations/lecture/randomwalktwo.py +++ b/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) diff --git a/simulations/lecture/simulations.tex b/simulations/lecture/simulations.tex index 4c7cdac..f2deb1d 100644 --- a/simulations/lecture/simulations.tex +++ b/simulations/lecture/simulations.tex @@ -212,31 +212,69 @@ the event does not occur. \end{exercise} -\subsection{Random walks} -A basic concept for stochastic models is the random walk. A walker -starts at some initial position $x_0$. At every iteration $n$ it takes -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 -left ($\delta_n = -1$) with probability $P_l = 1-P_r$. - -For a symmetric random walk (\figref{randomwalkonefig}), the -probabilities to step to the left or to the right are the same, $P_r = -P_l = \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} = 1$. Since -each step is independent of the steps before, the variances of each -step add up such that the total variance of the walker position at -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}$. +The gamma distribution (\code{gamrnd()}) phenomenologically describes +various types of interspike interval dsitributions +(chapter~\ref{pointprocesseschapter}). scale and shape. exercise. + +\subsection{Random integers} +\code{randi()} + + +\section{Random walks} +A basic concept for stochastic models is the \enterm{random walk}. A +walker starts at some initial position $x_0$. At every iteration $n$ +it takes 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$ +or to the left ($\delta_n = -1$) with probability $q = 1-p$. \begin{figure}[tb] \includegraphics{randomwalkone} \titlecaption{\label{randomwalkonefig} Random walk in one - dimension.}{Twelve symmetric random walks starting at $x_0=0$. The - ensemble average of the position of the random walkers is zero - (blue line) and the standard deviation increases with the square - root of the number of steps taken (shaded area).} + dimension.}{Twelve random walks starting at $x_0=0$. The ensemble + average of the position of the random walkers marked by the blue + line and the corresponding standard deviation by the shaded + 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 - iteration take a step either up, down, to the left, or to the - right.} + dimensions.}{Fivehundred steps of five random walkers starting at + the origin (black dot) and taking steps of distance one either up, + 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 -various types of interspike interval dsitributions -(chapter~\ref{pointprocesseschapter}). scale and shape. exercise. - -\subsection{Random integers} -\code{randi()} +\begin{exercise}{}{} + Random walk in 2D exercise. +\end{exercise} %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%