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

This commit is contained in:
Jan Grewe 2015-10-27 17:55:37 +01:00
commit 344ae6ba59
4 changed files with 95 additions and 76 deletions

View File

@ -1,9 +1,10 @@
function spikes = pifouspikes( trials, input, tmaxdt, D, outau ) function spikes = pifouspikes( trials, input, tmaxdt, D, outau )
% Generate spike times of a perfect integrate-and-fire neuron % Generate spike times of a perfect integrate-and-fire neuron
% with Ornstein-Uhlenbeck noise (colored noise).
% trials: the number of trials to be generated % trials: the number of trials to be generated
% input: the stimulus either as a single value or as a vector % 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 % tmaxdt: in case of a single value stimulus: the duration of a trial
% in case of a vector as a stimulus the time step % in case of a vector as a stimulus: the time step
% D: the strength of additive white noise % D: the strength of additive white noise
% outau: time constant of the colored noise % outau: time constant of the colored noise

View File

@ -99,10 +99,11 @@ jan.benda@uni-tuebingen.de}
der drei Neurone miteinander vergleichen. der drei Neurone miteinander vergleichen.
\begin{parts} \begin{parts}
\part Lade die Spiketrains aus den drei Dateien. Achte darauf, dass sie verschiedene \part Lade die Spiketrains aus den drei Dateien. Achte darauf, dass sie verschiedene
Variablennamen bekommen. Variablen\-namen bekommen.
\begin{solution} \begin{solution}
\begin{lstlisting} \begin{lstlisting}
clear all clear all
% not so good:
load poisson.mat load poisson.mat
whos whos
poissonspikes = spikes; poissonspikes = spikes;
@ -111,15 +112,23 @@ jan.benda@uni-tuebingen.de}
load lifadapt.mat; load lifadapt.mat;
lifadaptspikes = spikes; lifadaptspikes = spikes;
clear spikes; clear spikes;
% better:
clear all
x = load( 'poisson.mat' );
poissonspikes = x.spikes;
x = load( pifou.mat' );
pifouspikes = x.spikes;
x = load( 'lifadapt.mat' );
lifadaptspikes = x.spikes;
\end{lstlisting} \end{lstlisting}
\end{solution} \end{solution}
\part Schreibe eine Funktion, die die Spikezeiten der ersten \part Schreibe eine Funktion, die die Spikezeiten der ersten
\code{tmax} Sekunden in einem Rasterplot visualisiert. In jeder $t_{max}$ Sekunden in einem Rasterplot visualisiert. In jeder
Zeile des Rasterplots wird ein Spiketrain dargestellt. Jeder Zeile des Rasterplots wird ein Spiketrain dargestellt. Jeder
einzelne Spike wird als senkrechte Linie zu der Zeit des einzelne Spike wird als senkrechte Linie zu der Zeit des
Auftretens des Spikes geplottet. Benutze die Funktion, um die Auftretens des Spikes geplottet. Benutze die Funktion, um die
Spikeraster der ersten 1\,s der drei Neurone zu plotten. Spikeraster der ersten 1\,s der drei Neurone nebeneinander zu plotten.
\begin{solution} \begin{solution}
\lstinputlisting{../code/spikeraster.m} \lstinputlisting{../code/spikeraster.m}
\lstinputlisting{../code/plotspikeraster.m} \lstinputlisting{../code/plotspikeraster.m}
@ -127,22 +136,22 @@ jan.benda@uni-tuebingen.de}
\colorbox{white}{\includegraphics[width=1\textwidth]{spikeraster}} \colorbox{white}{\includegraphics[width=1\textwidth]{spikeraster}}
\end{solution} \end{solution}
\part Schreibe eine Funktion, die einen einzigen Vektor mit den Interspike-Intervallen \part Schreibe eine Funktion, die einen einzigen Vektor mit den
aller Trials von Spikezeiten zur\"uckgibt. Interspikeintervallen aller Trials von Spikezeiten zur\"uckgibt.
\begin{solution} \begin{solution}
\lstinputlisting{../code/isis.m} \lstinputlisting{../code/isis.m}
\end{solution} \end{solution}
\part Schreibe eine Funktion, die ein normiertes Histogramm aus \part Schreibe eine Funktion, die ein normiertes Histogramm aus
einem Vektor von Interspike-Intervallen, gegeben in Sekunden, einem Vektor von Interspikeintervallen, gegeben in Sekunden,
berechnet und dieses mit richtiger Achsenbeschriftung plottet. Die berechnet und dieses mit richtiger Achsenbeschriftung plottet.
Interspike-Intervalle sollen dabei in Millisekunden angegeben Die Interspikeintervalle sollen dabei in Millisekunden angegeben
werden. Die Funktion soll ausserdem den Mittelwert, die Standardabweichung, werden. Die Funktion soll zus\"atzlich den Mittelwert, die
und den Variationskoeffizienten der Interspike Intervalle berechnen Standardabweichung, und den Variationskoeffizienten der
und diese im Plot mit angeben. Interspikeintervalle berechnen und diese im Plot mit angeben.
Benutze diese und die vorherige Funktion, um die Interspike-Intervall Verteilung Benutze die vorherige und diese Funktion, um die
der drei Neurone zu vergleichen. Interspikeintervall Verteilung der drei Neurone zu vergleichen.
\begin{solution} \begin{solution}
\lstinputlisting{../code/isihist.m} \lstinputlisting{../code/isihist.m}
\lstinputlisting{../code/plotisih.m} \lstinputlisting{../code/plotisih.m}
@ -150,16 +159,19 @@ jan.benda@uni-tuebingen.de}
\colorbox{white}{\includegraphics[width=1\textwidth]{isihist}} \colorbox{white}{\includegraphics[width=1\textwidth]{isihist}}
\end{solution} \end{solution}
\part Schreibe eine Funktion, die die Seriellen Korrelationen der \part Schreibe eine Funktion, die die seriellen Korrelationen der
Interspike Intervalle f\"ur lags bis zu \code{maxlag} berechnet Interspikeintervalle f\"ur Lags bis zu \code{maxlag} berechnet und
und plottet. Die Seriellen Korrelationen $\rho_k$ f\"ur lag $k$ plottet. Die Seriellen Korrelationen $\rho_k$ f\"ur Lag $k$ der
der Interspike Intervalle $T_i$ sind wie folgt definiert: Interspikeintervalle $T_i$ sind die Korrelationskoeffizienten
zwischen den Interspikeintervallen $T_i$ und den um das Lag $k$
verschobenen Intervallen $T_{i+k}$:
\[ \rho_k = \frac{\langle (T_{i+k} - \langle T \rangle)(T_i - \[ \rho_k = \frac{\langle (T_{i+k} - \langle T \rangle)(T_i -
\langle T \rangle) \rangle}{\langle (T_i - \langle T \langle T \rangle) \rangle}{\langle (T_i - \langle T
\rangle)^2\rangle} = \frac{{\rm cov}(T_{i+k}, T_i)}{{\rm \rangle)^2\rangle} = \frac{{\rm cov}(T_{i+k}, T_i)}{{\rm
var}(T_i)} = {\rm corrcoef}(T_{i+k}, T_i) \] Benutze dies Funktion, var}(T_i)} = {\rm corr}(T_{i+k}, T_i) \]
um die Interspike Intervall Korrelationen der drei Neurone zu
vergleichen. Benutze diese Funktion, um die Interspikeintervall-Korrelationen
der drei Neurone zu vergleichen.
\begin{solution} \begin{solution}
\lstinputlisting{../code/isiserialcorr.m} \lstinputlisting{../code/isiserialcorr.m}
\lstinputlisting{../code/plotserialcorr.m} \lstinputlisting{../code/plotserialcorr.m}
@ -168,11 +180,10 @@ jan.benda@uni-tuebingen.de}
\part Schreibe eine Funktion, die aus Spikezeiten \part Schreibe eine Funktion, die aus Spikezeiten
Histogramme aus der Anzahl von Spikes, die in Fenstern gegebener L\"ange $W$ Histogramme aus der Anzahl von Spikes, die in Fenstern gegebener L\"ange $W$
gez\"ahlt werden, erzeugt und plottet. Zus\"atzlich soll die Funktion gez\"ahlt werden, erzeugt und plottet.
die Poisson-Verteilung
\[ P(k) = \frac{(\lambda W)^ke^{\lambda W}}{k!} \] mit der Rate Wende diese Funktion auf die drei
$\lambda$, die aus den Daten bestimmt werden kann, mit zu dem Datens\"atze an. Probiere verschiedene Fenstergr\"o{\ss}en $W$ aus.
Histogramm hineinzeichen.
\begin{solution} \begin{solution}
\lstinputlisting{../code/counthist.m} \lstinputlisting{../code/counthist.m}
\lstinputlisting{../code/plotcounthist.m} \lstinputlisting{../code/plotcounthist.m}

View File

@ -11,11 +11,11 @@
\usepackage[left=20mm,right=20mm,top=25mm,bottom=25mm]{geometry} \usepackage[left=20mm,right=20mm,top=25mm,bottom=25mm]{geometry}
\pagestyle{headandfoot} \pagestyle{headandfoot}
\ifprintanswers \ifprintanswers
\newcommand{\stitle}{: L\"osungen} \newcommand{\stitle}{L\"osungen}
\else \else
\newcommand{\stitle}{} \newcommand{\stitle}{\"Ubung}
\fi \fi
\header{{\bfseries\large \"Ubung 6\stitle}}{{\bfseries\large Statistik}}{{\bfseries\large 27. Oktober, 2015}} \header{{\bfseries\large \stitle}}{{\bfseries\large Punktprozesse 2}}{{\bfseries\large 27. Oktober, 2015}}
\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}{}
@ -89,26 +89,27 @@ jan.benda@uni-tuebingen.de}
\begin{questions} \begin{questions}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\question \qt{Homogeneous Poisson process} \question \qt{Homogener Poisson Prozess}
We use the Poisson process to generate spike trains on which we can test and imrpove some Wir wollen den homogenen Poisson Prozess benutzen um Spikes zu generieren,
standard analysis functions. mit denen wir die Analysfunktionen des vorherigen \"Ubungszettel \"uberpr\"ufen k\"onnen.
A homogeneous Poisson process of rate $\lambda$ (measured in Hertz) is a point process Ein homogener Poisson Prozess mit der Rate $\lambda$ (measured in Hertz) ist ein Punktprozess,
where the probability of an event is independent of time $t$ and independent of previous events. bei dem die Wahrschienlichkeit eines Ereignisses unabh\"angig von der Zeit $t$ und
The probability $P$ of an event within a bin of width $\Delta t$ is unabh\"angig von vorherigen Ereignissen ist.
Die Wahrscheinlichkeit $P$ eines Ereignisses innerhalb eines Bins der Breite $\Delta t$ ist
\[ P = \lambda \cdot \Delta t \] \[ P = \lambda \cdot \Delta t \]
for sufficiently small $\Delta t$. f\"ur gen\"ugend kleine $\Delta t$.
\begin{parts} \begin{parts}
\part Write a function that generates $n$ homogeneous Poisson spike trains of a given duration $T_{max}$ \part Schreibe eine Funktion die $n$ homogene Poisson Spiketrains
with rate $\lambda$. einer gegebenen Dauer $T_{max}$ mit rate $\lambda$ erzeugt.
\begin{solution} \begin{solution}
\lstinputlisting{hompoissonspikes.m} \lstinputlisting{hompoissonspikes.m}
\end{solution} \end{solution}
\part Using this function, generate a few trials and display them in a raster plot. \part Benutze diese Funktion um einige Trials von Spikes zu erzeugen
und plotte diese als Spikeraster.
\begin{solution} \begin{solution}
\lstinputlisting{../code/spikeraster.m}
\begin{lstlisting} \begin{lstlisting}
spikes = hompoissonspikes( 10, 100.0, 0.5 ); spikes = hompoissonspikes( 10, 100.0, 0.5 );
spikeraster( spikes ) spikeraster( spikes )
@ -117,41 +118,31 @@ for sufficiently small $\Delta t$.
\colorbox{white}{\includegraphics[width=0.7\textwidth]{poissonraster100hz}} \colorbox{white}{\includegraphics[width=0.7\textwidth]{poissonraster100hz}}
\end{solution} \end{solution}
\part Write a function that extracts a single vector of interspike intervals \part Berechne Histogramme aus den Interspikeintervallen von $n$
from the spike times returned by the first function. Poisson Spiketrains mit der Rate $\lambda=100$\,Hz. Ver\"andere
\"uber die Dauer $T_{max}$ der Spiketrains und die Anzahl $n$ der
Trials die Anzahl der Intervalle und ver\"andere auch die Binbreite
des Histograms (fange mit 1\,ms an). Wieviele Interspikeintervalle
werden ben\"otigt um ein ``sch\"ones'' Histogramm zu erhalten? Wie
lange m\"usste man also von dem Neuron ableiten?
\begin{solution} \begin{solution}
\lstinputlisting{../code/isis.m} About 5000 intervals for 25 bins. This corresponds to a $5000 /
100\,\hertz = 50\,\second$ recording of a neuron firing with
100\,\hertz.
\end{solution} \end{solution}
\part Write a function that plots the interspike-interval histogram \part Vergleiche das Histogramm mit der zu erwartenden Verteilung
from a vector of interspike intervals. The function should also der Intervalle $T$ des Poisson Prozesses
compute the mean, the standard deviation, and the CV of the intervals
and display the values in the plot.
\begin{solution}
\lstinputlisting{../code/isihist.m}
\end{solution}
\part Compute histograms for Poisson spike trains with rate
$\lambda=100$\,Hz. Play around with $T_{max}$ and $n$ and the bin width
(start with 1\,ms) of the histogram.
How many
interspike intervals do you approximately need to get a ``nice''
histogram? How long do you need to record from the neuron?
\begin{solution}
About 5000 intervals for 25 bins. This corresponds to a $5000 / 100\,\hertz = 50\,\second$ recording
of a neuron firing with 100\,\hertz.
\end{solution}
\part Compare the histogram with the true distribution of intervals $T$ of the Poisson process
\[ p(T) = \lambda e^{-\lambda T} \] \[ p(T) = \lambda e^{-\lambda T} \]
for various rates $\lambda$. mit rate $\lambda$.
\begin{solution} \begin{solution}
\lstinputlisting{hompoissonisih.m} \lstinputlisting{hompoissonisih.m}
\colorbox{white}{\includegraphics[width=0.48\textwidth]{poissonisih100hz}} \colorbox{white}{\includegraphics[width=0.48\textwidth]{poissonisih100hz}}
\colorbox{white}{\includegraphics[width=0.48\textwidth]{poissonisih20hz}} \colorbox{white}{\includegraphics[width=0.48\textwidth]{poissonisih20hz}}
\end{solution} \end{solution}
\part What happens if you make the bin width of the histogram smaller than $\Delta t$ \part \extra Was pasiert mit den Histogrammen, wenn die Binbreite der Histogramme kleiner
als das bei der Erzeugung der $\Delta t$ der
used for generating the Poisson spikes? used for generating the Poisson spikes?
\begin{solution} \begin{solution}
The bins between the discretization have zero entries. Therefore The bins between the discretization have zero entries. Therefore
@ -200,3 +191,11 @@ for sufficiently small $\Delta t$.
\end{questions} \end{questions}
\end{document} \end{document}
Zus\"atzlich soll die Funktion
die Poisson-Verteilung
\[ P(k) = \frac{(\lambda W)^ke^{\lambda W}}{k!} \] mit der Rate
$\lambda$, die aus den Daten bestimmt werden kann, mit zu dem
Histogramm hineinzeichen. Hinweis: es gibt eine \code{matlab} Funktion,
die die Fakult\"at $n!$ berechnet.

View File

@ -56,19 +56,28 @@ erzeugt. Zum Beispiel:
\item Diffusions Koeffizient $D_{ISI} = \frac{\sigma_{ISI}^2}{2\mu_{ISI}^3}$. \item Diffusions Koeffizient $D_{ISI} = \frac{\sigma_{ISI}^2}{2\mu_{ISI}^3}$.
\end{itemize} \end{itemize}
\subsection{Interval return maps} \subsection{Korrelationen der Intervalle}
Scatter plot von aufeinander folgenden Intervallen $(T_{i+k}, T_i)$ getrennt durch das ``lag'' $k$. In ``return maps'' werden die um das ``Lag'' $k$ verz\"ogerten
Intervalle $T_{i+k}$ gegen die Intervalle $T_i$ geplottet. Dies macht
m\"ogliche Abh\"angigkeiten von aufeinanderfolgenden Intervallen
sichtbar.
\begin{figure}[t] \begin{figure}[t]
\includegraphics[width=1\textwidth]{returnmapexamples} \includegraphics[width=1\textwidth]{returnmapexamples}
\includegraphics[width=1\textwidth]{serialcorrexamples} \includegraphics[width=1\textwidth]{serialcorrexamples}
\caption{\label{returnmapfig}Interspike-Intervall return maps and serial correlations.} \caption{\label{returnmapfig}Interspike-Intervall return maps und
serielle Korrelationen zwischen aufeinander folgenden Intervallen
im Abstand des Lags $k$.}
\end{figure} \end{figure}
\subsection{Serielle Korrelationen der Intervalle} Solche Ab\"angigkeiten werden durch die serielle Korrelation der
Korrelationskoeffizient zwischen aufeinander folgenden Intervallen getrennt durch ``lag'' $k$: Intervalle quantifiziert. Das ist der Korrelationskoeffizient
\[ \rho_k = \frac{\langle (T_{i+k} - \langle T \rangle)(T_i - \langle T \rangle) \rangle}{\langle (T_i - \langle T \rangle)^2\rangle} = \frac{{\rm cov}(T_{i+k}, T_i)}{{\rm var}(T_i)} \] zwischen aufeinander folgenden Intervallen getrennt durch ``Lag'' $k$:
$\rho_0=1$ (Korrelation jedes Intervalls mit sich selber). \[ \rho_k = \frac{\langle (T_{i+k} - \langle T \rangle)(T_i - \langle T \rangle) \rangle}{\langle (T_i - \langle T \rangle)^2\rangle} = \frac{{\rm cov}(T_{i+k}, T_i)}{{\rm var}(T_i)}
= {\rm corr}(T_{i+k}, T_i) \]
\"Ublicherweise wird die Korrelation $\rho_k$ gegen den Lag $k$
aufgetragen (\figref{returnmapfig}). $\rho_0=1$ (Korrelation jedes
Intervalls mit sich selber).
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
@ -87,7 +96,6 @@ Statistik der Anzahl der Ereignisse $N_i$ innerhalb von Beobachtungsfenstern $i$
\item Varianz der Anzahl: $\sigma_N^2 = \langle (N - \langle N \rangle)^2 \rangle$. \item Varianz der Anzahl: $\sigma_N^2 = \langle (N - \langle N \rangle)^2 \rangle$.
\item Fano Faktor (Varianz geteilt durch Mittelwert): $F = \frac{\sigma_N^2}{\mu_N}$. \item Fano Faktor (Varianz geteilt durch Mittelwert): $F = \frac{\sigma_N^2}{\mu_N}$.
\end{itemize} \end{itemize}
Insbesondere ist die mittlere Rate der Ereignisse $r$ (``Spikes pro Zeit'', Feuerrate) gemessen in Hertz Insbesondere ist die mittlere Rate der Ereignisse $r$ (``Spikes pro Zeit'', Feuerrate) gemessen in Hertz
\[ r = \frac{\langle N \rangle}{W} \; . \] \[ r = \frac{\langle N \rangle}{W} \; . \]