% Copyright (c) 2013, Jaakko Luttinen <jaakko.luttinen@aalto.fi>
% License CC-BY-SA any version or GNU GPL any version

\documentclass{article}
\usepackage{amsmath}
\usepackage{amsfonts}
\usepackage{graphicx}
%\usepackage{subfig}
% %\usepackage{upgreek}
\usepackage{bm}
\usepackage[authoryear]{natbib}

\newcommand{\mx}[1]{\mathbf{#1}}
\newcommand{\mxg}[1]{\boldsymbol{#1}}
\newcommand{\T}{^{\operatorname{T}}}
\newcommand{\invT}{^{\operatorname{-T}}}
\newcommand{\inv}{^{-1}}
\newcommand{\half}{\frac{1}{2}}
\newcommand{\der}[2]{\frac{\partial #1}{\partial #2}}
%\newcommand{\vec}{\mathoperator{vec}}
\DeclareMathOperator{\diag}{diag}
\DeclareMathOperator{\chol}{chol}
\DeclareMathOperator{\const}{const}
\DeclareMathOperator{\tr}{tr}

\title{Extra Material for \emph{Fast Variational Bayesian Linear State-Space Model}}
\author{Jaakko Luttinen}


\begin{document}

\maketitle

\begin{abstract}
  This document contains gradients for the rotation speed-ups presented in the
  paper \emph{Fast Variational Bayesian Linear State-Space Model} (Luttinen,
  ECML 2013). This document is available under the CC-BY-SA license (any
  version) or the GNU GPL license (any version).
\end{abstract}

\section{Useful matrix formulas}
\begin{align*}
  d\mx{R}\inv &= -\mx{R}\inv (d\mx{R}) \mx{R}\inv
  \\
  d\log|\mx{R}| &= \tr(\mx{R}\inv d\mx{R})
  \\
  d\tr(\mx{AR}) &= \tr(\mx{A}d\mx{R})
  \\
  d\tr(\mx{ARBR}\T) &= \tr\left((\mx{ARB} + \mx{A}\T\mx{R}\mx{B}\T)d\mx{R}\right)
\end{align*}

\section{Rotations}

\subsection{Simple Gaussian}

Prior:
\begin{align*}
  p(\mx{X}) &= \prod^N_{n=1} \mathcal{N} \left( \mx{x}_n | \mx{0},
    \mx{\Lambda}\inv \right)
\end{align*}

Transformed posterior:
\begin{align*}
  q(\mx{X}) &= \prod^N_{n=1}
  \mathcal{N} \left( \mx{x}_n | \mx{R}\mxg{\mu}_n, \mx{R}\mx{\Sigma}_n\mx{R}\T
  \right)
  \\
  \langle \mx{X} \rangle_R &= \mx{R}\langle\mx{X}\rangle
  \\
  \langle \mx{XX\T} \rangle_R &= \mx{R} \langle\mx{XX}\T\rangle \mx{R}\T
\end{align*}

Lower bound terms:
\begin{align*}
  \langle \log q(\mx{X}) \rangle_R &= - N \log|\mx{R}| + \const
  \\
  \langle \log p(\mx{X}) \rangle_R &= -\half \tr\left(
    \langle{\mx{XX}\T}\rangle_R \mx{\Lambda} \right) + \const
\end{align*}

Gradient terms:
\begin{align*}
  \der{}{\mx{R}} \log|\mx{R}| &= \mx{R}\invT
  \\
  \der{}{\mx{R}} \tr\left( \langle{\mx{XX}\T}\rangle_R \mx{\Lambda} \right) &= 2
  \mx{\Lambda} \mx{R} \langle\mx{XX}\T\rangle
\end{align*}


\subsection{Gaussian with ARD prior}

Prior:
\begin{align*}
  p(\mx{X}|\mxg{\alpha}) &= \prod^N_{n=1} \mathcal{N} \left( \mx{x}_n | \mx{0},
    \diag(\mxg{\alpha})\inv \right)
  \\
  p(\mxg{\alpha}) &= \prod^D_{d=1} \mathcal{G} ( \alpha_d | a_0, b_0)
\end{align*}

Transformed posteriors:
\begin{align*}
  q(\mx{X}) &= \prod^N_{n=1}
  \mathcal{N} \left( \mx{x}_n | \mx{R}\mxg{\mu}_n, \mx{R}\mx{\Sigma}_n\mx{R}\T
  \right)
  \\
  \langle \mx{XX\T} \rangle_R &= \mx{R} \sum^N_{n=1} (\mxg{\mu}_n\mxg{\mu}_n\T +
  \mx{\Sigma}_n) \mx{R}\T
\end{align*}
\begin{align*}
  q(\mxg{\alpha}) &= \prod^D_{d=1} \mathcal{G}( \alpha_d | a, b_d )
  \\
  a &= a_0 + \frac{N}{2}
  \\
  b_d &= b_0 + \half\left[ \langle\mx{XX}\T\rangle_R \right]_{dd}
  \\
  \langle \alpha_d \rangle_R &= \frac {a} { b_d}
  \\
  \langle \log\alpha_d \rangle_R &= \psi(a) - \log(b_d)
\end{align*}

Lower bound terms:
\begin{align*}
  \langle \log q(\mx{X}) \rangle_R &= - N \log|\mx{R}| + \const
  \\
  \langle \log p(\mx{X}|\mxg{\alpha}) \rangle_R &= -\half \tr\left(
    \langle{\mx{XX}\T}\rangle_R \langle\diag(\mxg{\alpha})\rangle_R \right) +
  \frac{N}{2} \sum^D_{d=1} \langle\log\alpha_d\rangle_R + \const
  \\
  \langle \log q(\mxg{\alpha}) \rangle_R &= \sum^D_{d=1} \langle \log \alpha_d
  \rangle_R + \const
  % \langle \log q(\mxg{\alpha}) \rangle_R &= - \sum^D_{d=1} \log \left( b_0 +
  %   \half\left[ \langle\mx{XX}\T\rangle_R \right]_{dd} \right) + \const
  \\
  \langle \log p(\mxg{\alpha}) \rangle_R &= \sum^D_{d=1} \langle
  \mathcal{G}(\alpha_d | a_0, b_0) \rangle_R
  \\
  &= (a_0-1) \sum^D_{d=1}\langle \log\alpha_d \rangle_R - b_0
  \sum^D_{d=1}\langle\alpha_d\rangle_d + \const
\end{align*}

Gradient terms:
\begin{align*}
  \der{}{\mx{R}} \log|\mx{R}| &= \mx{R}\invT
  \\
  \der{}{\mx{R}} \sum_{d=1}^D \langle \alpha_d \rangle_R &= -
  \begin{bmatrix}
    \frac{\langle\alpha_1\rangle_R}{b_1}
    \\
    &\ddots&
    \\
    && \frac{\langle\alpha_D\rangle_R}{b_D}
  \end{bmatrix}
  \mx{R} \langle\mx{XX}\T\rangle
  \\
  \der{}{\mx{R}} \sum_{d=1}^D \langle \log\alpha_d \rangle_R &= -
  \begin{bmatrix}
    \frac{1}{b_1} & &
    \\
    & \ddots &
    \\ & & \frac{1}{b_D}
  \end{bmatrix}
  \mx{R} \langle \mx{XX}\T \rangle
  \\
  \der{}{\mx{R}} \tr\left( \langle{\mx{XX}\T}\rangle_R
    \langle\diag(\mxg{\alpha})\rangle_R \right) &= 2
  \begin{bmatrix}
    \langle\alpha_1\rangle_R &&
    \\
    &\ddots&
    \\
    && \langle\alpha_D\rangle_R
  \end{bmatrix}
  \mx{R} \langle\mx{XX}\T\rangle
  \\
  &\quad -
  \begin{bmatrix}
    \frac{\langle\alpha_1\rangle_R}{b_1}
    [\mx{R}]_{1:}\langle\mx{XX}\T\rangle[\mx{R}]_{1:}\T && \\ &\ddots& \\ &&
    \frac{\langle\alpha_D\rangle_R}{b_D}
    [\mx{R}]_{D:}\langle\mx{XX}\T\rangle[\mx{R}]_{D:}\T
  \end{bmatrix}
  \mx{R} \langle\mx{XX}\T\rangle
\end{align*}



\subsection{Gaussian Markov chain}

Prior:
\begin{align*}
  p(\mx{x}_0) &= \mathcal{N}(\mx{x}_0|\mxg{\mu}, \mx{\Lambda}\inv)
  \\
  p(\mx{x}_n|\mx{x}_{n-1},\mx{A}) &= \mathcal{N}(\mx{x}_n|\mx{Ax}_{n-1},
  \mx{I})
  %\diag(\mxg{\alpha})\inv)
  \\
  &\text{some gaussian model for }\mx{A}
\end{align*}

Transformed posteriors:
\begin{align*}
  q_R(\mx{X}) &= \mathcal{N}\left([\mx{X}]_: \left|
      (\mx{I}\otimes\mx{R})\mxg{\mu},
      (\mx{I}\otimes\mx{R})\mx{\Sigma}(\mx{I}\otimes\mx{R})\T \right. \right)
  \\
  % q_R(\mx{A}) &= \prod^D_{d=1} \mathcal{N}\left( \mx{A}_{d:} |
  %   \mx{R}_{d:}\mx{A}\mx{R}\inv, \left(\sum^D_{k=1}\mx{R}_{kd}\right)^2
  %   \mx{R}\invT \mx{\Sigma}_A^{(d)} \mx{R}\inv \right)
  % \\
  q_R(\mx{A}) &= \prod^D_{d=1} \mathcal{N}\left( \mx{A}_{d:} |
    \mx{R}_{d:}\mx{A}\mx{R}\inv, \left(\sum^D_{k=1}\mx{R}_{kd}\right)^2
    \mx{R}\invT \mx{\Sigma}_A^{(d)} \mx{R}\inv \right)
\end{align*}

Lowerbound terms:
\begin{align*}
  \langle \log q(\mx{X}) \rangle_R &= -(N+1) \log|\mx{R}| + \const
  \\
  \langle \log p(\mx{X}|\mx{A}) \rangle_R &= \tr\Bigg( -\half \mx{\Lambda}
  \langle\mx{x}_0\mx{x}_0\T\rangle_R + \mx{\Lambda} \mxg{\mu}
  \langle\mx{x}_0\rangle_R\T + \sum^N_{n=1}\bigg[
  -\half\langle\mx{x}_n\mx{x}_n\T\rangle_R -
  \\
  &\quad \half \langle\mx{A}\T\mx{A}\rangle_R
  \langle\mx{x}_{n-1}\mx{x}_{n-1}\T\rangle_R + 
  \langle\mx{A}\rangle_R \langle\mx{x}_{n-1}\mx{x}_n\T\rangle_R\bigg] \Bigg)
\end{align*}

Relevant moments:
\begin{align*}
  \mx{c}\T \langle \mx{x}_n \rangle_R &= \mx{c}\T \mx{R} \langle \mx{x}_n
  \rangle
  \\
  \langle \mx{x}_n\mx{x}_n\T \rangle_R &= \mx{R} \langle \mx{x}_n\mx{x}_n\T
  \rangle \mx{R}\T
  % \\
  % \langle \mx{x}_{n-1}\mx{x}_n\T \rangle_R &= \mx{R} \langle
  % \mx{x}_{n-1}\mx{x}_n\T \rangle \mx{R}\T
  % \\
  % \langle \mx{A} \rangle_R &= \mx{R} \langle\mx{A}\rangle \mx{R}\inv
  % \\
  % \langle \mx{A}\T \mx{A} \rangle_R &= \mx{R}\invT \left(\sum^D_{d=1}
  %   \langle\mx{A}_{d:}\T\mx{A}_{d:}\rangle (\mx{R}_{:d}\T\mx{R}_{:d}) \right)
  % \mx{R}\inv
  \\
  \tr\left( \langle\mx{A}\rangle_R \langle\mx{x}_{n-1}\mx{x}_n\T\rangle_R
  \right) &= \tr\left( \mx{R}\T\mx{R} \langle\mx{A}\rangle
    \langle\mx{x}_{n-1}\mx{x}_{n}\T\rangle \right)
  \\
  \tr\left( \langle\mx{A}\T\mx{A}\rangle_R
    \langle\mx{x}_{n-1}\mx{x}_{n-1}\T\rangle_R \right) &= \tr\left(
    \langle\mx{A}\rangle\T\mx{R}\T\mx{R} \langle\mx{A}\rangle
    \langle\mx{x}_{n-1}\mx{x}_{n-1}\T\rangle \right) + \sum^D_{d=1}
  \left(\sum^D_{k=1}\mx{R}_{kd}\right)^2 \tr\left( \mx{\Sigma}_A^{(d)}
    \langle\mx{x}_{n-1}\mx{x}_{n-1}\T\rangle \right)
  % \\
  % \langle \mx{A}\T  \mx{A} \rangle_R &= \mx{R}\invT
  % \langle\mx{A}\rangle\T \mx{R}\T  \mx{R}
  % \langle\mx{A}\rangle \mx{R}\inv + \mx{R}\invT \left(\sum^D_{d=1}
  %   \mx{\Sigma}_A^{(d)} \sum^D_{k=1} \alpha_k r_{kd}^2 \right) \mx{R}\inv
\end{align*}

Gradient terms:
\begin{align*}
  \der{}{\mx{R}} \log|\mx{R}| &= \mx{R}\invT
  \\
  \der{}{\mx{R}} \mx{c}\T \langle \mx{x}_n \rangle_R &= \mx{c}
  \langle\mx{x}_n\rangle\T
  \\
  \der{}{\mx{R}} \tr\left( \mx{C} \langle \mx{x}_n\mx{x}_n\T \rangle_R \right)
  &= 2 \mx{C} \mx{R} \langle \mx{x}_n\mx{x}_n\T \rangle
  \\
  \der{}{\mx{R}} \tr\left( \langle\mx{A}\rangle_R
    \langle\mx{x}_{n-1}\mx{x}_n\T\rangle_R\right) &= \mx{R} \left( \langle
    \mx{x}_{n-1}\mx{x}_n\T \rangle\T \langle \mx{A} \rangle\T + \langle \mx{A}
    \rangle \langle \mx{x}_{n-1}\mx{x}_n\T \rangle \right)
  \\
  \der{}{\mx{R}} \tr\left( \langle\mx{A}\T\mx{A}\rangle_R
    \langle\mx{x}_{n-1}\mx{x}_{n-1}\T\rangle_R \right) &= 2\mx{R}
  \langle\mx{A}\rangle \langle\mx{x}_{n-1}\mx{x}_{n-1}\T\rangle
  \langle\mx{A}\rangle\T +
  \\
  &\quad
  2
  \begin{bmatrix}
    1 & \cdots & 1
    \\
    \vdots &  \ddots & \vdots
    \\
    1 & \cdots & 1
  \end{bmatrix}
  \mx{R}
  %\left( \langle\mx{A}\rangle \langle\mx{x}_{n-1}\mx{x}_{n-1}\T\rangle
  %  \langle\mx{A}\rangle\T +
    \begin{bmatrix}
      \tr\left( \mx{\Sigma}_A^{(1)}
        \langle\mx{x}_{n-1}\mx{x}_{n-1}\T\rangle \right) &&
      \\
      &\ddots&
      \\
      && \tr\left( \mx{\Sigma}_A^{(D)}
        \langle\mx{x}_{n-1}\mx{x}_{n-1}\T\rangle \right)
    \end{bmatrix}
  %\right)
\end{align*}

\subsection{Left-right rotation of Gaussian with ARD prior}

Prior:
\begin{align*}
  p(\mx{X}|\mxg{\alpha}) &= \prod^N_{n=1} \mathcal{N} \left( \mx{x}_n | \mx{0},
    \diag(\mxg{\alpha})\inv \right)
  \\
  p(\mxg{\alpha}) &= \prod^D_{d=1} \mathcal{G} ( \alpha_d | a_0, b_0)
\end{align*}

Transformed posteriors, inspired by (but not exactly) $\mx{RXQ}\T$:
\begin{align*}
  q_*(\mx{X}) &= \prod^N_{n=1} \mathcal{N} \left( \mx{x}_n | \sum_{m=1}^Nq_{nm}
    \mx{R}\mxg{\mu}_n, \left(\sum^N_{m=1}q_{mn}\right)^2
    \mx{R}\mx{\Sigma}_n\mx{R}\T \right)
  \\
  q_*(\mxg{\alpha}) &= \prod^D_{d=1} \mathcal{G}( \alpha_d | a, b_d )
  \\
  a &= a_0 + \frac{N}{2}
  \\
  b_d &= b_0 + \half\left[ \langle\mx{XX}\T\rangle_* \right]_{dd}
\end{align*}

Lower bound terms:
\begin{align*}
  \langle \log q(\mx{X}) \rangle_* &= - N \log|\mx{R}| - D
  \sum^N_{n=1}\log\left| \sum^N_{m=1} q_{mn} \right| + \const
  \\
  \langle \log p(\mx{X}|\mxg{\alpha}) \rangle_* &= -\half \tr\left(
    \langle{\mx{XX}\T}\rangle_* \langle\diag(\mxg{\alpha})\rangle_* \right) +
  \frac{N}{2} \sum^D_{d=1} \langle\log\alpha_d\rangle_* + \const
  \\
  \langle \log q(\mxg{\alpha}) \rangle_* &= \sum^D_{d=1} \langle \log \alpha_d
  \rangle_* + \const
  % \langle \log q(\mxg{\alpha}) \rangle_R &= - \sum^D_{d=1} \log \left( b_0 +
  %   \half\left[ \langle\mx{XX}\T\rangle_R \right]_{dd} \right) + \const
  \\
  \langle \log p(\mxg{\alpha}) \rangle_* &= (a_0-1) \sum^D_{d=1}\langle
  \log\alpha_d \rangle_* - b_0 \sum^D_{d=1}\langle\alpha_d\rangle_* + \const
\end{align*}

Relevant moments:
\begin{align*}
  \langle \alpha_d \rangle_* &= \frac {a} { b_d}
  \\
  \langle \log\alpha_d \rangle_* &= \psi(a) - \log(b_d)
  \\
  \langle \mx{XX\T} \rangle_Q &=\langle\mx{X}\rangle\mx{Q}\T \mx{Q}
  \langle\mx{X}\rangle\T + \sum^N_{n=1} \left(\sum^N_{m=1}q_{mn}\right)^2
  \mx{\Sigma}_n
  \\
  \langle \mx{XX\T} \rangle_* &= \mx{R} \langle \mx{XX\T} \rangle_Q \mx{R}\T
  % \\
  % \langle \mx{XX\T} \rangle_* &= \mx{R}\langle\mx{X}\rangle\mx{Q}\T \mx{Q}
  % \langle\mx{X}\rangle\T \mx{R}\T + \mx{R} \sum^N_{n=1}
  % \left(\sum^N_{m=1}q_{mn}\right)^2 \mx{\Sigma}_n\mx{R}\T
  % \\
  % \tr\left( \langle{\mx{XX}\T}\rangle_* \langle\diag(\mxg{\alpha})\rangle_*
  % \right) &=
\end{align*}

Gradient terms with respect to $\mx{R}$:
\begin{align*}
  \der{}{\mx{R}} \log|\mx{R}| &= \mx{R}\invT
  \\
  \der{}{\mx{R}} \sum_{d=1}^D \langle \alpha_d \rangle_* &= -
  \begin{bmatrix}
    \frac{\langle\alpha_1\rangle_*}{b_1}
    \\
    &\ddots&
    \\
    && \frac{\langle\alpha_D\rangle_*}{b_D}
  \end{bmatrix}
  \mx{R} \langle\mx{XX}\T\rangle_Q
  \\
  \der{}{\mx{R}} \sum_{d=1}^D \langle \log\alpha_d \rangle_* &= -
  \begin{bmatrix}
    \frac{1}{b_1} & &
    \\
    & \ddots &
    \\ & & \frac{1}{b_D}
  \end{bmatrix}
  \mx{R} \langle \mx{XX}\T \rangle_Q
  \\
  \der{}{\mx{R}} \tr\left( \langle{\mx{XX}\T}\rangle_*
    \langle\diag(\mxg{\alpha})\rangle_* \right) &= 2
  \begin{bmatrix}
    \langle\alpha_1\rangle_* &&
    \\
    &\ddots&
    \\
    && \langle\alpha_D\rangle_*
  \end{bmatrix}
  \mx{R} \langle\mx{XX}\T\rangle_Q
  \\
  &\quad -
  \begin{bmatrix}
    \frac{\langle\alpha_1\rangle_*}{b_1}
    [\mx{R}]_{1:}\langle\mx{XX}\T\rangle[\mx{R}]_{1:}\T && \\ &\ddots& \\ &&
    \frac{\langle\alpha_D\rangle_*}{b_D}
    [\mx{R}]_{D:}\langle\mx{XX}\T\rangle[\mx{R}]_{D:}\T
  \end{bmatrix}
  \mx{R} \langle\mx{XX}\T\rangle_Q
\end{align*}

Let us define a helpful function:
\begin{align*}
  \mx{\Psi}(\mx{v}) &= \der{}{\mx{Q}} \frac{1}{2} \tr\left( \diag(\mx{v})
    \langle\mx{XX}\T\rangle_* \right)
  \\
  &= \mx{Q} \langle\mx{X}\rangle\T \mx{R}\T \diag(\mx{v}) \mx{R}
  \langle\mx{X}\rangle
  \\
  &\quad +
  \begin{bmatrix}
    1 & \cdots & 1
    \\
    \vdots & \ddots & \vdots
    \\
    1 & \cdots & 1
  \end{bmatrix}
  \mx{Q}
  \begin{bmatrix}
    \tr\left( \mx{R}\T \diag(\mx{v}) \mx{R} \mx{\Sigma}_1 \right) & &
    \\
    & \ddots &
    \\
    &&\tr\left( \mx{R}\T \diag(\mx{v}) \mx{R} \mx{\Sigma}_N \right)
  \end{bmatrix}
\end{align*}

Gradient terms with respect to $\mx{Q}$:
\begin{align*}
  \der{}{\mx{Q}} \sum^N_{n=1}\log\left| \sum^N_{m=1} q_{mn} \right| &=
  \begin{bmatrix}
    1 \\ \vdots \\ 1
  \end{bmatrix}
  \begin{bmatrix}
    \left(\sum^N_{m=1} q_{m1}\right)\inv & \cdots & \left(\sum^N_{m=1}
      q_{mN}\right)\inv
  \end{bmatrix}
  \\
  \der{}{\mx{Q}} \sum_{d=1}^D \langle \alpha_d \rangle_* &= - \mx{\Psi} \left(
    \frac {\langle\mxg{\alpha}\rangle_*} {\mx{b}}\right)
  \\
  \der{}{\mx{Q}} \sum_{d=1}^D \langle \log\alpha_d \rangle_* &= - \mx{\Psi}
  \left( \frac {1} {\mx{b}} \right)
  \\
  \der{}{\mx{Q}} \tr\left( \langle{\mx{XX}\T}\rangle_*
    \langle\diag(\mxg{\alpha})\rangle_* \right) &= 2
  \mx{\Psi}\left(\langle\mxg{\alpha}\rangle_*\right) -
  \mx{\Psi}\left(\diag\langle\mx{XX}\T\rangle_* \circ \frac
    {\langle\mxg{\alpha}\rangle_*} {\mx{b}} \right)
\end{align*}


\end{document}

%%% Local Variables: 
%%% mode: latex
%%% TeX-PDF-mode: t
%%% TeX-master: t
%%% End: 
