Merge branch 'master' of whale.am28.uni-tuebingen.de:scientificComputing

This commit is contained in:
Jan Grewe 2017-01-21 16:55:02 +01:00
commit a32a22b35b
29 changed files with 5399 additions and 58 deletions

View File

@ -52,10 +52,10 @@ endobj
<?adobe-xap-filters esc="CRLF"?> <?adobe-xap-filters esc="CRLF"?>
<x:xmpmeta xmlns:x='adobe:ns:meta/' x:xmptk='XMP toolkit 2.9.1-13, framework 1.6'> <x:xmpmeta xmlns:x='adobe:ns:meta/' x:xmptk='XMP toolkit 2.9.1-13, framework 1.6'>
<rdf:RDF xmlns:rdf='http://www.w3.org/1999/02/22-rdf-syntax-ns#' xmlns:iX='http://ns.adobe.com/iX/1.0/'> <rdf:RDF xmlns:rdf='http://www.w3.org/1999/02/22-rdf-syntax-ns#' xmlns:iX='http://ns.adobe.com/iX/1.0/'>
<rdf:Description rdf:about='94827694-b0d9-11f0-0000-86f60cc553dd' xmlns:pdf='http://ns.adobe.com/pdf/1.3/' pdf:Producer='Artifex Ghostscript 8.54'/> <rdf:Description rdf:about='525cf0ae-142c-11f2-0000-86f60cc553dd' xmlns:pdf='http://ns.adobe.com/pdf/1.3/' pdf:Producer='Artifex Ghostscript 8.54'/>
<rdf:Description rdf:about='94827694-b0d9-11f0-0000-86f60cc553dd' xmlns:xap='http://ns.adobe.com/xap/1.0/' xap:ModifyDate='2015-10-22' xap:CreateDate='2015-10-22'><xap:CreatorTool>Artifex Ghostscript 8.54 PDF Writer</xap:CreatorTool></rdf:Description> <rdf:Description rdf:about='525cf0ae-142c-11f2-0000-86f60cc553dd' xmlns:xap='http://ns.adobe.com/xap/1.0/' xap:ModifyDate='2017-01-16' xap:CreateDate='2017-01-16'><xap:CreatorTool>Artifex Ghostscript 8.54 PDF Writer</xap:CreatorTool></rdf:Description>
<rdf:Description rdf:about='94827694-b0d9-11f0-0000-86f60cc553dd' xmlns:xapMM='http://ns.adobe.com/xap/1.0/mm/' xapMM:DocumentID='94827694-b0d9-11f0-0000-86f60cc553dd'/> <rdf:Description rdf:about='525cf0ae-142c-11f2-0000-86f60cc553dd' xmlns:xapMM='http://ns.adobe.com/xap/1.0/mm/' xapMM:DocumentID='525cf0ae-142c-11f2-0000-86f60cc553dd'/>
<rdf:Description rdf:about='94827694-b0d9-11f0-0000-86f60cc553dd' xmlns:dc='http://purl.org/dc/elements/1.1/' dc:format='application/pdf'><dc:title><rdf:Alt><rdf:li xml:lang='x-default'>/tmp/tpd0b45dc9_ff5a_4aa8_90bd_50aa8e8237b6.ps</rdf:li></rdf:Alt></dc:title></rdf:Description> <rdf:Description rdf:about='525cf0ae-142c-11f2-0000-86f60cc553dd' xmlns:dc='http://purl.org/dc/elements/1.1/' dc:format='application/pdf'><dc:title><rdf:Alt><rdf:li xml:lang='x-default'>/tmp/tp6e6add58_0ada_4a6c_ba7e_63eeba1cbd7c.ps</rdf:li></rdf:Alt></dc:title></rdf:Description>
</rdf:RDF> </rdf:RDF>
</x:xmpmeta> </x:xmpmeta>
@ -65,10 +65,10 @@ endstream
endobj endobj
2 0 obj 2 0 obj
<</Producer(Artifex Ghostscript 8.54) <</Producer(Artifex Ghostscript 8.54)
/CreationDate(D:20151022150138) /CreationDate(D:20170116181818)
/ModDate(D:20151022150138) /ModDate(D:20170116181818)
/Creator(MATLAB, The MathWorks, Inc. Version 8.3.0.532 \(R2014a\). Operating System: Linux 3.13.0-24-generic #47-Ubuntu SMP Fri May 2 23:30:00 UTC 2014 x86_64.) /Creator(MATLAB, The MathWorks, Inc. Version 8.3.0.532 \(R2014a\). Operating System: Linux 3.13.0-24-generic #47-Ubuntu SMP Fri May 2 23:30:00 UTC 2014 x86_64.)
/Title(/tmp/tpd0b45dc9_ff5a_4aa8_90bd_50aa8e8237b6.ps)>>endobj /Title(/tmp/tp6e6add58_0ada_4a6c_ba7e_63eeba1cbd7c.ps)>>endobj
xref xref
0 12 0 12
0000000000 65535 f 0000000000 65535 f
@ -85,7 +85,7 @@ xref
0000001375 00000 n 0000001375 00000 n
trailer trailer
<< /Size 12 /Root 1 0 R /Info 2 0 R << /Size 12 /Root 1 0 R /Info 2 0 R
/ID [<8FE161FFCC6D1C11BAD3CE59BA0E3F32><8FE161FFCC6D1C11BAD3CE59BA0E3F32>] /ID [<F0B160FA5B0C72A4E47D89F3F547B465><F0B160FA5B0C72A4E47D89F3F547B465>]
>> >>
startxref startxref
3070 3070

View File

@ -50,7 +50,7 @@ hold on;
bar(b, h, 'facecolor', 'b'); bar(b, h, 'facecolor', 'b');
bar(b(b>=rq), h(b>=rq), 'facecolor', 'r'); bar(b(b>=rq), h(b>=rq), 'facecolor', 'r');
plot( [rd rd], [0 4], 'r', 'linewidth', 2 ); plot( [rd rd], [0 4], 'r', 'linewidth', 2 );
xlim([-0.2 0.2]) xlim([-0.25 0.25])
xlabel('Correlation coefficient'); xlabel('Correlation coefficient');
ylabel('Probability density of H0'); ylabel('Probability density of H0');
hold off; hold off;

View File

@ -15,7 +15,7 @@
\else \else
\newcommand{\stitle}{} \newcommand{\stitle}{}
\fi \fi
\header{{\bfseries\large \"Ubung 3\stitle}}{{\bfseries\large Statistik}}{{\bfseries\large 21. Oktober, 2015}} \header{{\bfseries\large \"Ubung\stitle}}{{\bfseries\large Bootstrap}}{{\bfseries\large 17. Januar, 2017}}
\firstpagefooter{Prof. Dr. Jan Benda}{Phone: 29 74573}{Email: \firstpagefooter{Prof. Dr. Jan Benda}{Phone: 29 74573}{Email:
jan.benda@uni-tuebingen.de} jan.benda@uni-tuebingen.de}
\runningfooter{}{\thepage}{} \runningfooter{}{\thepage}{}
@ -84,10 +84,13 @@ jan.benda@uni-tuebingen.de}
\input{instructions} \input{instructions}
\begin{questions} \begin{questions}
\question \qt{Bootstrap des Standardfehlers} \question \qt{Bootstrap des Standardfehlers}
Wir wollen den Standardfehler, die Standardabweichung des Mittelwerts,
eines Datensatze mit Hilfe der Bootstrapmethode berechnen und mit der
Formel ``Standardabweichung geteilt durch Wurzel aus $n$''
vergleichen.
\begin{parts} \begin{parts}
\part Lade von Ilias die Datei \code{thymusglandweights.dat} herunter. \part Lade von Ilias die Datei \code{thymusglandweights.dat} herunter.
Darin befindet sich ein Datensatz vom Gewicht der Thymus Dr\"use in 14-Tage alten Darin befindet sich ein Datensatz vom Gewicht der Thymus Dr\"use in 14-Tage alten
@ -112,8 +115,11 @@ jan.benda@uni-tuebingen.de}
\end{solution} \end{solution}
\continue
\question \qt{Student t-Verteilung} \question \qt{Student t-Verteilung}
Durch Standardabweichungen normierte Mittelwerte sind nicht Gaussverteilt,
wenn beide aus Normalverteilten Daten abgesch\"atzt werden.
Die Verteilung von $t=\bar x/(\sigma_x/\sqrt{m})$ folgt vielmehr
der Student t-Verteilung.
\begin{parts} \begin{parts}
\part Erzeuge 100000 normalverteilte Zufallszahlen. \part Erzeuge 100000 normalverteilte Zufallszahlen.
\part Ziehe daraus 1000 Stichproben vom Umfang $m=3$, 5, 10, oder 50. \part Ziehe daraus 1000 Stichproben vom Umfang $m=3$, 5, 10, oder 50.
@ -123,6 +129,7 @@ dieser Mittelwerte.
\part Berechne ausserdem die Gr\"o{\ss}e $t=\bar x/(\sigma_x/\sqrt{m})$ \part Berechne ausserdem die Gr\"o{\ss}e $t=\bar x/(\sigma_x/\sqrt{m})$
(Standardabweichung $\sigma_x$) und vergleiche diese mit der Normalverteilung mit Standardabweichung Eins. Ist $t$ normalverteilt, bzw. unter welchen Bedingungen ist $t$ normalverteilt? (Standardabweichung $\sigma_x$) und vergleiche diese mit der Normalverteilung mit Standardabweichung Eins. Ist $t$ normalverteilt, bzw. unter welchen Bedingungen ist $t$ normalverteilt?
\end{parts} \end{parts}
\newsolutionpage
\begin{solution} \begin{solution}
\lstinputlisting{tdistribution.m} \lstinputlisting{tdistribution.m}
\includegraphics[width=1\textwidth]{tdistribution-n03}\\ \includegraphics[width=1\textwidth]{tdistribution-n03}\\
@ -132,7 +139,10 @@ dieser Mittelwerte.
\end{solution} \end{solution}
\question \qt{Korrelationen} \continue
\question \qt{Permutationstest}
Wir wollen die Signifikanz einer Korrelation durch einen
Permutationstest bestimmen.
\begin{parts} \begin{parts}
\part Erzeuge 1000 korrelierte Zufallszahlen $x$, $y$ durch \part Erzeuge 1000 korrelierte Zufallszahlen $x$, $y$ durch
\begin{verbatim} \begin{verbatim}

View File

@ -4,38 +4,3 @@
{\large Jan Grewe, Jan Benda}\\[-3ex] {\large Jan Grewe, Jan Benda}\\[-3ex]
Abteilung Neuroethologie \hfill --- \hfill Institut f\"ur Neurobiologie \hfill --- \hfill \includegraphics[width=0.28\textwidth]{UT_WBMW_Black_RGB} \\ Abteilung Neuroethologie \hfill --- \hfill Institut f\"ur Neurobiologie \hfill --- \hfill \includegraphics[width=0.28\textwidth]{UT_WBMW_Black_RGB} \\
\end{center} \end{center}
\ifprintanswers%
\else
% Die folgenden Aufgaben dienen der Wiederholung, \"Ubung und
% Selbstkontrolle und sollten eigenst\"andig bearbeitet und gel\"ost
% werden. Die L\"osung soll in Form eines einzelnen Skriptes (m-files)
% im ILIAS hochgeladen werden. Jede Aufgabe sollte in einer eigenen
% ``Zelle'' gel\"ost sein. Die Zellen \textbf{m\"ussen} unabh\"angig
% voneinander ausf\"uhrbar sein. Das Skript sollte nach dem Muster:
% ``variablen\_datentypen\_\{nachname\}.m'' benannt werden
% (z.B. variablen\_datentypen\_mueller.m).
\begin{itemize}
\item \"Uberzeuge dich von jeder einzelnen Zeile deines Codes, dass
sie auch wirklich das macht, was sie machen soll! Teste dies mit
kleinen Beispielen direkt in der Kommandozeile.
\item Versuche die L\"osungen der Aufgaben m\"oglichst in
sinnvolle kleine Funktionen herunterzubrechen.
Sobald etwas \"ahnliches mehr als einmal berechnet werden soll,
lohnt es sich eine Funktion daraus zu schreiben!
\item Teste rechenintensive \code{for} Schleifen, Vektoren, Matrizen
zuerst mit einer kleinen Anzahl von Wiederholungen oder kleiner
Gr\"o{\ss}e, und benutze erst am Ende, wenn alles \"uberpr\"uft
ist, eine gro{\ss}e Anzahl von Wiederholungen oder Elementen, um eine gute
Statistik zu bekommen.
\item Benutze die Hilfsfunktion von \code{matlab} (\code{help
commando} oder \code{doc commando}) und das Internet, um
herauszufinden, wie bestimmte \code{matlab} Funktionen zu verwenden
sind und was f\"ur M\"oglichkeiten sie bieten.
Auch zu inhaltlichen Konzepten bietet das Internet oft viele
Antworten!
\end{itemize}
\fi

View File

@ -4,10 +4,13 @@ all : pdf
include ../../chapter.mk include ../../chapter.mk
# script: include ../../slides.mk
pdf : chapter
clean : cleanchapter # script and slides:
pdf : chapter slides
clean : cleanchapter cleanslides
cleanall : clean cleanallchapter cleanallslides
cleanall : clean cleanchapter

File diff suppressed because one or more lines are too long

View File

@ -0,0 +1,61 @@
% Copyright 2007 by Till Tantau
%
% This file may be distributed and/or modified
%
% 1. under the LaTeX Project Public License and/or
% 2. under the GNU Public License.
%
% See the file doc/licenses/LICENSE for more details.
\usepackage{color}
\definecolor{karminrot}{RGB}{165,30,55}
\definecolor{gold}{RGB}{180,160,105}
\definecolor{anthrazit}{RGB}{50 ,65 ,75 }
\mode<presentation>
\setbeamercolor*{normal text}{fg=anthrazit,bg=white}
\setbeamercolor*{alerted text}{fg=anthrazit}
\setbeamercolor*{example text}{fg=anthrazit}
\setbeamercolor*{structure}{fg=gold,bg=karminrot}
\providecommand*{\beamer@bftext@only}{%
\relax
\ifmmode
\expandafter\beamer@bftext@warning
\else
\expandafter\bfseries
\fi
}
\providecommand*{\beamer@bftext@warning}{%
\ClassWarning{beamer}
{Cannot use bold for alerted text in math mode}%
}
\setbeamerfont{alerted text}{series=\beamer@bftext@only}
\setbeamercolor{palette primary}{fg=karminrot,bg=white}
\setbeamercolor{palette secondary}{fg=gold,bg=white}
\setbeamercolor{palette tertiary}{fg=anthrazit,bg=white}
\setbeamercolor{palette quaternary}{fg=black,bg=white}
\setbeamercolor{sidebar}{bg=karminrot!100}
\setbeamercolor{palette sidebar primary}{fg=karminrot}
\setbeamercolor{palette sidebar secondary}{fg=karminrot}
\setbeamercolor{palette sidebar tertiary}{fg=karminrot}
\setbeamercolor{palette sidebar quaternary}{fg=karminrot}
\setbeamercolor{item projected}{fg=black,bg=black!20}
\setbeamercolor*{block body}{}
\setbeamercolor*{block body alerted}{}
\setbeamercolor*{block body example}{}
\setbeamercolor*{block title}{parent=structure}
\setbeamercolor*{block title alerted}{parent=alerted text}
\setbeamercolor*{block title example}{parent=example text}
\setbeamercolor*{titlelike}{parent=structure}
\mode
<all>

View File

@ -0,0 +1,110 @@
\documentclass{beamer}
%%%%% title %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\title[]{Scientific Computing --- Bootstrap methods}
\author[]{Jan Benda}
\institute[]{Neuroethology}
\date[]{WS 16/17}
\titlegraphic{\includegraphics[width=0.3\textwidth]{UT_WBMW_Rot_RGB}}
%%%%% beamer %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\mode<presentation>
{
\usetheme{Singapore}
\setbeamercovered{opaque}
\usecolortheme{tuebingen}
\setbeamertemplate{navigation symbols}{}
\usefonttheme{default}
\useoutertheme{infolines}
% \useoutertheme{miniframes}
}
%\AtBeginSection[]
%{
% \begin{frame}<beamer>
% \begin{center}
% \Huge \insertsectionhead
% \end{center}
% \end{frame}
%}
\setbeamertemplate{blocks}[rounded][shadow=true]
\setcounter{tocdepth}{1}
%%%%% packages %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\usepackage[english]{babel}
\usepackage{amsmath}
\usepackage{bm}
\usepackage{pslatex} % nice font for pdf file
%\usepackage{multimedia}
\usepackage{dsfont}
\newcommand{\naZ}{\mathds{N}}
\newcommand{\gaZ}{\mathds{Z}}
\newcommand{\raZ}{\mathds{Q}}
\newcommand{\reZ}{\mathds{R}}
\newcommand{\reZp}{\mathds{R^+}}
\newcommand{\reZpN}{\mathds{R^+_0}}
\newcommand{\koZ}{\mathds{C}}
%%%% graphics %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\usepackage{graphicx}
\newcommand{\texpicture}[1]{{\sffamily\small\input{#1.tex}}}
%%%%% listings %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\usepackage{listings}
\lstset{
basicstyle=\ttfamily,
numbers=left,
showstringspaces=false,
language=Matlab,
commentstyle=\itshape\color{darkgray},
keywordstyle=\color{blue},
stringstyle=\color{green},
backgroundcolor=\color{blue!10},
breaklines=true,
breakautoindent=true,
columns=flexible,
frame=single,
captionpos=b,
xleftmargin=1em,
xrightmargin=1em,
aboveskip=10pt
}
\graphicspath{{figures/}}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\begin{document}
\begin{frame}[plain]
\frametitle{}
\vspace{-1cm}
\titlepage % erzeugt Titelseite
\end{frame}
\subsection{Bootstrap distribution}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\begin{frame}
\frametitle{Sampling distribution}
\includegraphics[width=1\textwidth]{srs1}
\end{frame}
\begin{frame}
\frametitle{Theoretical sampling distribution}
\includegraphics[width=1\textwidth]{srs2}
\end{frame}
\begin{frame}
\frametitle{Bootstrap distribution}
\includegraphics[width=1\textwidth]{srs3}\\[1ex]
\begin{itemize}
\item Resampling with replacement.
\item Size of resample equals size of original sample.
\end{itemize}
\end{frame}
\end{document}

View File

@ -1 +1 @@
printfaculty(5); printfactorial(5);

View File

@ -1,4 +1,4 @@
function x = mufactorial(n) function x = myfactorial(n)
% return the factorial of n % return the factorial of n
x = 1; x = 1;
for i = 1:n for i = 1:n

View File

@ -1,5 +1,5 @@
function printfaculty(n) function printfactorial(n)
% compute the faculty of n and print it % compute the factorial of n and print it
x = 1; x = 1;
for i = 1:n for i = 1:n
x = x * i; x = x * i;

View File

@ -26,6 +26,7 @@
\item Matrizen erstmal 2-D und nur kurz n-D \item Matrizen erstmal 2-D und nur kurz n-D
\item Zusammenfassung von Fehlermeldungen bezueglich Matrizen und Vektoren \item Zusammenfassung von Fehlermeldungen bezueglich Matrizen und Vektoren
\item Random-walk behutsam mit einer Schleife, dann zwei Schleifen. Plus Bildchen. \item Random-walk behutsam mit einer Schleife, dann zwei Schleifen. Plus Bildchen.
\item Random-Walk with drift is model for decision making! See resoures/GoldShadlen2007 paper, Figure 2. Das wird auch in Systems Neurosciene \texttt{9 - Kognitive Kontrolle II.pdf} behandelt.
\item Doppelte for-Schleife \item Doppelte for-Schleife
\item File output and input (load, save, fprintf, scanf) (extra chapter?) \item File output and input (load, save, fprintf, scanf) (extra chapter?)
\item help() und doc() \item help() und doc()

View File

@ -0,0 +1,11 @@
function [counts, cbins] = counthist(spikes, tmin, tmax, T, cmax)
tbins = tmin+T/2:T:tmax;
cbins = 0.5:cmax;
counts = zeros(1, length(cbins));
for k = 1:length(spikes)
times = spikes{k};
n = hist(times((times>=tmin)&(times<=tmax)), tbins);
counts = counts + hist(n, cbins);
end
counts = counts / sum(counts);
end

View File

@ -0,0 +1,23 @@
function [d, thresholds, true1s, false1s, true2s, false2s, pratio] = discriminability(spikes1, spikes2, tmax, T, cmax)
[c1, b1] = counthist(spikes1, 0.0, tmax, T, cmax);
[c2, b2] = counthist(spikes2, 0.0, tmax, T, cmax);
thresholds = 0:cmax;
true1s = zeros(length(thresholds), 1);
true2s = zeros(length(thresholds), 1);
false1s = zeros(length(thresholds), 1);
false2s = zeros(length(thresholds), 1);
for k = 1:length(thresholds)
th = thresholds(k);
t1 = sum(c1(b1<=th));
f1 = sum(c1(b1>th));
t2 = sum(c2(b2>=th));
f2 = sum(c2(b2<th));
true1s(k) = t1;
true2s(k) = t2;
false1s(k) = f1;
false2s(k) = f2;
end
%pratio = (true1s + true2s)./(false1s+false2s);
pratio = (true1s + true2s)/2;
d = max(pratio);
end

View File

@ -0,0 +1,71 @@
%% general settings for the model neuron:
trials = 10;
tmax = 50.0;
D = 0.01;
%% generate and plot spiketrains for two inputs:
I1 = 14.0;
I2 = 15.0;
spikes1 = lifspikes(trials, I1, tmax, D);
spikes2 = lifspikes(trials, I2, tmax, D);
subplot(1, 2, 1);
tmin = 10.0;
spikeraster(spikes1, tmin, tmin+2.0);
title(sprintf('I_1=%g', I1))
subplot(1, 2, 2);
spikeraster(spikes2, tmin, tmin+2.0);
title(sprintf('I_2=%g', I2))
%savefigpdf(gcf(), 'spikeraster.pdf')
%% spike count histograms:
Ts = [0.01 0.1 0.5 1.0];
cmax = 100;
figure()
for k = 1:length(Ts)
T = Ts(k);
[c1, b1] = counthist(spikes1, 0.0, tmax, T, cmax);
[c2, b2] = counthist(spikes2, 0.0, tmax, T, cmax);
subplot(2, 2, k)
bar(b1, c1, 'r');
hold on;
bar(b2, c2, 'b');
xlim([0 cmax])
title(sprintf('T=%gms', 1000.0*T))
hold off;
end
%% discrimination measure:
T = 0.1;
cmax = 20;
[d, thresholds, true1s, false1s, true2s, false2s, pratio] = discriminability(spikes1, spikes2, tmax, T, cmax);
figure()
subplot(1, 3, 1);
plot(thresholds, true1s, 'b');
hold on;
plot(thresholds, true2s, 'b');
plot(thresholds, false1s, 'r');
plot(thresholds, false2s, 'r');
hold off;
% Ratio:
subplot(1, 3, 2);
fprintf('discriminability = %g\n', d);
plot(thresholds, pratio);
% ROC:
subplot(1, 3, 3);
plot(false2s, true1s);
%% discriminability:
Ts = 0.01:0.01:1.0;
cmax = 100;
ds = zeros(length(Ts), 1)
for k = 1:length(Ts)
T = Ts(k);
[c1, b1] = counthist(spikes1, 0.0, tmax, T, cmax);
[c2, b2] = counthist(spikes2, 0.0, tmax, T, cmax);
[d, thresholds, true1s, false1s, true2s, false2s, pratio] = discriminability(spikes1, spikes2, tmax, T, cmax);
ds(k) = d;
end
figure()
plot(Ts, ds)

View File

@ -0,0 +1,39 @@
function spikes = lifspikes(trials, input, tmaxdt, D)
% Generate spike times of a leaky integrate-and-fire neuron.
% trials: the number of trials to be generated
% input: the stimulus either as a single value or as a vector
% tmaxdt: in case of a single value stimulus the duration of a trial
% in case of a vector as a stimulus the time step
% D: the strength of additive white noise
tau = 0.01;
if nargin < 4
D = 1e0;
end
vreset = 0.0;
vthresh = 10.0;
dt = 5e-5;
if max(size(input)) == 1
input = input * ones(ceil(tmaxdt/dt), 1);
else
dt = tmaxdt;
end
spikes = cell(trials, 1);
for k=1:trials
times = [];
j = 1;
v = vreset + (vthresh-vreset)*rand(1);
noise = sqrt(2.0*D)*randn(length(input), 1)/sqrt(dt);
for i=1:length(noise)
v = v + (- v + noise(i) + input(i))*dt/tau;
if v >= vthresh
v = vreset;
spiketime = i*dt;
times(j) = spiketime;
j = j + 1;
end
end
spikes{k} = times;
end
end

View File

@ -0,0 +1,28 @@
function savefigpdf(fig, name, width, height)
% Saves figure fig in pdf file name.pdf with appropriately set page size
% and fonts
% default width:
if nargin < 3
width = 11.7;
end
% default height:
if nargin < 4
height = 9.0;
end
% paper:
set(fig, 'PaperUnits', 'centimeters');
set(fig, 'PaperSize', [width height]);
set(fig, 'PaperPosition', [0.0 0.0 width height]);
set(fig, 'Color', 'white')
% font:
set(findall(fig, 'type', 'axes'), 'FontSize', 12)
set(findall(fig, 'type', 'text'), 'FontSize', 12)
% save:
saveas(fig, name, 'pdf')
end

View File

@ -0,0 +1,30 @@
function spikeraster(spikes, tmin, tmax)
% Display a spike raster of the spike times given in spikes.
%
% spikeraster(spikes, tmax)
% spikes: a cell array of vectors of spike times in seconds
% tmin: plot spike raster starting at tmin seconds
% tmax: plot spike raster upto tmax seconds
ntrials = length(spikes);
for k = 1:ntrials
times = spikes{k};
times = times((times>=tmin) & (times<=tmax));
if tmax < 1.5
times = 1000.0*times; % conversion to ms
end
for i = 1:length( times )
line([times(i) times(i)],[k-0.4 k+0.4], 'Color', 'k');
end
end
if (tmax-tmin) < 1.5
xlabel('Time [ms]');
xlim([1000.0*tmin 1000.0*tmax]);
else
xlabel('Time [s]');
xlim([tmin tmax]);
end
ylabel('Trials');
ylim([0.3 ntrials+0.7 ]);
end

File diff suppressed because one or more lines are too long

15
slides.mk Normal file
View File

@ -0,0 +1,15 @@
# slides:
slides : $(BASENAME)-slides.pdf
$(BASENAME)-slides.pdf : $(BASENAME)-slides.tex $(PYPDFFILES) $(GPTTEXFILES)
pdflatex -interaction=scrollmode $< | tee /dev/stderr | fgrep -q "Rerun to get cross-references right" && pdflatex -interaction=scrollmode $< || true
watchslides :
while true; do ! make -q slides && make slides; sleep 0.5; done
cleanslides :
rm -f *~
rm -f $(BASENAME)-slides.aux $(BASENAME)-slides.log $(BASENAME)-slides.out $(BASENAME)-slides.idx
cleanallslides : cleanslides
rm -f $(BASENAME)-slides.pdf

View File

@ -0,0 +1,34 @@
import numpy as np
import matplotlib.pyplot as plt
# first passage time of random walker
# the average of the first passage time grows
# the more experiments you average (n)
# eventually diverges to infinity
# see https://www.math.ucdavis.edu/~tracy/courses/math135A/UsefullCourseMaterial/firstPassage.pdf last paragraph
p = 0.5
thresh = 1.0
rep = 1
rn = 10000
n = 2
meantimes = np.zeros(rep)
for i in range(rep):
times = np.zeros(n)
for k in range(n):
xx = 0.0
xxinx = 0
while True:
r = np.random.rand(rn)
x = np.zeros(len(r))
x[r < p] = 1.0
x[r >= p]= -1.0
y = xx + np.cumsum(x)
inx = np.where(y > thresh)[0]
if len(inx) > 0:
times[k] = xxinx + inx[0]
break
xx = y[-1]
xxinx += len(x)
meantimes[i] = np.mean(times)
print np.mean(meantimes)