defense_Yuanjun.tex

\documentclass[16pt,presentation]{beamer}
%\mode<presentation>{\usetheme{default}}
%\setbeamersize{text margin left=0.5cm} 
%\setbeamersize{text margin right=0.5cm} 

\usepackage{color}
\usepackage{rotating}
%\usepackage{movie15}
\usepackage{graphicx}
\usepackage[english]{babel}
\usepackage[latin1]{inputenc}
\usepackage{helvet}
\usepackage[T1]{fontenc}
\usepackage{units}
\usepackage{bm}
\usepackage{amssymb}
\usepackage{amsthm}
\usepackage{amsmath}
\usepackage{algorithm}
\usepackage{algorithmicx}
\usepackage{rotating}
\usepackage{hyperref}
\usepackage{verbatim}
%\usepackage{multimedia}
\usepackage{tikz}
\usepackage{bbm}
\usepackage{smartdiagram}
\usepackage{animate}
\usesmartdiagramlibrary{additions}
\usetikzlibrary{fit,positioning}
\usepackage[backend=bibtex,style=authoryear]{biblatex}
\usepackage{pifont}% http://ctan.org/pkg/pifont
\newcommand{\cmark}{\ding{51}}%
\newcommand{\xmark}{\ding{55}}%
\addbibresource{refs.bib}

\input{texdefs.tex}


%TODO
% start of the presentation
% notation in calcium imaging
% video for calcium imaing?

%\usetheme{Warsaw}
\usetheme{metropolis}

\setbeamerfont{footline}{size=\fontsize{6}{10}\selectfont}

%%%% title etc.
\author[Yuanjun Gao]{Yuanjun Gao}
\institute{\small Department of Statistics\\ Columbia University}
\title[Department of Statistics, Columbia University]{\large Statistical Machine Learning Methods for High-dimensional Neural Population Data Analysis}
\date{}

\begin{document}

\defbeamertemplate{subsection in toc}{bullets}{%
  \leavevmode
  \parbox[t]{2em}{\hspace{1em}\textbullet\hfill}%
  \parbox[t]{\dimexpr\textwidth-1em\relax}{\inserttocsubsection}\par}
\defbeamertemplate{section in toc}{sections numbered roman}{%
  \leavevmode%
  \MakeUppercase{\romannumeral\inserttocsectionnumber}.\ %
  \inserttocsection\par}

\beamertemplatetransparentcovereddynamicmedium
\setbeamertemplate{section in toc}[sections numbered roman]
\setbeamertemplate{subsection in toc}[bullets]
\setbeamercolor{background canvas}{bg=white}


%%%%%%%%%%%%%%%%%%% Title page %%%%%%%%%%%%%%%%%
\begin{frame}\vspace*{1.5cm}
\maketitle
\begin{center}
\vspace{-1cm}
\end{center}
\end{frame}


%%%%%%%%%%%%%%%%%%%%%%%% Overview %%%%%%%%%%%%%%%%%%

\begin{frame}
\frametitle{Overview}
%\begin{itemize}
%\item Brain, an \alert{amazing} organ!
%\end{itemize}
\begin{center}
\includegraphics[width = 0.45\textwidth]{./figs/overview/543px-Brain_human_sagittal_section.png}
\end{center}
%https://commons.wikimedia.org/wiki/File%3ABrain_human_sagittal_section.svg
\end{frame}




\begin{frame}
\frametitle{Overview}
\begin{center}

\begin{tabular}{ccc}

\begin{minipage}{0.5\textwidth}
\begin{tabular}[t]{c}
\visible<2->{\includegraphics[width=1\textwidth, height=0.1\textwidth,clip = true]{figs/ROI_cartoon/fig_cartoon_spike_1.pdf}}\\
\visible<2->{\includegraphics[width=1\textwidth,height=0.1\textwidth, clip = true]{figs/ROI_cartoon/fig_cartoon_spike_2.pdf}}\\
\visible<1->{\includegraphics[width=1\textwidth,height=0.1\textwidth,clip = true]{figs/ROI_cartoon/fig_cartoon_spike_3.pdf}}\\
\scriptsize{Time}
\end{tabular}
\end{minipage}
&
\visible<2->{\Large $\Leftarrow$}
&
\visible<2->{
\begin{minipage}{0.30\textwidth}
%\begin{center}
\includegraphics[width = 0.8\textwidth]{./figs/overview/monkey.jpg}
%\end{center}
\end{minipage}
}
\\
\visible<3->{
{\Large $\Downarrow$}&&\\
\begin{minipage}{0.5\textwidth}
\includegraphics[width = 1.0\textwidth]{./figs/overview/spike_train_data.png}
\end{minipage}}
&
\multicolumn{2}{l}{
\visible<4->{
\begin{minipage}{0.3\textwidth}
\begin{itemize}
\item Rich structure
\item Large data
\item \alert{Opportunities!}
\end{itemize}
\end{minipage}}}

%$\Rightarrow$
%&\\
%\includegraphics[width = 0.25\textwidth]{./figs/overview/Zebrafish_embryos.png}
\end{tabular}
%}
\end{center}
\end{frame}




\begin{frame}
\frametitle{Overview}
%{\scripsize
\begin{minipage}[t][3.5cm]{0.8\textwidth}
\begin{center}
\smartdiagramset{
  uniform color list=blue!30!white for 4 items,
  back arrow disabled=false,
  additions={
   additional item offset=0.85cm,
   additional item border color=red,
   additional arrow color=red,
   additional arrow tip=stealth,
   additional arrow line width=1pt,
   additional arrow style=]-latex,
   }
 }
\smartdiagramadd[flow diagram:horizontal]{%
 {\scriptsize Data collection}, {\scriptsize Pre-processing}, {\scriptsize Modeling}, {\scriptsize Hypothesis testing}%
 }{%
 below of module2/Part II, below of module3/Part I, below of module4/Part III%
 }
\smartdiagramconnect{{]-latex}}{additional-module1/module2}
\smartdiagramconnect{{]-latex}}{additional-module2/module3}
\smartdiagramconnect{{]-latex}}{additional-module3/module4}
\end{center}
%}
\end{minipage}
\end{frame}

\begin{frame}
\frametitle{Table of Contents}
\tableofcontents
\end{frame}


\AtBeginSection[]
{
  \begin{frame}
    \frametitle{Table of Contents}
    \tableofcontents[
    currentsection, currentsubsection
]
  \end{frame}
}

\AtBeginSubsection[]
{
  \begin{frame}
    \frametitle{Table of Contents}
    \tableofcontents[
    currentsection, currentsubsection
]
  \end{frame}
}

%%%%%%%%%%%%%%%%%%%%%%% latent variable model for neural population data %%%%%%%%
\section[]{Neural Population Data Analysis with Latent Variable Models}

\begin{frame}
\frametitle{State space models}
\begin{figure}
\centering
{\small
\begin{tikzpicture}
\tikzstyle{main}=[circle, minimum size = 11mm, thick, draw =black!80, node distance = 6mm]
\tikzstyle{connect}=[-latex, thick]
\tikzstyle{box}=[rectangle, draw=black!100]
  \node[main, fill = white!100] (z1) [] { $\vz_{t-1}$};
  \node[main] (z2) [right=of z1] {$\vz_{t}$ };
  \node[main] (z3) [right=of z2] {$\vz_{t+1}$};
  \node[main, fill = blue!20] (x1) [above=of z1] { $\vx_{t-1}$};
  \node[main, fill = blue!20] (x2) [above=of z2] {$\vx_{t}$ };
  \node[main, fill = blue!20] (x3) [above=of z3] {$\vx_{t+1}$};
  \node (z0) [left=of z1] {$\cdots$};
  \node (zT) [right=of z3] {$\cdots$};
  \path (z1) edge [connect] (z2)
        (z2) edge [connect] (z3)
        (z1) edge [connect] (x1)
        (z2) edge [connect] (x2)
        (z3) edge [connect] (x3)
        (z0) edge [connect] (z1)
        (z3) edge [connect] (zT);
\end{tikzpicture}
}
\end{figure}
\begin{itemize}
\item $\vx_t \in \mathbb{N}^n$: spike counts; $\vz_t \in \mathbb{R}^m$: latent variables
%\item $p(\vx, \vz) = p(\vz_1) \prod_{t=1}^{T-1}p(\vz_{t+1} | \vz_t) \prod_{t=1}^T p(\vx_t | \vz_t)$
%\item Joint distribution
%\[\log p(\vx, \vz) = \underbrace{\log p(\vz_1)}_{\text{Initial distribution}} + 
%\underbrace{\sum_{t=1}^{T-1}\log p(\vz_{t+1} | \vz_t)}_{\text{Transition model}} + 
%\underbrace{\sum_{t=1}^T \log p(\vx_t | \vz_t)}_{\text{Observation model}}\]
\item Joint distribution
\[p(\vx, \vz) = \underbrace{p(\vz_1)}_{\text{Initial distribution}} 
\underbrace{\prod_{t=1}^{T-1}p(\vz_{t+1} | \vz_t)}_{\text{Transition model}} 
\underbrace{\prod_{t=1}^T p(\vx_t | \vz_t)}_{\text{Observation model}}\]
%\item Interpretation of latent variables: Unobserved neuron, intention; Dynamical view of motor data %(TODO: elaborate this line)
\end{itemize}
\end{frame}


\begin{frame}
\frametitle{State space models: multiple trials}
\begin{figure}
\centering
{\footnotesize
\begin{tikzpicture}
\tikzstyle{main}=[circle, minimum size = 11mm, text width = 8mm, align=center, thick, draw =black!80, node distance = 6mm]
\tikzstyle{connect}=[-latex, thick]
\tikzstyle{box}=[rectangle, draw=black!100]
  \node[main, fill = white!100] (z1) [] { $\vz_{r(t-1)}$};
  \node[main] (z2) [right=of z1] {$\vz_{rt}$ };
  \node[main] (z3) [right=of z2] {$\vz_{r(t+1)}$};
  \node[main, fill = blue!20] (x1) [above=of z1] { $\vx_{r(t-1)}$};
  \node[main, fill = blue!20] (x2) [above=of z2] {$\vx_{rt}$ };
  \node[main, fill = blue!20] (x3) [above=of z3] {$\vx_{r(t+1)}$};
  \node (z0) [left=of z1] {$\cdots$};
  \node (zT) [right=of z3] {$\cdots$};
  \path (z1) edge [connect] (z2)
        (z2) edge [connect] (z3)
        (z1) edge [connect] (x1)
        (z2) edge [connect] (x2)
        (z3) edge [connect] (x3)
        (z0) edge [connect] (z1)
        (z3) edge [connect] (zT);
\end{tikzpicture}
}
\end{figure}
\begin{itemize}
\item $r=1,...,R$: trial number
\item $\vx_{rt} \in \mathbb{N}^n$: spike counts; $\vz_{rt} \in \mathbb{R}^m$: latent variables
%\item $p(\vx, \vz) = p(\vz_1) \prod_{t=1}^{T-1}p(\vz_{t+1} | \vz_t) \prod_{t=1}^T p(\vx_t | \vz_t)$
%\item Joint distribution
%\[\log p(\vx, \vz) = \underbrace{\log p(\vz_1)}_{\text{Initial distribution}} + 
%\underbrace{\sum_{t=1}^{T-1}\log p(\vz_{t+1} | \vz_t)}_{\text{Transition model}} + 
%\underbrace{\sum_{t=1}^T \log p(\vx_t | \vz_t)}_{\text{Observation model}}\]
\item Joint distribution
\[p(\vx, \vz) = \prod_{r=1}^R \left[ \underbrace{p(\vz_{r1})}_{\text{Initial distribution}} 
\underbrace{\prod_{t=1}^{T-1}p(\vz_{r(t+1)} | \vz_{rt})}_{\text{Transition model}} 
\underbrace{\prod_{t=1}^T p(\vx_{rt} | \vz_{rt})}_{\text{Observation model}} \right]\]
%\item Common input; Dynamical view of motor data (TODO: elaborate this line)
\end{itemize}
\end{frame}



\begin{frame}
\frametitle{Common parameterization and our extensions}
\begin{itemize}
\item Common assumptions for latent dynamics: Gaussian linear dynamical system (LDS)
 \[\begin{split}
 \vz_{r1} &\sim \N(\mu_1, Q_1)\\
 \vz_{r(t+1)} | \vz_{rt} &\sim \N(A \vz_{rt}, Q)
  \end{split}\]
\item Common observation models:
 \[\vx_{rt} | \vz_{rt} \sim \underbrace{\underbrace{\N(C \vz_{rt} + d, \Sigma)}_{\text{model mismatch}}
 \text{ or }
 \underbrace{\text{Poisson}\left(\exp(C \vz_{rt} + d)\right)}_{\text{equal dispersion}}}_\text{stringent assumptions}\]
\item Our extensions for observation model:
\begin{itemize}
\item Generalized count distribution (GCLDS) \parencite{Gao2015}
\item Flexible nonlinear observation (fLDS) \parencite{gao2016linear}
\end{itemize}
\end{itemize}
\end{frame}


%%%%%%%%%%%%%%%%%%%%%%% GCLDS %%%%%%%%
\subsection[]{Generalized count linear dynamical system}
\begin{frame}
\frametitle{Motivation}
\begin{itemize}
\item Doubly stochastic Poisson model implies \alert{overdispersion}
%\[ \vx \sim \text{Poisson}(f(\vz))  , \vz \sim p(\vz) \rightarrow \text{var}(\vx) \geq E(\vx)\]
\[\left. \begin{array}{ll} \vz &\sim p(\vz) \\ \vx | \vz &\sim \text{Poisson}(f(\vz)) \end{array} \right\} \Rightarrow 
\begin{array}{ll} \text{var}(\vx | \vz) &= E(\vx | \vz) \\ \text{var}(\vx) &\geq E(\vx)  \end{array}\]
%\item Need a more flexible distribution to separate \alert{firing rate variability} with \alert{noise variability}.
%\[\text{var}(\vx) = \underbrace{\text{var}\left(E(\vx | \vz)\right)}_\text{firing rate variability} + \underbrace{E\left(\text{var}(\vx | \vz)\right)}_\text{noise variability}\]
\end{itemize}
\end{frame}


\begin{frame}
\frametitle{Generalized count distribution family}
\begin{itemize}
\item Generalized count (GC) distribution family
\[\begin{split}
p_{\text{Poisson}}(x; \lambda) \propto& \frac{\exp\left\{\log{\lambda} \cdot x\right\}}{x!},~~~x \in \mathbb{N}\\
\Downarrow&\\
p_{\mathcal{GC}}(x; \theta, g(\cdot)) \propto& \frac{\exp(\theta \cdot x + g(x) )}{x!}, ~~~x \in \mathbb{N}
\end{split}\]
where $\theta \in \mathbb{R}$, $g(\cdot): \mathbb{N} \rightarrow \mathbb{R}$.
\item Parameterizes \alert{all} the count distributions \alert{redundantly}.
\item $\theta$ acts as ``rate'' parameter.
\item $g(\cdot)$ acts as ``shape'' parameter. %Convex/concave/linear $g(\cdot)$ implies overdispersed/underdispered/Poisson distribution.
\end{itemize}
\end{frame}


\begin{frame}
\frametitle{Generalized count distribution family}
\begin{center}
\includegraphics[width = 0.60\textwidth]{./figs/gclds/fig_var_GPoiss.pdf}
\end{center}
\begin{itemize}
\item Probability mass function for GC distribution family
\[p_{\mathcal{GC}}(x; \theta, g(\cdot)) \propto \frac{\exp(\theta \cdot x + g(x) )}{x!}, ~~~x \in \mathbb{N}\]
\end{itemize}
\end{frame}


\begin{frame}
\frametitle{Model formulation}
 \begin{itemize}
 \item Linear dynamical systems with generalized count observation 
 \[\begin{split}
 \vz_{r1} &\sim \N(\mu_1, Q_1)\\
 \vz_{r(t+1)} | \vz_{rt} &\sim \N(A \vz_{rt}, Q)\\
 x_{rti} | \vz_{rt} & \sim \mathcal{GC}(c_i^T \vz_{rt}, g_i(\cdot))%, i = 1,...,n
 \end{split}\]
 \item Practical considerations
 \begin{itemize}
 %\item Set $g_i(0) = 0$ without loss of generality;
 \item Set $g_i(k) = -\infty$ for $k > K$ to facilitate computation.
 \item Ridge penalty on the $2^{\text{nd}}$ difference of $g_i(\cdot)$ to avoid overfitting; $\text{penalty} = \lambda \sum_{k=1}^{K-1} (g_i(k-1) - 2 g_i(k) + g_i(k+1))^2$.
 %\item Set $g_i(0) = 0$ without loss of generality.
 \end{itemize}
 \end{itemize}
\end{frame}


\begin{frame}
\frametitle{Variational Bayes Expectation Maximization (VBEM)}
\begin{itemize}
\item $\vx$: data, $\vz$: latent variables, $\theta$: model parameters, 
\item Often hard to compute $p_\theta(\vx) = \int p_\theta(\vx, \vz) d \vz$ and $p_\theta(\vz | \vx)$.
\item Approximate the posterior by a \alert{tractable} distribution family.
\[p_{\theta}(\vz | \vx) \approx q(\vz) \in \mathcal{Q}\]
\item Optimize a \alert{lower bound of log likelihood}, or ELBO %w.r.t. both the model parameters and the variational distribution.
\[\begin{split}
&\text{ELBO}(\theta, q) = \int \left[\log \frac{p_{\theta}(\vx, \vz)}{q(\vz)}\right] q(\vz) d\vz \leq \log p_{\theta}(\vx) \\
%&= \log p_{\theta}(\vx) - \text{KL}(q(\vz) || p_\theta(\vz | \vx)) \leq \log p_{\theta}(\vx) % \log p_\theta(\vx) =& \log \int p_{\theta}(\vx, \vz) d \vz\\ 
\end{split}\]
%\[\begin{split}
%\log p_\theta(\vx) =& \log \int p_{\theta}(\vx, \vz) d \vz\\
%\geq& \int \left[\log p_{\theta}(\vx, \vz) - \log q(\vz)\right] q(\vz) d\vz
%\end{split}\]
\end{itemize}
\end{frame}

\begin{frame}
\frametitle{Variational Bayes Expectation Maximization (VBEM)}
\begin{itemize}
\item VBEM: Optimize $\text{ELBO}(\theta, q) \leq \log p_\theta(\vx)$ iteratively %by \alert{coordinate ascent}
\begin{itemize}
\item E-step: For a fixed $\theta$, optimize $q$ %(such that $q(\vz) \approx p_{\theta}(\vz | \vx)$)
\item M-step: For a fixed $q$, optimize $\theta$
%\item Iterate between E-step and M-step until convergence %, each iteration gives a higher ELBO and hopefully a higher log likelihood
\end{itemize}
\item VBEM for GCLDS 
\begin{itemize}
\item We set $q$ to be multivariate Gaussian
\item We derive a looser but tractable ELBO%The vanilla ELBO is still intractable, we use Jensen's inequality to get a tractable but looser bound. %(looser bound $\Rightarrow$ worse performance?)
%\item Fast initialization by Laplace approximation in .
\item E-step: fast Laplace approximation initialization + dual optimization \parencite{emtiyaz2013fast}
\item M-step: convex optimization + analytical solution
\end{itemize}
\end{itemize}
\end{frame}


\begin{frame}
\frametitle{Experiments}
\begin{itemize}
\item For both simulated and real dataset, we compare GCLDS with PLDS (Poisson observation model)
%\XSolid \Checkmark 
%\item On both simulated and real data
%GCLDS captures both mean and variance of the data well, PLDS capture 
\begin{center}
\begin{tabular}{ cccc } 
 \hline
  & Mean & \alert{Variance} & \alert{Likelihood} \\
 \hline
 PLDS & \cmark & \xmark& \xmark \\ 
 \alert{GCLDS} &\cmark &\cmark &\cmark \\ 
 \hline
\end{tabular}
\end{center}
\end{itemize}
\end{frame}


\begin{frame}
\frametitle{Real data analysis: data}
\begin{tabular}{ cc } 
{\small Data} & {\small Variance and mean of spike counts}\\
%\includegraphics[width = 0.40\textwidth]{./figs/gclds/fig_monkey.png}&
\includegraphics[width = 0.40\textwidth]{./figs/overview/monkey.jpg}&
\includegraphics[width = 0.45\textwidth]{./figs/gclds/fig_var_obs_Move_seq14.pdf}
\end{tabular}
\begin{itemize}
\item Center-out reaching experiments
\item Multi-electrode array recording
\item Strong \alert{under-dispersion}
\end{itemize}
\end{frame}


\begin{frame}
\frametitle{Real data analysis: algorithms}
\begin{itemize}
\item Main algorithms to be compared 
\begin{itemize}
\item \alert{PLDS}: Poisson observation \parencite{Macke2015}
\item \alert{GCLDS-full}: Generalized count observation, individual $g(\cdot)$ across neurons
\end{itemize}
\item Two control cases for GCLDS
\begin{itemize}
\item \alert{GCLDS-linear}: truncated linear $g(\cdot)$ (truncated Poisson)
\item \alert{GCLDS-simple}: $g(\cdot)$ shared across neurons (up to a linear function)
\end{itemize}
\end{itemize}
\end{frame}


\begin{frame}
\frametitle{Real data analysis: single neuron fit}
\begin{tabular}[t]{ccc}
{\small Fitted $g(\cdot)$} & {\small Fitted mean} & {\small Fitted variance}\\
		\raisebox{-0.85\totalheight}{\includegraphics[scale=0.5,clip = true]{figs/gclds/fig_962346_seq14_xDim8_neu1_g.pdf}}&
		\raisebox{-0.85\totalheight}{\includegraphics[scale=0.5,clip = true]{figs/gclds/fig_962346_seq14_xDim8_neu1_mean.pdf}}&
		\raisebox{-0.85\totalheight}{\includegraphics[scale=0.5,clip = true]{figs/gclds/fig_962346_seq14_xDim8_neu1_var.pdf}}\\

		%B)&\raisebox{-0.85\totalheight}{\includegraphics[width=0.3\textwidth,clip = true]{figs/fig_gclds/fig_1625526_seq14_xDim8_neu2_g.pdf}}&
		%\raisebox{-0.85\totalheight}{\includegraphics[width=0.3\textwidth,clip = true]{figs/fig_gclds/fig_1625526_seq14_xDim8_neu2_mean.pdf}}&
		%\raisebox{-0.85\totalheight}{\includegraphics[width=0.3\textwidth,clip = true]{figs/fig_gclds/fig_1625526_seq14_xDim8_neu2_var.pdf}}\\
		\raisebox{-0.85\totalheight}{\includegraphics[scale=0.5,clip = true]{figs/gclds/fig_962346_seq14_xDim8_neu3_g.pdf}}&
		\raisebox{-0.85\totalheight}{\includegraphics[scale=0.5,clip = true]{figs/gclds/fig_962346_seq14_xDim8_neu3_mean.pdf}}&
		\raisebox{-0.85\totalheight}{\includegraphics[scale=0.5,clip = true]{figs/gclds/fig_962346_seq14_xDim8_neu3_var.pdf}}\\
				%&{\graphFont Dates}&\\
\end{tabular}
\end{frame}

\begin{frame}
\frametitle{Real data analysis: population fit}
%\begin{centering}
\begin{itemize}
\item Leave-one-neuron-out prediction
%\begin{itemize}
%\item Separate trials into training and testing
%\item Use training data to learn parameters
%\item For test data, drop one neuron and use other neurons to predict its firing rate
%\end{itemize}
\end{itemize}
\begin{center}
\begin{tabular}[t]{cc}
{\small MSE reduction} & {\small NLL reduction} \\
\includegraphics[scale=0.5,clip = true]{figs/gclds/fig_MSE_band_George_Move_NULL.pdf}&
\includegraphics[scale=0.5,clip = true]{figs/gclds/fig_likelihood_band_George_Move_NULL.pdf}%&
%\includegraphics[scale=0.43,clip = true]{figs/gclds/fig_var_JOB962346_seq14_xDim8_plot.pdf}
\end{tabular}
\end{center}
%\end{centering}
\end{frame}



\begin{frame}
\frametitle{Conclusion and discussion}
\begin{itemize}
\item Summary
\begin{itemize}
\item Incorporated generalized count family into state space models.
\item Developed VBEM algorithm.
\item Observed superior fitted results on real neural data.
\end{itemize}
%\item Future work
%\begin{itemize}
%\item Time-varying dispersion structure.
%\item Hierarchical model that share information of $g(\cdot)$ across neurons.
%\item Generative models for under-dispersion.
%\end{itemize}
\item \alert{Gao Y}, Buesing L, Shenoy KV, Cunningham JP (2015) High-dimensional neural spike train analysis with generalized count linear dynamical systems. NIPS 2015.
\end{itemize}
\end{frame}



%%%%%%%%%%%%%%%%%%%%%%% fLDS %%%%%%%%
\subsection[]{Linear dynamical neural population models through nonlinear embeddings}

\begin{frame}
\frametitle{Motivation}
% TODO: Explain this data better? Or not using this figure?
%\begin{tabular}{cl}
%\parbox{0.6\textwidth}{
\begin{itemize}
\item Neural activities lie in a low-dimensional \alert{nonlinear manifold} rather than a \alert{linear subspace}
\item Flexible observation model makes the state space model more expressive
\end{itemize}
%}
%&
%\hspace{-1cm}
%\begin{tabular}{c}
%\includegraphics[scale=0.5,clip = true]{figs/flds/fig_V1_single_3concat_onlyobservation.pdf}
%\end{tabular}
%\end{tabular}
\end{frame}


\begin{frame}
\frametitle{Model formulation: fLDS}
\begin{itemize}
 \item Linear dynamical systems with \alert{nonlinear link} and count observation
 \[\begin{split}
 \vz_{r1} \sim& \N(\mu_1, Q_1)\\
 \vz_{r(t+1)} | \vz_{rt} \sim& \N(A \vz_{rt}, Q)\\
 x_{rti} \sim& \text{Poisson}(\alert{f}_i(\vz_{rt})) \text{ (PfLDS) } \\
 &\text{   or }\mathcal{GC}(\alert{f}_i(\vz_{rt}), g_i(\cdot)) \text{ (GCfLDS)}%, i = 1,...,n
 \end{split}\]
 where $f_i$ is a nonlinear function parameterized by a neural network.
%\item Linear latent dynamics: simple, tractable, interpretable
%\item Nonlinear observation: flexible
\end{itemize}
\end{frame}

\begin{frame}
\frametitle{Inference algorithm: AEVB (high level idea)}
%\begin{itemize}
%\item Auto-encoding Variational Bayes (AEVB)
\begin{itemize}
\item Auto-encoding Variational Bayes (AEVB)
\item Learn a \alert{mapping (recognition model)} from data to the approximate posterior distribution of latent variable. 
\item Jointly optimize the generative model parameters and recognition model parameters.
\item Naturally incorporate stochastic optimization to handle large datasets.
\item Tractable for a large class of graphical models
\end{itemize}
%\end{itemize}
\end{frame}

\begin{frame}
\frametitle{Inference algorithm: AEVB (algorithm)}
\begin{itemize}
\item Decompose ELBO by trials
\[\begin{split}
\text{ELBO}(\theta, q) %=& \int \left[\log p_{\theta}(\vx, \vz) - \log q(\vz)\right] q(\vz) d\vz\\
=& \sum_{r=1}^{R} \int \left[\log \frac{p_{\theta}(\vx_r, \vz_r)}{q(\vz_r)}\right] q(\vz_r) d\vz_r 
\end{split}\]
\item Map data $\vx_r$ to $q(\vz_r)$ by a parameterized function
\[q(\vz_r) = q_\phi(\vz_r | \vx_r) = \N\left(\mu_\phi(\vx_r), \Sigma_\phi(\vx_r)\right)\]
\item Learn both $\theta$ and $\phi$ by \alert{stochastic} gradient descent on ELBO
\[\begin{split}
\grad\text{ELBO}(\theta, q_\phi) \approx& R \times \grad \int \left[\log \frac{p_{\theta}(\vx_r, \vz_r)}{q_\phi(\vz_r | \vx_r)}\right] q_\phi(\vz_r| \vx_r) d\vz_r\\
\approx& R \times \text{an unbiased estimator of gradient} %\\
%&\text{{\small \alert{(Reparameterization trick)}}}
\end{split}\]
%\item Do stochastic optimization with gradient of a single trial
\end{itemize}
\end{frame}

\begin{frame}
\frametitle{Structure of the recognition model}
\begin{itemize}
\item Generative model:
\[p_\theta(\vz_r | \vx_r) \propto \underbrace{p_\theta(\vz_{r1})
\prod_{t=1}^{T-1}p_\theta(\vz_{t(t+1)} | \vz_{rt})}_{Gaussian} 
\prod_{t=1}^T \underbrace{p_\theta(\vx_{rt} | \vz_{rt})}_{\alert{Complicated}}
\]
%No analytical form, hard to sample.
\item Recognition model, product-of-Gaussian form:
\[q_\phi(\vz_r | \vx_r) \propto \underbrace{q_\phi(\vz_{r1}) 
\prod_{t=1}^{T-1}q_\phi(\vz_{r(t+1)} | \vz_{rt}) }_{Gaussian}
\prod_{t=1}^T \underbrace{q_\phi(\vz_{rt} | \vx_{rt})}_{\alert{Gaussian}}
\]
%Approximates a complicated factor with a Gaussian factor dependent on the data in a complicated way.
\item Multivariate Gaussian distribution with block tri-diagonal precision matrix. %Maintaining the Markovian structure.
\end{itemize}
\end{frame}




\begin{frame}
\frametitle{Experiments}
\begin{center}
\begin{tabular}{ ccccc } 
 \hline
  & Mean & Variance & Likelihood & \alert{Concise representation}\\
 \hline
 PLDS & \cmark & \xmark& \xmark &  \xmark  \\ 
 GCLDS &\cmark &\cmark &\cmark &  \xmark  \\ 
 \alert{PfLDS} & \cmark & \xmark& \xmark &   \cmark\\ 
 \alert{GCfLDS} &\cmark &\cmark &\cmark &  \cmark \\ 
 \hline
\end{tabular}
\end{center}
\end{frame}

\begin{frame}
\frametitle{Real data analysis: primate visual cortex}
%\hspace{-1cm}
\begin{tabular}{cc}
\small{Firing rate} & \\
\begin{tabular}{c}
\includegraphics[scale=0.45,clip = true]{figs/flds/fig_V1_single_3concat_onlyobservation.pdf} 
\end{tabular}
&
\begin{tabular}{c}
%PLDS \\
%\includegraphics[scale=0.5,clip = true]{figs/flds/fig_V1_mat_trajectory_3concat.pdf}\\
%PfLDS \\
%\includegraphics[scale=0.5,clip = true]{figs/flds/fig_V1_py_trajectory_3concat.pdf}
\end{tabular}
\end{tabular}
\end{frame}



\begin{frame}
\frametitle{Real data analysis: primate visual cortex}
%\hspace{-1.2cm}
\begin{tabular}{ccc}
{\small \visible<1->Firing rate} & \visible<2->{\hspace{-0.6cm}\small{Latent projection}} & \visible<3->{\hspace{-0.6cm}\small{1-step-ahead prediction}}\\
%\hline
\begin{tabular}{c}
\visible<1->{\includegraphics[scale=0.45,clip = true]{figs/flds/fig_V1_single_3concat.pdf}}
%\visible<2->{\includegraphics[scale=0.5,clip = true]{figs/flds/fig_V1_single_3concat.pdf} }
\end{tabular}
&
\hspace{-0.6cm}
\begin{tabular}{c}
\visible<2->{\scriptsize{PLDS} \\
\includegraphics[scale=0.4,clip = true]{figs/flds/fig_V1_mat_trajectory_3concat.pdf}\\
\scriptsize{PfLDS} \\
\includegraphics[scale=0.4,clip = true]{figs/flds/fig_V1_py_trajectory_3concat.pdf}}
\end{tabular}
&
\hspace{-0.6cm}
\visible<3->{\begin{tabular}{c}
\scriptsize{MSE reduction}\\
\includegraphics[scale=0.40,clip = true, trim = 0cm 0cm 0cm 0.2cm]{figs/flds/fig_V1_MSE_3concat.pdf}\\
\scriptsize{NLL reduction}\\
\includegraphics[scale=0.40,clip = true, trim = 0cm 0cm 0cm 0.2cm]{figs/flds/fig_V1_NLL_3concat.pdf}\\
\end{tabular}}
\end{tabular}
\end{frame}


\begin{frame}
\frametitle{Real data analysis: Primate motor cortex}
\begin{tabular}{ cc } 
{\small Data} & {\small Reaching trajectory} \\
%\includegraphics[width = 0.45\textwidth]{./figs/gclds/fig_monkey.png}&
\includegraphics[width = 0.45\textwidth]{./figs/overview/monkey.jpg}&
\includegraphics[width = 0.40\textwidth]{./figs/flds/fig_GeorgeMove_trajectory.pdf}
\end{tabular}
\end{frame}


\begin{frame}
\frametitle{Real data analysis: Primate motor cortex}
\begin{itemize}
\item Latent projection with $2$ latent dimensions
\end{itemize}
\begin{tabular}{ ccc } 
{\small Reaching trajectory} & {\small PLDS} & {\small PfLDS} \\
\includegraphics[width = 0.30\textwidth]{./figs/flds/fig_GeorgeMove_trajectory.pdf}&
\includegraphics[width = 0.30\textwidth]{./figs/flds/fig_GeorgeMove_2dvisual_mat_PLDS.pdf}&
\includegraphics[width = 0.30\textwidth]{./figs/flds/fig_GeorgeMove_2dvisual_py_PLDS.pdf}
\end{tabular}
\begin{itemize}
\item Recovering the behavior structure from neural data in a \alert{unsupervised} fashion.
\item \alert{Concise} and \alert{informative} latent dimensions.
\end{itemize}
\end{frame}

\begin{frame}
\frametitle{Real data analysis: Primate motor cortex}
\begin{itemize}
\item One-step-ahead predictive performance
\end{itemize}
\begin{center}
\begin{tabular}{ cc} 
{\small MSE reduction} & {\small NLL reduction}\\
\includegraphics[scale = 0.6]{./figs/flds/fig_GeorgeMove_MSE_GC.pdf}&
\includegraphics[scale = 0.6]{./figs/flds/fig_GeorgeMove_NLL_GC.pdf}
\end{tabular}
\end{center}
\end{frame}


\begin{frame}
\frametitle{Conclusion and discussion}
\begin{itemize}
\item Summary
\begin{itemize}
\item Incorporated nonlinear observation into state space models.
\item Developed AEVB algorithm (flexible and scalable).
\item Obtained concise latent representations.
\end{itemize}
%\item Future work
%\begin{itemize}
%\item Better stochastic optimization scheme.
%\item Interpretable nonlinearity.
%\item Application on more complex datasets.
%\end{itemize}
\item \alert{Gao Y}*, Archer E*, Paninski L, Cunningham JP (2016) Linear dynamical neural population models through nonlinear embeddings. NIPS 2016. (* = equal contribution)
\end{itemize}
\end{frame}


%%%%%%%%%%%%%%%%%%%%%% ROI detection %%%%%%%%%%%%%%%%%%%%%%%%%
\section[]{Region of Interest Detection for Calcium Imaging Data}

\begin{frame}
\frametitle{Introduction: calcium imaging data}
\begin{itemize}
\item Basic principle: the \alert{spiking} activity of a neuron induces a transient increase in \alert{calcium concentration}, which can be indirectly observed by recording the \alert{fluorescent properties} of certain calcium indicators.
%\item Allows simultaneous recording from hundreds of thousands of neurons.
\end{itemize}
\begin{center}
{\scriptsize
\begin{tabular}[t]{c}
  Spike train (cartoon) \\
\includegraphics[scale=0.4,clip = true]{figs/ROI_cartoon/fig_cartoon_spike_1.pdf}\\
Calcium trace (cartoon) \\
\includegraphics[scale=0.4,clip = true]{figs/ROI_cartoon/fig_cartoon_calcium_1.pdf}
\end{tabular}
}
\end{center}
\end{frame}


\begin{frame}
\frametitle{Introduction: calcium imaging data}
\begin{itemize}
\item Calcium imaging enables \alert{simultaneous} recording of \alert{many} neurons.
\end{itemize}
%{\scriptsize
\begin{tabular}[t]{ccccc}
Neuron shape & $\times$ & Calcium trace  & $+$  Noise $=$ &Obseravtion\\
\hspace{-1cm}\begin{minipage}{0.2\textwidth}
\begin{tabular}[t]{ccc}
\includegraphics[scale=0.2,clip = true]{figs/ROI_cartoon/fig_cartoon_neuron_1.pdf}&
\hspace{-0.3cm}\includegraphics[scale=0.2,clip = true]{figs/ROI_cartoon/fig_cartoon_neuron_2.pdf}&
\hspace{-0.3cm}\includegraphics[scale=0.2,clip = true]{figs/ROI_cartoon/fig_cartoon_neuron_3.pdf}
\end{tabular}
\end{minipage}
&&
\begin{minipage}{0.2\textwidth}
%\begin{tabular}[t]{c}
\includegraphics[scale=0.18,clip = true]{figs/ROI_cartoon/fig_cartoon_calcium_1.pdf}\\
\includegraphics[scale=0.18,clip = true]{figs/ROI_cartoon/fig_cartoon_calcium_2.pdf}\\
\includegraphics[scale=0.18,clip = true]{figs/ROI_cartoon/fig_cartoon_calcium_3.pdf}
%\end{tabular}
\end{minipage}
&&
\begin{minipage}{0.06\textwidth}
\begin{center}
\animategraphics[width=1.3cm, height=2.6cm, loop]{12}{figs/sim_vid/fig_cartoon_spike-}{1}{50}
%\movie[width=1.3cm, height=2.6cm, borderwidth=0pt, autostart, loop]{}{video_cal_simulation_3neurons_autumn.avi}
  \end{center}
 \end{minipage}
\end{tabular}
\begin{itemize}
\item \alert{Goal}: recover the neuron shape and calcium trace given the observation.
\end{itemize}
%}
\end{frame}


\begin{comment}
\begin{frame}
\frametitle{Introduction: calcium imaging data}
\begin{itemize}
\item Calcium imaging enables \alert{simultaneous} recording of \alert{many} neurons.
\end{itemize}
%{\scriptsize
\begin{tabular}[t]{ccccc}
Neuron shape & $\times$ & Calcium trace  & $+$  Noise $=$ &Obseravtion\\
\hspace{-1cm}\begin{minipage}{0.2\textwidth}
\begin{tabular}[t]{ccc}
\includegraphics[scale=0.2,clip = true]{figs/ROI_cartoon/fig_cartoon_neuron_1.pdf}&
\hspace{-0.3cm}\includegraphics[scale=0.2,clip = true]{figs/ROI_cartoon/fig_cartoon_neuron_2.pdf}&
\hspace{-0.3cm}\includegraphics[scale=0.2,clip = true]{figs/ROI_cartoon/fig_cartoon_neuron_3.pdf}
\end{tabular}
\end{minipage}
&&
\begin{minipage}{0.2\textwidth}
%\begin{tabular}[t]{c}
\includegraphics[scale=0.18,clip = true]{figs/ROI_cartoon/fig_cartoon_calcium_1.pdf}\\
\includegraphics[scale=0.18,clip = true]{figs/ROI_cartoon/fig_cartoon_calcium_2.pdf}\\
\includegraphics[scale=0.18,clip = true]{figs/ROI_cartoon/fig_cartoon_calcium_3.pdf}
%\end{tabular}
\end{minipage}
&&
\begin{minipage}{0.06\textwidth}
\begin{center}
 \includegraphics[scale=0.2,clip = true]{figs/ROI_cartoon/fig_cartoon_frame_297.pdf}
  \end{center}
 \end{minipage}
\end{tabular}
\begin{itemize}
\item \alert{Goal}: recover the neuron shape and calcium trace given the observation.
\end{itemize}
%}

\end{frame}
\end{comment}


\begin{frame}
\frametitle{Model formulation}
\begin{itemize}
\item  $X \in \mathbb{R}^{N \times T}$ represents the calcium imaging data, where each column is a (vectorized) frame that contains $N$ pixels
\item Decompose $X$ into a product of $K$ \alert{spatial component} and \alert{temporal component}% (neural acitivities)
\[X = D A^T + \text{noise}\]
\begin{itemize}
\item $D = [D_1,...,D_K] \in \mathbb{R}^{N \times K}$ represents the neuron shapes
\item $A = [A_1,...,A_K] \in \mathbb{R}^{T \times K}$ is the neural activities%, $K$ is the number of neurons.
\end{itemize}
\item Further exploit structure of the components (localized neuron shapes)
\end{itemize}
\end{frame}

\begin{frame}
\frametitle{Model formulation: objective}
\begin{itemize}
\item Structured matrix factorization
\[
%\label{equ:opt}
\begin{aligned}
& \underset{D, A}{\text{minimize}}
& & \| X - D A^T \|_2^2 + f_D(D), \\
& \text{subject to}
& & D_k \in \mathcal{D}_{w}^+; k = 1, \ldots, K,\\
& 
& & \|A_k\|_2 \leq c_k,
\end{aligned}
\]
\item $\mathcal{D}_{w}^+$: non-negative vectors whose nonzero values is within a $w \times w$ window
\item $f_D(D)$ regularizes the neuron shape (will discuss later)
\item $\|A_k\|_2 \leq c_k$ avoids degenerate solution
\end{itemize}
\end{frame}


\begin{frame}
\frametitle{Greedy algorithm}
\begin{itemize}
\item Identify ROIs one at a time, using the residual un-explained by existing signals.
\item At iteration $k$, given the current residue (unexplained by existing ROIs):
\begin{itemize}
\item \alert{Greedy identification}: Identify the location $p_k$ where the Gaussian kernel explains most of the data (across time).
\item \alert{Shape fine tuning}: Locally optimize the spatial and temporal component.
\item \alert{Residue update}: Subtract the newly identified signal.
\end{itemize}
\end{itemize}
\end{frame}


\begin{frame}
\frametitle{Shape fine tuning}
\begin{itemize}
\item Given current residue $R$, an identified center pixel $p_k$, denote $S_k$ as a $w \times w$ window centered at $p_k$
\[
\begin{aligned}
& \underset{D_k, A_k}{\text{minimize}}
& & \| R - D_k A_k^T \|^2 + f(D_k), \\%\sum_{i = 1}^3 \lambda_i f_i(D_k),\\
& \text{subject to}
& & D_{kp} \geq 0, p \in S_k,\\
&
& & D_{kp} = 0, p \notin S_k,\\
& 
& & \|A_k\|_2 \leq c_k,
\end{aligned}
\]
%\begin{itemize}
\item $f(D_k) = \sum_{i=1}^3 \lambda_i f_i(D_k)$
%where
\begin{itemize}
\item $f_1(D_k) = \sum_p \tau_{(p, p_k)} | D_{kp} |$ encourages sparsity
\item $f_2(D_k) = \sum_p (D_{kp} - G_{p_k})^2$ encourages Gaussian shape
\item $f_3(D_k) = \sum_{\text{$p_1$ and $p_2$ are neighbors}} (D_{kp_1} - D_{kp_2})^2$ encourages smoothness
\end{itemize}
\item Optimize by block coordinate descent
%\item Need to tune regularization parameters, not necessarily good since those would shrink the estimated shape
\end{itemize}
%\end{itemize}
\end{frame}




\begin{frame}
\frametitle{Real data analysis: sample patch, no shape regularization}
{\small
\begin{tabular}[t]{cc}
\parbox{0.3\textwidth}{Mean image with fitted ROI locations}& Fitted ROI shape 
\vspace{0.5cm}
\\
%\multicolumn{2}{c}{Gaussian kernel shape}  \\
\includegraphics[scale=0.35,clip = true]{figs/ROI/fig_Misha_plain_comp.pdf}&
\includegraphics[scale=0.35,clip = true]{figs/ROI/fig_Misha_plain_shape.pdf}\\
%\multicolumn{2}{c}{Ring shape} \\
%\includegraphics[scale=0.45,clip = true]{figs/fig_ROI/fig_simulation_comp_Ring.pdf}&
%\includegraphics[scale=0.45,clip = true]{figs/fig_ROI/fig_simulation_shape_Ring.pdf}\\
%\\
\end{tabular}
}
\end{frame}

\begin{frame}
\frametitle{Real data analysis: sample patch, with shape regularization}
{\small
\begin{tabular}[t]{cc}
\parbox{0.3\textwidth}{Mean image with fitted ROI locations}& Fitted ROI shape 
\vspace{0.5cm}
\\

%\multicolumn{2}{c}{Gaussian kernel shape}  \\
\includegraphics[scale=0.35,clip = true]{figs/ROI/fig_Misha_smooth_comp.pdf}&
\includegraphics[scale=0.35,clip = true]{figs/ROI/fig_Misha_smooth_shape.pdf}\\
%\multicolumn{2}{c}{Ring shape} \\
%\includegraphics[scale=0.45,clip = true]{figs/fig_ROI/fig_simulation_comp_Ring.pdf}&
%\includegraphics[scale=0.45,clip = true]{figs/fig_ROI/fig_simulation_shape_Ring.pdf}\\
%\\
\end{tabular}
}
\end{frame}


\begin{frame}
\frametitle{Real data analysis: video}
\begin{center}
%\movie[width=1\textwidth, height=0.75\textwidth, borderwidth=0pt, autostart, repeat]{}{video_misha_result.avi}
\animategraphics[width=1\textwidth, height=0.75\textwidth, loop]{12}{figs/res_vid/fig_frame-}{1}{100}
\end{center}
\end{frame}



\begin{frame}
\frametitle{Conclusion and discussion}
\begin{itemize}
\item Summary
\begin{itemize}
\item Formulated calcium imaging ROI detection as a structured matrix factorization problem.
\item Developed a fast greedy algorithm.
%\item Fast ROI detection algorithm
\end{itemize}
%\item Future work
%\begin{itemize}
%\item More spatial and temporal structure.
%\item Overlapping neuron.
%\item Online ROI detection.
%\item Motion correction and background elimination.
%\end{itemize}
\item {\scriptsize Pnevmatikakis EA, Soudry D, \alert{Gao Y}, Machado TA, Merel J, Pfau D, Reardon T, Mu Y, Lacefield C, Yang W, Ahrens M, Bruno R, Jessell TM, Peterka DS, Yuste R, Paninski L (2016) Simultaneous denoising, deconvolution, and demixing of calcium imaging data. Neuron, 89(2), 285-299.}
\end{itemize}
\end{frame}


%%%%%%%%%%%%%%%%%%%%% MEFN %%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\section[]{Maximum Entropy Flow Networks}

\begin{frame}
\frametitle{Maximum entropy principle}
\begin{itemize}
\item \alert{Entropy}: for a continuous distribution with density $p(\vz)$ where $\vz \in \mathbb{R}^d$, the entropy is defined as
\[H(p) = - \int p(\vz) \log p(\vz)  d \vz = E_{\mZ \sim p} \left[ -\log p(\mZ) \right].\]
\item A popular measure of diversity or information content.
\item \alert{Maximum entropy (ME) principle}: Subject to some given prior knowledge (usually in the form of expectation or support constraints), the distribution that makes \alert{minimal additional assumptions} is that which has the \alert{largest entropy} of any distribution obeying those constraints
\end{itemize}
\end{frame}

\begin{frame}
\frametitle{Applications of maximum entropy}
\begin{itemize}
\item \alert{Neuroscience}: generate a distribution of neural activity by specifying a set of features (pairwise correlation etc.) for hypothesis testing.
\item \alert{Texture modeling}: generate an image with a certain texture by specifying expected value of features relevant to texture.
%\item \alert{Ecology}: Fit species distribution with certain feature constraints (altitude, temperature etc.).
%\item \alert{Finance}: Fit the risk-neutral distribution of an asset given a list of option prices.
\item Ecology, natural language processing, finance...
\end{itemize}
\end{frame}

\begin{frame}
\frametitle{Maximum entropy problem}
\begin{itemize}
\item Maximum entropy (ME) problem
\begin{eqnarray*}
%\substack{E_{\mZ\sim p}[T(\mZ)]=0 \\ \supp(p)=\mathcal{X}}
p^{*} &=& \underset{p}{\text{maximize}}~~~ H(p) \\
&& \text{subject to} ~~~E_{\mZ\sim p}[T(\mZ)]=0 \nonumber\\
&& ~~~~~~~~~~~~~~~~~~\supp(p)=\mathcal{Z} \nonumber,
\end{eqnarray*}
where $T(\vz) = (T_1(\vz),...,T_m(\vz)):\mathcal{Z}\rightarrow \mathbb{R}^{m}$ is the vector of known statistics, and $\mathcal{Z}$ is the given support.
%\item Useful for model formulation and hypothesis testing.
%\item \alert{Problem}: hard to learn and sample from maximum entropy distribution when constraints are complicated
\end{itemize}
\end{frame}

\begin{frame}
\frametitle{Gibbs distribution}
\begin{itemize}
\item Lagrange multipliers argument $\Rightarrow$ an exponential family form for ME distribution (Gibbs distribution):
\[\label{curmet}
p^*(\vz) \propto e^{<\eta, T(\vz)>} \mathbbm{1}(\vz \in \mathcal{Z})
\]
For some $\eta \in  \mathbb{R}^m$ such that $E_{\mZ\sim p^*}[T(\mZ)]=0$.
\item Hard to Identify $\eta$.
\item Hard to sample from the distribution.
\item \alert{Question}: is there a better way to do this?
\end{itemize}
\end{frame}


\begin{frame}
\frametitle{Maximum entropy flow network (MEFN)}
\begin{itemize}
\item Identify a transformation $f_{\phi^*} \in \mathcal{F}=\{f_\phi:  \mathbb{R}^d \rightarrow \mathbb{R}^d, \phi \in \mathbb{R}^q\}$ that transforms a simple distribution $\mZ_0 \sim p_0$ (Gaussian) to $\mZ = f_{\phi^*}(\mZ_0) \sim p_{\phi^*}$ that approximates the ME distribution.
\begin{eqnarray*}
%\substack{E_{\mZ\sim p}[T(\mZ)]=0 \\ \supp(p)=\mathcal{X}}
p^{*} &=& \underset{p}{\text{maximize}}~~~ H(p) \\
&& \text{subject to} ~~~E_{\mZ\sim p}[T(\mZ)]=0 \nonumber\\
&& ~~~~~~~~~~~~~~~~~~\supp(p)=\mathcal{Z} \nonumber,\\
&&~~~~~~~~\Downarrow\\
%\substack{E_{\mZ\sim p}[T(\mZ)]=0 \\ \supp(p)=\mathcal{X}}
\phi^{*} &=& \underset{\phi}{\text{maximize}}~~~ H(p_\phi) \\
&& \text{subject to} ~~~R(\phi) = E_{\mZ_0\sim p_0}[T(f_\phi(\mZ_{0}))]=0 \nonumber\\
&& ~~~~~~~~~~~~~~~~~~\supp(p_\phi)=\mathcal{Z}. \nonumber
\end{eqnarray*}
%where $p_\phi$ is the distribution of $\mZ = f_\phi(\mZ_0)$ for $\mZ_0 \sim p_0$.
%\[p_\phi(\vz) = p_0( f_\phi^{-1}(\vz) ) | \det (J_\phi (\vz)) |^{-1}\]
\end{itemize}
\end{frame}


\begin{frame}
\frametitle{Normalizing flow}
\begin{itemize}
\item \alert{Normalizing flow}: A \alert{smooth and invertible} transformation. %\[\mathcal{F}=\{f_\phi:  \mathbb{R}^d \rightarrow \mathbb{R}^d, \phi \in \mathbb{R}^q\}.\]
\item Explicit probability density form for $Z = f_\phi(Z_0) \sim p_\phi$ by \alert{change-of-variable theorem}
\[p_\phi(\vz) = p_0\left( f_\phi^{-1}(\vz) \right) | \det \left(J_\phi (\vz)\right) |^{-1}.\]
\item Constructing a flexible family of normalizing flows by composing simple normalizing flows (\alert{deep learning})
\[\mZ = f_{k}\circ f_{k-1}\circ\dots \circ f_{1}(\mZ_{0}).\]
\end{itemize}
\end{frame}


\begin{frame}
\frametitle{Augmented Lagrangian method}
\begin{itemize}
\item Augmented Lagrangian method minimizes the objective
\[L(\phi; \lambda, c) = - H(p_\phi) + \lambda^{\top} R(\phi) + \dfrac{c}{2}||R(\phi)||^2.\]
for a sequence of $\lambda \in \mathbb{R}^m$ and $c \geq 0$.
\item Update rule: at iteration $k$, given $\lambda_k$ and $c_k$.
\begin{itemize}
\item Find $\phi_k$ that optimizes $L(\phi; \lambda_k, c_k)$.
\item Update $\lambda$ and $c$ by
\[\begin{split}
\lambda_{k+1} =& \lambda_k + c_k R(\phi_k)\\
c_{k+1} =&\begin{cases}
\beta c_{k} & ||R(\phi_{k})||>\gamma ||R(\phi_{k-1})|| \\
c_k & \text{otherwise}
\end{cases}
\end{split}\]
for some $\gamma \in (0,1)$, $\beta > 1$
\end{itemize}
\end{itemize}
\end{frame}

\begin{frame}
\frametitle{Augmented Lagrangian method in stochastic setting}
\begin{itemize}
\item Both $H(p_\phi)$ and $R(\phi)$ are intractable, but we can approximate with an unbiased estimation,
\[\begin{split}
H(p_\phi) \approx& \frac{1}{n} \sum_{i=1}^n -\log p_\phi(\vz^{(i)}),\\
R(\phi) \approx& \frac{1}{n} \sum_{i=1}^nT(f_\phi( \vz^{(i)})),
\end{split}\]
where $\vz^{(i)} \sim p_0$. We can then optimize the objective by stochastic gradient descent.
\end{itemize}
\end{frame}

\begin{frame}
\frametitle{Texture modeling (overview)}
\begin{itemize}
\item Goal: construct a distribution on images that mimics a given texture (from a training image).
\begin{itemize}
\item Authenticity: the samples from the distribution should mimic the given texture (constraints).
\item Sample diversity: the distribution should generate a variety of images with certain texture (entropy).
\end{itemize}
\item Natural application of maximum entropy modeling.
\end{itemize}
\end{frame}


%http://www.texturex.com/Stone-Textures/
\begin{frame}
\frametitle{Texture modeling (formulation)}
\begin{itemize}
\item Sample from the space of images $\vz \in [0,1]^{d = w \times h \times c}$ where $w \times h$ is the image size, $c$ is the number of channels ($c = 3$ for RGB representation). 
\item Texture loss: $T(\vz): [0,1]^{w \times h \times c} \rightarrow \mathbb{R}$, a complicated form proposed in \textcite{ulyanov2016texture}.
\item \textcite{ulyanov2016texture} proposes texture net, which solves the problem \alert{without considering entropy term.}
\[\min E_{\mZ \sim p_0}\left[ T(f_\phi(\mZ)) \right] \]
%\item We compare our formulation with texture net. 
\item We build MEFN using $T(\vz)$ as the expectation constraint. 
%\item We use real-nvp \parencite{dinh2016density} as the normalizing flow structure.
\end{itemize}
\end{frame}

\begin{frame}
\frametitle{Experiment: Texture modeling (result)}
\centering
{\scriptsize
\begin{tabular}[t]{ccc}
Input & \hspace{-0.3cm} \parbox{0.3\textwidth}{\centering Texture net \\\parencite{ulyanov2016texture}} & \hspace{-0.3cm} MEFN (ours) \\
%\vspace{-0cm}
\includegraphics[scale=0.3,clip = true]{figs/MEFN/fig_texture_raw.jpg}&
\hspace{-0.3cm}\includegraphics[scale=0.3,clip = true]{figs/MEFN/fig_texture_noent_0.jpg}&
\hspace{-0.3cm}\includegraphics[scale=0.3,clip = true, trim=0cm 0cm 0cm 0cm]{figs/MEFN/fig_texture_ent_0.jpg}
\vspace{0.6cm}
\\
Texture loss & \hspace{-0.3cm} Negative entropy & \hspace{-0.3cm} RGB histogram\\
\includegraphics[scale=0.33,clip = true]{figs/MEFN/fig_texture_cost.pdf}&
\hspace{-0.3cm}\includegraphics[scale=0.33,clip = true]{figs/MEFN/fig_texture_entropy.pdf}&
\hspace{-0.3cm}\includegraphics[scale=0.33,clip = true]{figs/MEFN/fig_texture_RGB_hist.pdf} \\
%\\
\end{tabular}
}
\end{frame}


\begin{frame}
\frametitle{Experiment: Texture modeling (diversity measure)}
%\begin{itemize}
%\item We provide two numerical measures for assessing sample diversity given a set of images $\{\vz^{(1)},...,\vz^{(n)}\}$.
%\end{itemize}
\centering
\begin{table}[htpb]
%\caption{Quantitative measure of image diversity using $20$ randomly sampled images}
%\label{sample-table}
\begin{center}
\begin{tabular}{l | l | lll}
Method  & $d_{L^2}$& SST & SSW & SSB \\
\hline 
Texture net  & 0.077 & 0.043 & 0.036 & 0.006\\
MEFN   & \alert{0.113} & 0.058 & \alert{0.054} & \alert{0.005} \\
\end{tabular}
\end{center}
\begin{itemize}
\item For $n=20$ randomly sampled image $\vz^{(1)},...,\vz^{(20)}$
\item $d_{L^2} = \text{mean}_{i \neq j} \frac{1}{d}\|\vz^{(i)} - \vz^{(j)} \|_2^2$: Mean Euclidean distance.
\item ANOVA: $\text{SST} = \text{SSW} + \text{SSB}$, $\bar{\vz} = \frac{1}{n} \sum_{i} \vz^{(i)}$.
\begin{itemize}
\item $\text{SST} = \frac{1}{nd}\sum_{i,k} (z_k^{(i)} - \bar{z})^2$: Total var.
\item $\text{SSW} = \frac{1}{nd}\sum_{i,k} (z_k^{(i)} - \bar{z}_k)^2$: residue var, larger $\Rightarrow$ better.
\item $\text{SSB} = \frac{1}{n}\sum_{k} (\bar{z}_k - \bar{z})^2$: mean image var, smaller $\Rightarrow$ better.
\end{itemize}
\end{itemize}
\label{tab:texture}
\end{table}
\end{frame}

\begin{frame}
\frametitle{Conclusion and discussion}
\begin{itemize}
\item Summary
\begin{itemize}
%\item Bridging information theory and deep learning
\item Solved maximum entropy problem by optimizing a normalizing flow.
\item Combined augmented Lagrangian optimization with stochastic optimization.
\item Obtained promising texture modeling result.
\end{itemize}
%\item Future work
%\begin{itemize}
%\item Better normalizing flow structure
%\item Better constrained stochastic optimization algorithm
%\end{itemize}
\item Loaiza G*, \alert{Gao Y}*, Cunningham JP (2017) Maximum entropy flow networks. ICLR 2017. (*=equal contribution) 
\end{itemize}
\end{frame}

\begin{frame}
\frametitle{Acknowledgement}
\begin{center}
\Large{Thank you!}
\end{center}
\end{frame}


\begin{frame}
\frametitle{Table of Contents}
\tableofcontents
\begin{center}
\alert{End of the talk.\\ Questions?}
\end{center}
\end{frame}

\begin{frame}
\begin{center}
\Large{Backup slides}
\end{center}
\end{frame}

\begin{frame}
\frametitle{GCLDS: Supervised case, GCGLM}
\begin{itemize}
\item Given count data $x_i \in \mathbb{N}$, and associated covariates $z_i \in \mathbb{R}^p$, build Generalized Count Generalized Linear Model (GCGLM).
\[x_i \sim \mathcal{GC}(\theta(z_i), g(\cdot)), ~~~~\text{where} ~~~ \theta(z_i) =  z_i \beta.\]
\item Reminder: generalized count distribution
\[p_{\mathcal{GC}}(x; \theta, g(\cdot)) \propto \frac{\exp(\theta \cdot x + g(x) )}{x!}, ~~~x \in \mathbb{N}\]
\end{itemize}
\end{frame}

\begin{frame}
\frametitle{GCLDS: special cases of GCGLM}
\begin{centering}
{\scriptsize
\begin{tabular}{ lll}
\hline
 \multicolumn{1}{c}{\bf Model Name} &  \multicolumn{1}{c}{\bf Typical Parameterization} &  \multicolumn{1}{c}{\bf GCGLM Parametrization} \\
\hline
\parbox{0.2\textwidth}{Logistic regression } & $\begin{aligned} P(x=k) = \frac{\exp\left(k (\alpha + z\beta)\right)}{1 + \exp(\alpha + x\beta)}  \end{aligned}$ & $\begin{aligned} g(k) = \alpha k; k = 0,1\end{aligned}$\\
 \hline
   \parbox{0.2\textwidth}{Poisson regression} & $\begin{aligned} P(x=k) =& \frac{\lambda^{k}}{k!}\exp(-\lambda);\\ \lambda =& \exp(\alpha + z\beta)\end{aligned}$& $\begin{aligned} g(k) = \alpha k \end{aligned}$\\
  \hline

  %Adjacent Category Regression & $\frac{P(y = k)}{P(y = k-1)} = \logit^{-1}(\alpha_k + x\beta); k = 1,...,K$ & $g(k) = \sum_{i = 1}^{k} (\alpha_{i} + \log i); k = 0,1,...,K$\\
  \parbox{0.2\textwidth}{Adjacent category regression }  
   & $\begin{aligned}\frac{P(x = k+1)}{P(x = k)} =& \exp(\alpha_k + z\beta)\end{aligned}$ & $\begin{aligned} g(k) =& \sum_{i = 1}^{k} (\alpha_{i - 1} + \log{i});\\ k =& 0,1,...,K \end{aligned}$\\
    \hline
%     Binomial Regression & $\begin{aligned}P(y=k) =& \binom{N}{K} p^k(1-p)^{N-k}\\ k =& 0,1,...,N\\ p =& \logit^{-1}(\alpha + x \beta)\end{aligned}$ & $\begin{aligned}g(k) =&
%\alpha k + \sum\limits_{j = n - k +1}^n \log j;\\ k =& 0, 1,...,N \end{aligned}$\\
%\hline
  \parbox{0.2\textwidth}{Negative binomial regression}&  $\begin{aligned}P(x = k) =& \frac{(k+r-1)!}{k!(r-1)!}(1-p)^r p^k \\ p =& \exp (\alpha + z\beta) \end{aligned}$ & $\begin{aligned} g(k) =& \alpha k + \log{(k+r-1)!} \end{aligned}$\\
\hline
\parbox{0.2\textwidth}{COM-Poisson regression} & $\begin{aligned} P(x = k) =& \frac{\lambda^k}{(k!)^{\nu}} / \sum_{j=1}^{+\infty} \frac{\lambda^j}{(j!)^{\nu}} \\ \lambda =& \exp(\alpha + z \beta) \end{aligned}$ &
$\begin{aligned} g(k) = \alpha k + (1-\nu) \log{k!} \end{aligned}$ \\
%$\begin{aligned} &\text{Negative-binomial Regression}\\ &\text{(log link)} \end{aligned}$ & $E(Y) =   \exp (\alpha + x\beta)$ & Not Available\\
  \hline
\end{tabular}
}
\end{centering}
\end{frame}


\begin{frame}
\frametitle{GCLDS: VBEM detail}
\end{frame}



\begin{frame}
\frametitle{MEFN: Normalizing flow structures}
\begin{itemize}
\item \textcite{Rezende2015} Proposes two specific families of transformations for variational inference
\item Planar flow%: cut the space by half and distort the relations
\[f_i(\vz) = \vz + \vu_i h(\vw_i^{T} \vz+ b_i),\]
where $b_i \in \mathbb{R}$, $\vu_i, \vw_i \in \mathbb{R}^d$ and $h$ is an activation function. %in the planar case
\item Circular flow
\[f_i(\vz) = \vz + \beta_i h(\alpha_i, r_i)(\vz - \vz'_i),\]
$\beta_i \in \mathbb{R}$, $\alpha_i >0$, $\vz'_i \in \mathbb{R}^d$ , $h(\alpha,r)=1/(\alpha+r)$ and $r_i = ||\vz-\vz'_i||$.
\end{itemize}
\end{frame}

\begin{frame}
\frametitle{MEFN: real-nvp normalizing flow structure}
\begin{itemize}
\item Used in \textcite{dinh2016density} for image denstiy estimation.
\item Affine coupling layer: split variable $\vz \in \mathbb{R}^D$ into $\vz_1 \in \mathbb{R}^d$ and $\vz_2 \in \mathbb{R}^{D-d}$. linearly transform $\vz_2$ given $\vz_1$. 
\[f\left(\begin{array}{c} \vz_1\\ \vz_2 \end{array} \right) = 
\left(\begin{array}{c} \vz_1 \\ \vz_2 \odot \exp \left( s(\vz_1) \right) + t(\vz_1) \end{array}\right) \]
where $\odot$ is element-wise product. $s(\vz_1), t(\vz_1): \mathbb{R}^d \rightarrow \mathbb{R}^{D-d}$.
\item The transformation family is flexible because
\begin{itemize}
\item $s$ and $t$ can be complex while maintaining tractability (inversion and Jacobian computation). %without hurting the tractability of the Jacobian computation
\item Partitioning $\vz$ into $(\vz_1, \vz_2)$ can be chosen arbitrarily.
\end{itemize}
\end{itemize}
\end{frame}

\begin{frame}
\frametitle{Experiment: Dirichlet}
\begin{itemize}
\item Dirichlet distribution is the ME distribution on a simplex $\mathcal{S}=\{\vz = (z_1,\dots,z_{d-1}):z_i\geq 0\text{ and }\sum_{k=1}^{d-1}z_k\leq 1\}$ with expetation on the log of each coordinate $E[\log Z_k] = \kappa_k (k = 1,...,d)$, where $Z_d = 1 - \sum_{k=1}^{d-1} Z_k$.
\item $10$ layers of planar flow \parencite{Rezende2015}.
\end{itemize}
\begin{tabular}[t]{ccc}
\vspace{-0cm}
\includegraphics[scale=0.3,clip = true, trim=1cm 1.5cm 1cm 1.5cm]{figs/MEFN/normal_contours.pdf}&
\includegraphics[scale=0.3,clip = true, trim=1cm 1.5cm 1cm 1.5cm]{figs/MEFN/network_contours_moments.pdf}&
\includegraphics[scale=0.3,clip = true, trim=1cm 1.5cm 1cm 1.5cm]{figs/MEFN/true_contours.pdf}\\
{\scriptsize Initial distribution $p_0$}& \parbox{0.3\textwidth}{\scriptsize MEFN result $p_{\phi^*}$ (Control case: moment matching)}&{\scriptsize Ground truth $p^*$}\\
%& {\small (Control case: moment matching network)} & \\
%\\
\end{tabular}
\end{frame}


\end{document}

%%%%%%%%%%%% deprecated slides %%%%%%%%%%%%%%
\begin{frame}
\frametitle{fLDS: AEVB form}
\begin{itemize}
\item Generative model, for a single trial:
\[p_\theta(\vz | \vx) \propto \underbrace{p_\theta(\vz_{1})
\prod_{t=1}^{T-1}p_\theta(\vz_{r(t+1)} | \vz_{rt})}_{Gaussian} 
\underbrace{\prod_{t=1}^T p_\theta(\vx_{rt} | \vz_{rt})}_{\alert{Complicated}}
\]
%No analytical form, hard to sample.
\item Recognition model, for a single trial:
\[q_\phi(\vz | \vx) \propto \underbrace{q_\phi(\vz_{1}) 
\prod_{t=1}^{T-1}q_\phi(\vz_{r(t+1)}|\vz_{rt}) }_{Gaussian}
\underbrace{\prod_{t=1}^T q_\phi(\vz_{rt} | \vx_{rt})}_{\alert{Gaussian}}
\]
Approximates a complicated factor with a Gaussian factor dependent on the data in a complicated way.
\item Jointly Gaussian distribution with block tri-diagonal precision matrix.
\end{itemize}
\end{frame}

\begin{frame}
\frametitle{Inference algorithm: AEVB (important details)}
\begin{itemize}
\item Specific parameterization of the recognition model
\[q(\vz_r) = q_\phi(\vz_r; \vx_r) = \N\left(\mu_\phi(\vx_r), \Sigma_\phi(\vx_r)\right)\]
\vspace{-0.5cm}
\begin{itemize}
\item Block tri-diagonal precision matrix that agrees with Markovian structure
\item Potentially useful to perform filtering in an online fashion
\end{itemize}
\item Reparameterization trick for stochastic optimization
\begin{itemize}
\item Easy implementation
\item Low variance
\end{itemize}
\end{itemize}
\end{frame}

%%%%%%% OLD calcium ROI
\begin{frame}
\frametitle{Model formulation}
\begin{itemize}
\item  $X \in \mathbb{R}^{N \times T}$ represents the calcium imaging data, where each column is a (vectorized) frame that contains $N$ pixels
\item Decompose $X$ into a product of $K$ \alert{spatial components} and \alert{temporal components}% (neural acitivities)
\[X = D A^T + \text{noise}\]
\begin{itemize}
\item $D = [D_1,...,D_K] \in \mathbb{R}^{N \times K}$ represents the neuron shapes
\item $A = [A_1,...,A_K] \in \mathbb{R}^{T \times K}$ is the neural activities%, $K$ is the number of neurons.
\end{itemize}
\item Further exploit structure of the components (localized neuron shapes)
\end{itemize}
\end{frame}

\begin{frame}
\frametitle{Model formulation: objective}
\begin{itemize}
\item Structured matrix factorization
\[
%\label{equ:opt}
\begin{aligned}
& \underset{D, A}{\text{minimize}}
& & \| X - D A^T \|_2^2 + f_D(D), \\
& \text{subject to}
& & D_k \in \mathcal{D}_{w}^+; k = 1, \ldots, K,\\
& 
& & \|A_k\|_2 \leq c_k,
\end{aligned}
\]
\item $\mathcal{D}_{w}^+$: non-negative vectors whose nonzero values is within a $w \times w$ window
\item $f_D(D)$ regularizes the neuron shape (will discuss later)
\item $\|A_k\|_2 \leq c_k$ avoids degenerate solution
\end{itemize}
\end{frame}


\begin{frame}
\frametitle{Greedy algorithm}
\begin{itemize}
\item Identify ROIs one at a time, using the residual un-explained by existing signals.
\item At iteration $i$, given the current residue (unexplained by existing ROIs).
\begin{itemize}
\item \alert{Greedy identification}: Identify the location $p_k$ where the Gaussian kernel explains the largest variance of the current residue across time.
\item \alert{Shape fine tuning}: Locally optimize the spatial and temporal component.
\item \alert{Residue update}: Subtract the newly identified signal.
\end{itemize}
\end{itemize}
\end{frame}


\begin{frame}
\frametitle{Shape fine tuning}
\begin{itemize}
\item Given current residue $R$, an identified center pixel $p_k$, denote $S_k$ as a $w \times w$ window centered at $p_k$
\[
\begin{aligned}
& \underset{D_k, A_k}{\text{minimize}}
& & \| R - D_k A_k^T \|^2 + f(D_k), \\%\sum_{i = 1}^3 \lambda_i f_i(D_k),\\
& \text{subject to}
& & D_{kp} \geq 0, p \in S_k,\\
&
& & D_{kp} = 0, p \notin S_k,\\
& 
& & \|A_k\|_2 \leq c_k,
\end{aligned}
\]
%\begin{itemize}
\item $f(D_k) = \sum_{i=1}^3 \lambda_i f_i(D_k)$
%where
\begin{itemize}
\item $f_1(D_k) = \sum_p \tau_{(p, p_k)} | D_{kp} |$ encourages sparsity
\item $f_2(D_k) = \sum_p (D_{kp} - G_{p_k})^2$ encourages Gaussian shape
\item $f_3(D_k) = \sum_{\text{$p_1$ and $p_2$ are neighbors}} (D_{kp_1} - D_{kp_2})^2$ encourages smoothness
\end{itemize}
\item Optimize by block coordinate descent
%\item Need to tune regularization parameters, not necessarily good since those would shrink the estimated shape
\end{itemize}
%\end{itemize}
\end{frame}


%%%% Calcium imaging ROI, consistent notation

\begin{frame}
\frametitle{Model formulation}
\begin{itemize}
\item  $X \in \mathbb{R}^{d \times T}$ represents the calcium imaging data, where each column is a (vectorized) frame that contains $d$ pixels
\item Decompose $X$ into a product of $n$ \alert{spatial components} and \alert{temporal components}% (neural acitivities)
\[X = D A^T + \text{noise}\]
%\begin{itemize}
\item $D = [D_1,...,D_n] \in \mathbb{R}^{d \times n}$ represents the neuron shapes
\item $A = [A_1,...,A_n] \in \mathbb{R}^{T \times n}$ is the neural activities%, $K$ is the number of neurons.
%\end{itemize}
\item Further exploit structure of the components (localized neuron shapes)
\end{itemize}
\end{frame}

\begin{frame}
\frametitle{Model formulation: objective function}
\begin{itemize}
\item Structured matrix factorization
\[
%\label{equ:opt}
\begin{aligned}
& \underset{D, A}{\text{minimize}}
& & \| X - D A^T \|_2^2 + f_D(D), \\
& \text{subject to}
& & D_i \in \mathcal{D}_{w}^+; i = 1, \ldots, n,\\
& 
& & \|A_i\|_2 \leq c_i,
\end{aligned}
\]
\item $\mathcal{D}_{w}^+$: non-negative vectors whose nonzero values is within a $w \times w$ window
\item $f_D(D)$ regularizes the neuron shape (will discuss later)
\item $\|A_i\|_2 \leq c_i$ avoids degenerate solution
\end{itemize}
\end{frame}


\begin{frame}
\frametitle{Greedy algorithm}
\begin{itemize}
\item Identify ROIs one at a time, using the residual un-explained by existing signals.
\item At iteration $i$, given the current residue (unexplained by existing ROIs).
\begin{itemize}
\item \alert{Greedy identification}: Identify the location $p_i$ where the Gaussian kernel explains the largest variance of the current residue across time.
\item \alert{Shape fine tuning}: Locally optimize the spatial and temporal component.
\item \alert{Residue update}: Subtract the newly identified signal.
\end{itemize}
\end{itemize}
\end{frame}


\begin{frame}
\frametitle{Shape fine tuning}
\begin{itemize}
\item Given current residue $R$, an identified center pixel $p_i$, denote $S_i$ as a $w \times w$ window centered at $p_i$
\[
\begin{aligned}
& \underset{D_i, A_i}{\text{minimize}}
& & \| R - D_i A_i^T \|^2 + f(D_i), \\%\sum_{i = 1}^3 \lambda_i f_i(D_k),\\
& \text{subject to}
& & D_{ip} \geq 0, p \in S_i,\\
&
& & D_{ip} = 0, p \notin S_i,\\
& 
& & \|A_i\|_2 \leq c_i,
\end{aligned}
\]
%\begin{itemize}
\item $f(D_i) = \lambda_1 f_1(D_i) +  \lambda_2 f_2(D_i) +  \lambda_3 f_3(D_i)$
%where
\begin{itemize}
\item $f_1(D_i) = \sum_p \tau_{(p, p_i)} | D_{ip} |$ encourages sparsity
\item $f_2(D_i) = \sum_p (D_{ip} - G_{p_i})^2$ encourages Gaussian shape
\item $f_3(D_i) = \sum_{\text{$p_1$ and $p_2$ are neighbors}} (D_{ip_1} - D_{ip_2})^2$ encourages smoothness
\end{itemize}
\item Block coordinate descent
%\item Need to tune regularization parameters, not necessarily good since those would shrink the estimated shape
\end{itemize}
%\end{itemize}
\end{frame}


\begin{frame}
%% TODO: shorten things
\frametitle{Overview}
\centering
\includegraphics[width = 0.50\textwidth]{./figs/gclds/fig_monkey.png}
\begin{itemize}
\item Neuroscience + Big data = \alert{Opportunities}!
\end{itemize}
%https://web.stanford.edu/group/murmann_group/cgi-bin/mediawiki/index.php/Ross_Walker

\end{frame}




\begin{frame}
\frametitle{Real data analysis: video}
\begin{center}
\movie[width=1\textwidth, height=0.75\textwidth, borderwidth=0pt, autostart, repeat]{}{video_misha_result.avi}
\end{center}
\end{frame}