Merge branch 'master' of raven.am28.uni-tuebingen.de:scientificComputing
This commit is contained in:
commit
e64cc71623
BIN
programming/exercises/STA/p-unit_spike_times.mat
Normal file
BIN
programming/exercises/STA/p-unit_spike_times.mat
Normal file
Binary file not shown.
BIN
programming/exercises/STA/p-unit_stimulus.mat
Normal file
BIN
programming/exercises/STA/p-unit_stimulus.mat
Normal file
Binary file not shown.
18
programming/exercises/STA/sta.m
Normal file
18
programming/exercises/STA/sta.m
Normal file
@ -0,0 +1,18 @@
|
|||||||
|
function [st_avg, std_sta, valid_spikes]= sta(stimulus, spike_times, count, sampling_rate)
|
||||||
|
|
||||||
|
snippets = zeros(numel(spike_times), 2*count);
|
||||||
|
valid_spikes = 1;
|
||||||
|
for i = 1:numel(spike_times)
|
||||||
|
t = spike_times(i);
|
||||||
|
index = round(t*sampling_rate);
|
||||||
|
if index <= count || (index + count) > length(stimulus)
|
||||||
|
continue
|
||||||
|
end
|
||||||
|
snippets(valid_spikes,:) = stimulus(index-count:index+count-1);
|
||||||
|
valid_spikes = valid_spikes + 1;
|
||||||
|
end
|
||||||
|
|
||||||
|
snippets(end-(end-valid_spikes):end,:) = [];
|
||||||
|
|
||||||
|
st_avg = mean(snippets, 1);
|
||||||
|
std_sta = std(snippets,[],1);
|
67
programming/exercises/STA/sta_script.m
Normal file
67
programming/exercises/STA/sta_script.m
Normal file
@ -0,0 +1,67 @@
|
|||||||
|
clear
|
||||||
|
clc
|
||||||
|
|
||||||
|
%% load data
|
||||||
|
load p-unit_spike_times.mat
|
||||||
|
load p-unit_stimulus.mat
|
||||||
|
sample_rate = 20000;
|
||||||
|
|
||||||
|
%% calculate STA
|
||||||
|
all_times = [];
|
||||||
|
for i = 1:length(spike_times)
|
||||||
|
all_times = cat(1, all_times, spike_times{i});
|
||||||
|
end
|
||||||
|
|
||||||
|
[st_average, sta_sd, num] = sta(stimulus(:,2), all_times, 1000, sample_rate);
|
||||||
|
fig = figure();
|
||||||
|
set(fig, 'PaperUnits', 'centimeters');
|
||||||
|
set(fig, 'PaperSize', [11.7 9.0]);
|
||||||
|
set(fig, 'PaperPosition',[0.0 0.0 11.7 9.0]);
|
||||||
|
set(fig,'Color', 'white')
|
||||||
|
plot(((1:length(st_average))-1000)/sample_rate, st_average)
|
||||||
|
xlabel('time [s]')
|
||||||
|
ylabel('stimulus')
|
||||||
|
|
||||||
|
%% reverse reconstruction
|
||||||
|
|
||||||
|
% make binary representation of the spike times
|
||||||
|
binary_spikes = zeros(size(stimulus, 1), length(spike_times));
|
||||||
|
estimated_stims = zeros(size(binary_spikes));
|
||||||
|
for i = 1:length(spike_times)
|
||||||
|
binary_spikes(round(spike_times{i}*sample_rate), i) = 1;
|
||||||
|
estimated_stims(:,i) = conv(binary_spikes(:,i), st_average, 'same');
|
||||||
|
end
|
||||||
|
|
||||||
|
fig = figure();
|
||||||
|
set(fig, 'PaperUnits', 'centimeters');
|
||||||
|
set(fig, 'PaperSize', [11.7 9.0]);
|
||||||
|
set(fig, 'PaperPosition',[0.0 0.0 11.7 9.0]);
|
||||||
|
set(fig,'Color', 'white')
|
||||||
|
plot(stimulus(:,1), stimulus(:,2), 'displayname','original')
|
||||||
|
hold on
|
||||||
|
plot(stimulus(:,1), mean(estimated_stims,2), 'r', 'displayname', 'reconstruction')
|
||||||
|
xlabel('time [s]')
|
||||||
|
ylabel('stimulus')
|
||||||
|
legend show
|
||||||
|
|
||||||
|
%% calculate STC
|
||||||
|
|
||||||
|
% we need to downsample the data otherwise the covariance matrixs gets too
|
||||||
|
% large 20Khz to 1kHz
|
||||||
|
|
||||||
|
% downsampled_binary = zeros(size(stimulus, 1)/20, length(spike_times));
|
||||||
|
downsampled_stim = zeros(size(stimulus,1)/20,1);
|
||||||
|
|
||||||
|
% for i = 1:length(spike_times)
|
||||||
|
% indices = round(spike_times{i}.*1000);
|
||||||
|
% indices(indices < 1) = [];
|
||||||
|
% downsampled_binary(indices, i) = 1;
|
||||||
|
% end
|
||||||
|
for i = 1:length(downsampled_stim)
|
||||||
|
start_index = (i-1) * 20 + 1;
|
||||||
|
downsampled_stim(i) = mean(stimulus(start_index:start_index+19,2));
|
||||||
|
end
|
||||||
|
|
||||||
|
[st_average, ~, ~] = sta(downsampled_stim, all_times, 50, 1000);
|
||||||
|
|
||||||
|
|
380
programming/lectures/sta_stc.tex
Normal file
380
programming/lectures/sta_stc.tex
Normal file
@ -0,0 +1,380 @@
|
|||||||
|
\documentclass{beamer}
|
||||||
|
\usepackage{xcolor}
|
||||||
|
\usepackage{listings}
|
||||||
|
\usepackage{pgf}
|
||||||
|
%\usepackage{pgf,pgfarrows,pgfnodes,pgfautomata,pgfheaps,pgfshade}
|
||||||
|
%\usepackage{multimedia}
|
||||||
|
|
||||||
|
\usepackage[english]{babel}
|
||||||
|
\usepackage{movie15}
|
||||||
|
\usepackage[latin1]{inputenc}
|
||||||
|
\usepackage{times}
|
||||||
|
\usepackage{amsmath}
|
||||||
|
\usepackage{bm}
|
||||||
|
\usepackage[T1]{fontenc}
|
||||||
|
\usepackage[scaled=.90]{helvet}
|
||||||
|
\usepackage{scalefnt}
|
||||||
|
\usepackage{tikz}
|
||||||
|
\usepackage{ textcomp }
|
||||||
|
\usepackage{soul}
|
||||||
|
\usepackage{hyperref}
|
||||||
|
\definecolor{lightblue}{rgb}{.7,.7,1.}
|
||||||
|
\definecolor{mygreen}{rgb}{0,1.,0}
|
||||||
|
|
||||||
|
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
|
||||||
|
\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}
|
||||||
|
% \frametitle{\insertsectionhead}
|
||||||
|
% \tableofcontents[currentsection,hideothersubsections]
|
||||||
|
\end{frame}
|
||||||
|
}
|
||||||
|
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%5
|
||||||
|
|
||||||
|
\setbeamertemplate{blocks}[rounded][shadow=true]
|
||||||
|
|
||||||
|
\title[]{Introduction to Scientific Computing -- \\
|
||||||
|
Cross-Correlation, Spike--Triggered--Average and Reverse Reconstruction}
|
||||||
|
\author[]{Jan Grewe\\Abteilung f\"ur Neuroethologie\\
|
||||||
|
Universit\"at T\"ubingen}
|
||||||
|
|
||||||
|
\institute[Introduction to Scientific Computing]{}
|
||||||
|
\date{2014/10/13 - 2014/11/07}
|
||||||
|
%\logo{\pgfuseimage{../../resources/UT_BM_Rot_RGB.pdf}}
|
||||||
|
|
||||||
|
\subject{Einf\"uhrung in wissenschaftliche Datenverarbeitung}
|
||||||
|
\vspace{1em}
|
||||||
|
\titlegraphic{
|
||||||
|
\includegraphics[width=0.5\linewidth]{../../resources/UT_WBMW_Rot_RGB}
|
||||||
|
}
|
||||||
|
%%%%%%%%%% configuration for code
|
||||||
|
\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
|
||||||
|
}
|
||||||
|
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
|
||||||
|
\newcommand{\mycite}[1]{
|
||||||
|
\begin{flushright}
|
||||||
|
\tiny \color{black!80} #1
|
||||||
|
\end{flushright}
|
||||||
|
}
|
||||||
|
|
||||||
|
\input{../../latex/environments.tex}
|
||||||
|
\makeatother
|
||||||
|
|
||||||
|
\begin{document}
|
||||||
|
|
||||||
|
\begin{frame}[plain]
|
||||||
|
\frametitle{}
|
||||||
|
\vspace{-1cm}
|
||||||
|
\titlepage % erzeugt Titelseite
|
||||||
|
\end{frame}
|
||||||
|
|
||||||
|
|
||||||
|
\begin{frame}[plain]
|
||||||
|
\frametitle{Rekapitulation}
|
||||||
|
\begin{enumerate}
|
||||||
|
\item PSTH\pause
|
||||||
|
\end{enumerate}
|
||||||
|
\end{frame}
|
||||||
|
|
||||||
|
|
||||||
|
\begin{frame}
|
||||||
|
\frametitle{Introduction to scientific computing}
|
||||||
|
\frametitle{Menue}
|
||||||
|
\begin{enumerate}
|
||||||
|
\item Cross-correlation
|
||||||
|
\item Estimation of the Spike--Triggered--Average -- STA
|
||||||
|
\item Reverse reconstruction
|
||||||
|
\end{enumerate}
|
||||||
|
\end{frame}
|
||||||
|
|
||||||
|
|
||||||
|
\begin{frame}[plain]
|
||||||
|
\huge{1. Recapitulation: PSTH}
|
||||||
|
\end{frame}
|
||||||
|
|
||||||
|
|
||||||
|
\begin{frame}
|
||||||
|
\frametitle{Recapitulation}
|
||||||
|
\only<1>{
|
||||||
|
\framesubtitle{Displaying the neuronal response over time - Rasterplot}
|
||||||
|
\begin{figure}
|
||||||
|
\centering
|
||||||
|
\includegraphics[width=0.375\columnwidth]{images/rasterplot}
|
||||||
|
\end{figure}
|
||||||
|
}
|
||||||
|
\only<2>{
|
||||||
|
\framesubtitle{Displaying the neuronal response over time - PSTH}
|
||||||
|
\begin{figure}
|
||||||
|
\centering
|
||||||
|
\includegraphics[width=0.5\columnwidth]{images/conv}
|
||||||
|
\end{figure}
|
||||||
|
}
|
||||||
|
\end{frame}
|
||||||
|
|
||||||
|
|
||||||
|
\begin{frame}[plain]
|
||||||
|
\huge{2. Relating stimulus and response}
|
||||||
|
\end{frame}
|
||||||
|
|
||||||
|
|
||||||
|
\begin{frame}
|
||||||
|
\frametitle{Relating stimulus and response}
|
||||||
|
\framesubtitle{How can we relate the response to the stimulus?}
|
||||||
|
\begin{figure}
|
||||||
|
\centering
|
||||||
|
\includegraphics[height=0.9\textheight]{images/conv_stim}
|
||||||
|
\end{figure}
|
||||||
|
\end{frame}
|
||||||
|
|
||||||
|
|
||||||
|
\begin{frame}
|
||||||
|
\frametitle{Relating stimulus and response}
|
||||||
|
\framesubtitle{Cross--correlation}
|
||||||
|
The cross - correlation for (time) discrete data is defined as:
|
||||||
|
\begin{equation}
|
||||||
|
Q_{rs}(\tau) = \frac{1}{2\cdot N-|\tau|} \displaystyle\sum^{T}_{\tau=-T}{r(t) \cdot s(t-\tau)} \cdot \Delta \tau
|
||||||
|
\end{equation}
|
||||||
|
\pause
|
||||||
|
Don't worry, this is already implemented in MATLAB!
|
||||||
|
\end{frame}
|
||||||
|
|
||||||
|
|
||||||
|
\begin{frame}[fragile]
|
||||||
|
\frametitle{Relating stimulus and response}
|
||||||
|
\framesubtitle{Cross--correlation}
|
||||||
|
\footnotesize
|
||||||
|
\begin{lstlisting}
|
||||||
|
[c,lags] = xcorr(stimulus, response, 10000, 'unbiased');
|
||||||
|
fig = figure();
|
||||||
|
set(fig, 'PaperUnits', 'centimeters');
|
||||||
|
set(fig, 'PaperSize', [11.7 9.0]);
|
||||||
|
set(fig, 'PaperPosition',[0.0 0.0 11.7 9.0]);
|
||||||
|
set(fig,'Color', 'white')
|
||||||
|
plot(lags/sample_rate, c)
|
||||||
|
xlabel('lag [s]')
|
||||||
|
ylabel('correlation')
|
||||||
|
\end{lstlisting}
|
||||||
|
\end{frame}
|
||||||
|
|
||||||
|
|
||||||
|
\begin{frame}
|
||||||
|
\frametitle{Relating stimulus and response}
|
||||||
|
\framesubtitle{Cross--correlation}
|
||||||
|
\begin{figure}
|
||||||
|
\includegraphics[width=0.5\linewidth]{images/correlation}
|
||||||
|
\end{figure}
|
||||||
|
\end{frame}
|
||||||
|
|
||||||
|
|
||||||
|
|
||||||
|
\begin{frame}[fragile]
|
||||||
|
\frametitle{Relating stimulus and response}
|
||||||
|
\framesubtitle{Cross--correlation - Exercises}
|
||||||
|
\begin{enumerate}
|
||||||
|
\item calculate the cross-correlation between two vectors of random
|
||||||
|
numbers.
|
||||||
|
\item Calculate the cross-correlation between one of these vectors
|
||||||
|
and itself (auto-correlation).
|
||||||
|
\item Calculate the cross-correlation between one vector and a
|
||||||
|
time-shifted version of itself (use \verb+circshift+ to do this).
|
||||||
|
\end{enumerate}
|
||||||
|
\textbf{Note:} Select max\_lag to be less than 10\% of the length of
|
||||||
|
your vectors!
|
||||||
|
\begin{enumerate}
|
||||||
|
\item Create the cross correlation of the p-unit data and stimulus.
|
||||||
|
\item \textbf{Note:} you have to convert the spike\_times to a PSTH!
|
||||||
|
\item Find out the position of the correlation peak.
|
||||||
|
\item What does this tell you?
|
||||||
|
\end{enumerate}
|
||||||
|
\end{frame}
|
||||||
|
|
||||||
|
|
||||||
|
\begin{frame}[plain]
|
||||||
|
\huge{2. Spike--Triggered--Average --- STA}
|
||||||
|
\end{frame}
|
||||||
|
|
||||||
|
|
||||||
|
\begin{frame}[fragile]
|
||||||
|
\frametitle{Relating stimulus and response}
|
||||||
|
\framesubtitle{Spike--Triggered--Average --- STA}
|
||||||
|
|
||||||
|
The \textbf{STA} is the average stimulus that led to a spike in the
|
||||||
|
neuron.
|
||||||
|
|
||||||
|
\begin{equation}
|
||||||
|
STA(\tau) = \left\langle \frac{1}{n} \displaystyle\sum^{i=1}_{n}{s(t_i - \tau)} \right\rangle
|
||||||
|
\end{equation}
|
||||||
|
|
||||||
|
Sadly, we need to implement this ourselves.\newline\newline\pause
|
||||||
|
\vspace{1em}
|
||||||
|
\noindent
|
||||||
|
\textbf{Algorithm:}
|
||||||
|
\begin{enumerate}
|
||||||
|
\item For each spike a snippet is taken from the respective stimulus.
|
||||||
|
\item The snippets are averaged.
|
||||||
|
\end{enumerate}
|
||||||
|
\end{frame}
|
||||||
|
|
||||||
|
|
||||||
|
\begin{frame}[fragile]
|
||||||
|
\frametitle{Relating stimulus and response}
|
||||||
|
\framesubtitle{Spike--Triggered--Average -- STA}
|
||||||
|
\vspace{-1em}
|
||||||
|
\begin{figure}
|
||||||
|
\centering
|
||||||
|
\includegraphics[width=0.65\columnwidth]{images/sta}
|
||||||
|
\end{figure}
|
||||||
|
\pause
|
||||||
|
\vspace{-0.5em}
|
||||||
|
What does the \textbf{STA} tell us?
|
||||||
|
\begin{enumerate}
|
||||||
|
\item Is there a relation between stimulus and response? \pause
|
||||||
|
\item Is there a lag between them and how large is it? \pause
|
||||||
|
\item How far in the past does a neuron encode? \pause
|
||||||
|
\item Can it see into the future?
|
||||||
|
\end{enumerate}
|
||||||
|
\end{frame}
|
||||||
|
|
||||||
|
|
||||||
|
\begin{frame}[fragile]
|
||||||
|
\frametitle{Relating stimulus and response}
|
||||||
|
\framesubtitle{STA -- Exercises}
|
||||||
|
\textbf{Exercise:}
|
||||||
|
\begin{enumerate}
|
||||||
|
\item Write a function \verb+sta(x, y, count, sample_rate)+ that
|
||||||
|
takes the stimulus (x), the response (y, as spike times), the
|
||||||
|
number (count) of sampling points it should cut out from the
|
||||||
|
stimulus and the sampling\_rate to convert from times to
|
||||||
|
indices.
|
||||||
|
\item \textbf{Beware:} sometimes the spike\_time may be too close
|
||||||
|
to the beginning or the end of the stimulus to cut out enough
|
||||||
|
data.
|
||||||
|
\item Calculate the STA for P-unt data.
|
||||||
|
\end{enumerate}
|
||||||
|
\end{frame}
|
||||||
|
|
||||||
|
|
||||||
|
\begin{frame}[fragile]
|
||||||
|
\frametitle{Relating stimulus and response}
|
||||||
|
\framesubtitle{Spike--Triggered--Average -- STA}
|
||||||
|
What does the \textbf{STA} tell us?
|
||||||
|
\begin{figure}
|
||||||
|
\centering
|
||||||
|
\includegraphics[width=0.25\columnwidth]{images/sta}
|
||||||
|
\end{figure}
|
||||||
|
\begin{enumerate}
|
||||||
|
\item Is there a relation between stimulus and response?\pause
|
||||||
|
\item Is there a lag between them and how large is it?\pause
|
||||||
|
\item How far in the past does a neuron encode?\pause
|
||||||
|
\item Can it see into the future?
|
||||||
|
\end{enumerate}
|
||||||
|
\end{frame}
|
||||||
|
|
||||||
|
|
||||||
|
\begin{frame}[plain]
|
||||||
|
\huge{3. Reverse reconstruction using the \textbf{STA}}
|
||||||
|
\end{frame}
|
||||||
|
|
||||||
|
|
||||||
|
\begin{frame}[fragile]
|
||||||
|
\frametitle{Relating stimulus and response}
|
||||||
|
\framesubtitle{Reverse reconstruction using the \textbf{STA}}
|
||||||
|
\textbf{Basic idea:}\newline
|
||||||
|
\begin{itemize}
|
||||||
|
\item The \textbf{STA} is the average stimlus that led to a spike. \pause
|
||||||
|
\item The other way round: Whenever we observe a spike, the stimulus
|
||||||
|
was likely to have looked like the \textbf{STA}.\pause
|
||||||
|
\item Thus, we should be able to reconstruct how the stimulus looked
|
||||||
|
like to evoke the observed response.
|
||||||
|
\end{itemize}
|
||||||
|
\end{frame}
|
||||||
|
|
||||||
|
|
||||||
|
\begin{frame}[fragile]
|
||||||
|
\frametitle{Relating stimulus and response}
|
||||||
|
\framesubtitle{Reverse reconstruction using the \textbf{STA}}
|
||||||
|
\textbf{Algorithm:}\newline
|
||||||
|
\begin{enumerate}
|
||||||
|
\item Estimate the \textbf{STA}. \pause
|
||||||
|
\item Create a binary representation of the spike train (vector of
|
||||||
|
\verb+zeros+ in which each \textbf{1} signals the occurence of a
|
||||||
|
spike).\pause
|
||||||
|
\item Use the STA as a Kernel (compare to the spike convolution
|
||||||
|
method to create the PSTH) and convolve (\verb+conv+) it with the
|
||||||
|
\textbf{STA}.
|
||||||
|
\end{enumerate}
|
||||||
|
\textbf{Exercise:}
|
||||||
|
\begin{enumerate}
|
||||||
|
\item Reconstruct the stimulus from the P-unit data.
|
||||||
|
\item Plot original and reconstructed stimulus into the same plot.
|
||||||
|
\item Are they similar? In which instances are they different?
|
||||||
|
\end{enumerate}
|
||||||
|
\end{frame}
|
||||||
|
|
||||||
|
\end{document}
|
||||||
|
|
||||||
|
\begin{frame}[plain]
|
||||||
|
\huge{4. Spike--Triggered--Covariance \textbf{STC}}
|
||||||
|
\end{frame}
|
||||||
|
|
||||||
|
|
||||||
|
\begin{frame}
|
||||||
|
\frametitle{Relating stimulus and response}
|
||||||
|
\framesubtitle{Spike--Triggered--Covariance \textbf{STC}}
|
||||||
|
\textbf{Problem:}
|
||||||
|
\begin{itemize}
|
||||||
|
\item The \textbf{STA} is a linear filter.
|
||||||
|
\item It can only reveal linear relationships between stimulus and response.\pause
|
||||||
|
\end{itemize}
|
||||||
|
\vspace{1em}
|
||||||
|
Further relations can be recovered using the Spike-triggered-covariance.\pause
|
||||||
|
\begin{equation}
|
||||||
|
STC = \frac{1}{N-1} \cdot \displaystyle\sum_{n=1}^{N}[\overrightarrow{s}(t_n) - STA]
|
||||||
|
[\overrightarrow{s}(t_n) - STA]^T
|
||||||
|
\end{equation}
|
||||||
|
where:\\ $N$ is the total number of spikes,
|
||||||
|
$\overrightarrow{s}(t_n)$ the stimulus snippet for the $n^{th}$
|
||||||
|
spike, and $STA$ the Spike-Triggered-Average.
|
||||||
|
\end{frame}
|
||||||
|
|
||||||
|
|
||||||
|
\begin{frame}
|
||||||
|
\frametitle{Relating stimulus and response}
|
||||||
|
\framesubtitle{Spike--Triggered--Covariance}
|
||||||
|
\begin{itemize}
|
||||||
|
\item
|
||||||
|
\end{itemize}
|
||||||
|
\end{frame}
|
||||||
|
|
||||||
|
|
||||||
|
\end{document}
|
Reference in New Issue
Block a user