\chapter{\glstopfull{EM}}%
\label{cha:Expectation Maximization}
\gls{EM} ist ein Verfahren zur Optimierung,
welches auf alle \glspl{latent_variable_model} angewandt werden kann,
wobei jedes \gls{latent_variable_model} bekannte Variablen $\bm x$ und unbekannte Variablen $\bm z$ enthält.
Hierdurch ergibt sich
\begin{itemize}
    \item Parametric model:\tabto{4cm}$p_{\bm\theta}(\bm x,\bm z)$
    \item \glslink{marginal}{Marginal distribution}:\tabto{4cm}$\underbrace{p_{\bm\theta}(\bm x) = \sum_z p(\bm x,z)}_{\text{discrete latent variable}},
            \qquad \underbrace{p_{\bm\theta} = \int_{\bm z} p_{\bm\theta}(\bm x,\bm z) d\bm z}_{\text{continuous latent variable}} $
\end{itemize}
Das Ziel ist es die (\glslink{marginal}{marginal}) Log-Likelihood
\begin{equation} \label{eq:latent_variable_model_marginal_log-likelihood}
    L(\bm\theta) = \sum_{i=1}^N \log p_{\bm\theta}(\bm x_i) = \sum_{i=1}^N \log \left( \sum_z p_{\bm\theta}(\bm x_i,z) \right) 
\end{equation}
zu optimieren.
Zwei Verfahren für die Durchführung der \gls{EM} sind bekannt,
wobei in der Vorlesung nur auf die \dref{sec:Decomposition in lower-bound and KL-term} eingegangen wird.
Bei der \gls{EM} handelt es sich um ein Verfahren,
welches lediglich ein lokales Optimum findet,
weshalb die Initialisierung sehr wichtig ist.
Hierbei werden die Cluster-Mittelpunkte meist über den \dref{sec:K-Means Algorithm}us approximiert und ein fester Wert für die Kovarianzen gewählt.

\section{Decomposition in lower-bound and \glstopshort{KL}-term}%
\label{sec:Decomposition in lower-bound and KL-term}
Die \glslink{marginal}{Marginal} Log-Likelihood eines \gls{latent_variable_model}s lässt sich wie folgt zerlegen (Herleitung: {\color{red}Vorlesung 11 Folie 24})
\begin{align} \label{eq:decomposition_in_lower-bound_and_KL-term}
    \underbrace{\log p(\bm x|\bm\theta)}_{\text{\glslink{marginal}{marginal} log-like}} 
        &= \mathcal L(q|\bm\theta) + \nomeq{kl_divergence}(q(z)\|p(z|\bm x))\\
        &= \underbrace{\sum_z q(z)\log \frac{p(\bm x,z|\bm\theta)}{q(z)}}_{\text{Lower-Bound $\mathcal L(q,\bm\theta)$}}
            + \underbrace{\sum_z q(z)\log \frac{q(z)}{p(z|\bm x)}}_{\text{\noms{kl_divergence}: $\nomeq{kl_divergence}(q(z)\|p(z|\bm x))$}}
\end{align}
Nach dieser Zerlegung wird die \gls{EM} zu einem iterativen Verfahren mit 2 Schritten (\nameref{sub:Expectation-Step} und \nameref{sub:Maximization-Step}).
Dieses Verfahren funktioniert auch,
wenn nicht nur ein einzelner Datenpunkt $\bm x$ betrachtet wird.
Für einen Datensatz mit mehreren Datenpunkten definiert sich die Lower-Bound Funktion durch
\begin{equation} \label{eq:lower-bound_for_dataset}
    \mathcal L(q,\bm\theta) = \sum_i \left( \int_z q_i(z)\log p(\bm x_i,z|\bm\theta)dz - \int_z q_i(z)\log q_i(z)dz \right) 
\end{equation}
\begin{itemize}
    \item eine latente Variable $z_i$ für jeden Datenpunkt $\bm x_i$
    \item wenn $z$ diskret mit $K$ verschiedenen Werten ist, gilt $q_i(z=k) = p(z=k|\bm x_i,\bm\theta_{\text{old}})$.
        Dies kann als eine $N\times K$ Matrix dargestellt werden.
\end{itemize}

\subsection{Expectation-Step}%
\label{sub:Expectation-Step}
Der Expectation-Step ist der erste der 2 Schritte der \gls{EM} mit \nameref{sec:Decomposition in lower-bound and KL-term}.
Hierbei ist es das Ziel,
das $q(z)$ zu finden,
welches die \gls{KL} Divergenz (\cref{sec:KL-Divergenz}) minimiert.
\begin{equation} \label{eq:expectation_step}
    q(z) = \argmin_q \nomeq{kl_divergence}(q(z)\|p(z|\bm x))
\end{equation}
Für diskrete $z$ (z.B. bei \dref{sec:Mixture Models}) hat dieses Problem eine \gls{closed_form_solution}.
\begin{equation} \label{eq:expectation_step_closed_form_solution}
    q(z) = p(z|\bm x,\bm\theta_{\text{old}}) = \frac{p(\bm x,z|\bm\theta_{\text{old}})}{\sum_z p(\bm x,z|\bm\theta_{\text{old}})}
\end{equation}
Dies hat zur folge,
dass (für diskrete $z$) $\nomeq{kl_divergence}(q(z)\|p(z|\bm x)) = 0$ und daher wegfällt,
weshalb nach dem Expectation-Step $\log p(\bm x|\bm\theta) = \mathcal L(q,\bm\theta)$ gilt.

\subsection{Maximization-Step}%
\label{sub:Maximization-Step}
Nachdem im Expectation-Step die \gls{KL} Divergenz minimiert wurde (0 für diskrete $z$),
kann nun die Lower-Bound maximiert werden um im Gegenzug die gesamte Log-Likelihood zu erhöhen (maximieren (für diskrete $z$)).
\begin{equation} \label{eq:maximization-step}
    \bm\theta = \argmax_{\bm\theta}\mathcal L(q,\bm\theta) 
        = \argmax_{\bm\theta}\sum_z q(z)\log p(\bm x,z|\bm\theta) + \text{const}
\end{equation}
\begin{itemize}
    \item $q(z) = p(z|\bm x,\bm\theta_{\text{old}})$
\end{itemize}


\section{\glstoplong{EM} for \glstopshortpl{GMM}}%
\label{sec:Expectation Maximization for GMMs}
\begin{mybox}
    \textbf{\Large \glsxtrshort{EM} for \glsxtrshortpl{GMM}}\\ 
    \begin{enumerate}
        \item \textbf{Initialisierung:} Komponenten und Koeffizienten des \gls{GMM} (z.B. auf Basis von k-means) erstellen
        \item \textbf{Expectation"~Step:} Responsibilities errechnen
            \begin{equation*}
                q_{ik} = \frac{\pi_k\nomeq{gaussian_distribution}(\bm x_i|\bm\mu_k,\nomeq{covariance}_k)}{\sum_{j=1}^K \pi_j\nomeq{gaussian_distribution}(\bm x_i|\bm\mu_j,\nomeq{covariance}_j)}
            \end{equation*}
        \item \textbf{Maximization"~Step:} Koeffizienten, Komponenten"~Durchschnitte und Kom\-po\-nen\-ten"~Va\-ri\-an\-zen updaten
            (Herleitung: {\color{red} Vorlesung 11 Folie 16})
            \begin{equation*}
                \pi_k = \dfrac{\sum_i q_{ik}}{N}\qquad
                \bm\mu_k = \dfrac{\sum_i q_{ik}\bm x_i}{\sum_i q_{ik}}
                \qquad\nomeq{covariance}_k = \dfrac{\sum_i q_{ik}(\bm x_i - \bm\mu_k)(\bm x_i - \bm\mu_k)^T}{\sum_i q_{ik}}
            \end{equation*}
        \item Falls noch kein ausreichende Konvergenz erreicht ist, ab Schritt 2 wiederholen
    \end{enumerate} 
\end{mybox}
Aus \cref{sub:GMMs} ist bekannt,
dass die \glslink{marginal}{marginal} Log-Likelihood nicht gut für den \dref{cha:Gradient Descent} geeignet ist.
Das Problem ist vor allem deshalb schwer zu lösen,
da nicht bekannt ist,
zu welcher Teil-Komponente des \gls{GMM} ein einzelner Datenpunkt gehört.
Hierbei stellen die Datenpunkte $\bm x$ die bekannten Variablen (observed variables) dar
und die Zugehörigkeit der Datenpunkte zu den jeweiligen Teilkomponenten die latenten Variablen (latent variables) dar.
Daher ist jedes Mixture Model auch ein \gls{latent_variable_model},
auf welches das Verfahren des \glsxtrlong{EM} angewandt werden kann.
Wäre $p(\bm x_i,k|\bm\theta)$ bekannt wäre es sehr einfach eine Maximum Likelihood Abschätzung (\cref{sec:MLE}) durchzuführen:
\begin{itemize}
    \item $q_{ik} = \mathbb{I}(k,k_i) = \begin{cases} 1, &\text{ wenn der $i$-te Datenpunkt zur $k$-ten Komponente gehört }\\ 0, &\text{ sonst } \end{cases}$
    \item Koeffizienten:\tabto{3cm}$\pi_k = \dfrac{\sum_i q_{ik}}{N}$
    \item Durchschnitte:\tabto{3cm}$\bm\mu_k = \dfrac{\sum_i q_{ik}\bm x_i}{\sum_i q_{ik}}$
    \item Kovarianzen:\tabto{3cm}$\nomeq{covariance}_k = \dfrac{\sum_i q_{ik}(\bm x_i - \bm\mu_k)(\bm x_i - \bm\mu_k)^T}{\sum_i q_{ik}}$
\end{itemize}
Die Zugehörigkeit (cluster probability\slash\,responsibilites) $q_{ik}$ kann wiederum mithilfe des \gls{GMM} errechnet werden.
\begin{equation} \label{eq:responsibilites_GMM}
    q_{ik} = p(k_i=k|\bm x_i) = \frac{p(\bm x_i|k)p(k)}{\sum_j p(\bm x_i|j)p(j)} 
        = \frac{\pi_k\nomeq{gaussian_distribution}(\bm x_i|\bm\mu_k,\nomeq{covariance}_k)}{\sum_{j=1}^K \pi_j\nomeq{gaussian_distribution}(\bm x_i|\bm\mu_j,\nomeq{covariance}_j)}
\end{equation}
Hierdurch ergibt sich allerdings ein Henne"~Ei"~Problem.
Für die Ermittlung des \gls{GMM} brauchen wir die Responsibility $q_{ik}$ und umgekehrt.
Um dieses Problem zu umgehen wird das \gls{GMM} zunächst auf Basis einer Approximation (z.B. mittels k-means-Algorithmus) initialisiert
und anschließend der Expectation-Step und Maximization-Step solange wiederholt,
bis das \gls{GMM} konvergiert.
\begin{figure}[H]
    \centering
    \includegraphics[width=0.8\textwidth]{EM_for_GMM.png}
    \caption{Beispielhafte Anwendung von \texorpdfstring{\glsxtrlong{EM}}{\glsfmtlong{EM}} an einem \texorpdfstring{\glsxtrshort{GMM}}{\glsfmtshort{GMM}}}
    \label{fig:EM_for_GMM}
\end{figure}

\section{\glstopshort{EM} for Dimensionality Reduction}%
\label{sec:EM for Dimensionality Reduction}
\begin{mybox}
   \textbf{\Large \glsxtrshort{EM} for \nameref{cha:Dimensionality Reduction}} \\
   \begin{enumerate}
       \item \textbf{Initialisierung:} \noms{mean} von $\bm x$ für $\bm\mu$ verwenden, zufällige Matrix $\bm W$ generieren
       \item \textbf{Expectation-Step:} Posterior, \noms{mean} und \noms{covariance} errechnen
           \begin{equation} \label{eq:EM_for_dimensionality_reduction_expectation-step}
                \bm\mu_{\bm z|\bm x_i} = (\bm W^T\bm W + \nomeq{variance}\nomeq{identity_matrix})^{-1}\bm W^T(\bm x_i-\bm\mu),\qquad
                \nomeq{covariance}_{\bm z|\bm x_i} = \nomeq{variance}(\bm W^T\bm W + \nomeq{variance}\nomeq{identity_matrix})^{-1}
           \end{equation}
           \begin{itemize}
               \item $\bm z_i\sim\nomeq{gaussian_distribution}(\bm\mu_{\bm_z|\bm x_i},\nomeq{covariance}_{\bm q,\bm x_i})$
           \end{itemize}
       \item \textbf{Maximization-Step:} $\bm W$, $\bm\mu$ und \nomsym{variance} updaten
            \begin{align}
                \begin{bmatrix} \bm\mu\\\bm W \end{bmatrix} &= (\bm Z^T\bm Z)^{-1}\bm Z^T\bm X,\quad
                \text{ mit } \bm Z = \begin{bmatrix} 1 & \bm z_1^T \\ \vdots&\vdots \\ 1 & \bm z_n^T \end{bmatrix} 
                \text{ und } \bm X = \begin{bmatrix} \bm x_1^T \\\vdots\\ \bm x_n^T \end{bmatrix}
                \label{eq:EM_for_dimensionality_reduction_maximization-step_1} \\
                \nomeq{variance} &= \frac{1}{nd}\sum_{i=1}^n\sum_{k=1}^d(y_{ik} - x_{ik})^2,\quad\text{ mit } \bm y_i = \bm W\bm z_i + \bm\mu
                \label{eq:EM_for_dimensionality_reduction_maximization-step_2}
            \end{align}
   \end{enumerate}
\end{mybox}
Das \gls{EM} Verfahren kann auch zur \dref{cha:Dimensionality Reduction} verwendet werden.
Allerdings wird dieses Verfahren seit der Entdeckung der \dref{sub:PCA} kaum noch verwendet.

Die \dref{cha:Dimensionality Reduction} lässt sich als \gls{latent_variable_model} formulieren:
\begin{equation} \label{eq:dimensionality_reduction_latend_variable_model}
    \bm x = \bm W \bm z + \bm\mu + \bm\epsilon
\end{equation}
\begin{itemize}
    \item $\bm z$ ist die latente Variable (die Repräsentation in einer kleineren Dimension)
    \item $\bm W$ ist eine $D\times M$ Matrix, die eine Verbindung zwischen $\bm z\in\mathbb{R}^M$ und $\bm x\in\mathbb{R}^D$ herstellt ($M<D$)
    \item $\bm\mu$ ist ein konstanter Offset-Vektor
    \item $\bm\epsilon$ ist ein d-dimensionaler gaußscher Störvektor $\bm\epsilon \sim \nomeq{gaussian_distribution}(0,\nomeq{variance}\nomeq{identity_matrix})$
\end{itemize}
Die \glslink{marginal}{Marginal} Likelihood ist gegeben durch: 
\begin{equation} \label{eq:EM_for_dimensionality_reduction_marginal_likelihood}
    p(x|\bm\theta) = \int_{\bm z} p(\bm x|\bm z,\bm\theta)p(\bm z)d\bm z 
        = \int_{\bm z} \nomeq{gaussian_distribution}(\bm x|\bm W\bm z + \bm\mu,\nomeq{variance}\nomeq{identity_matrix})\nomeq{gaussian_distribution}(\bm z|0,\nomeq{identity_matrix}) d\bm z
\end{equation}
\begin{itemize}
    \item Sample Latent Variable $\bm z\sim\nomeq{gaussian_distribution}(0,\nomeq{identity_matrix})$
    \item Continuous Latent Variable $p(\bm z) = \nomeq{gaussian_distribution}(0,\nomeq{identity_matrix})$
    \item Observation Model $p(\bm x|\bm z,\bm\theta = \nomeq{gaussian_distribution}(\bm W\bm z) + \bm\mu,\nomeq{variance}\nomeq{identity_matrix})$
\end{itemize}

\subsection{Expectation-Step}%
\label{sub:EM for Dimensionality Reduction: Expectation-Step}
\includegraphics[scale=.6]{EM_for_Dimensionality_Reduction_expectation-step.png}

\subsection{Maximization-Step}%
\label{sub:EM for Dimensionality Reduction: Maximization-Step}
Die Lower-Bound kann in Bezug auf $\bm\theta$ bestimmt werden.
\begin{align} \label{eq:EM_for_dimensionality_reduction_lower_bound}
    \mathcal L(q,\bm\theta) &= \sum_i \left( \int_{\bm z} q_i(\bm z)\log p(\bm x_i,z|\bm\theta)dz - \int_{\bm z}q_i(\bm z)\log q_i(z)dz \right) \\
                            &= \dots\text{ Herleitung: {\color{red} Vorlesung 11 Folie 38} }\\
                            &= \sum_i \int_{\bm z} q_i(\bm z)\log p(\bm x_i|\bm z\bm\theta)d\bm z + \underbrace{\text{const}}_{\text{independent of $\bm\theta$}}
                            = \sum_i \mathbb{E}_{q_i(\bm z)}[\log p(\bm x_i|\bm z,\bm\theta)] + \underbrace{\text{const}}_{\text{independent of $\bm\theta$}}
\end{align}
Bei $\mathbb{E}_{q_i(\bm z)}$ handelt es sich um eine Approximation des Integrals auf Basis von Sample-Datenpunkten.
Hierfür wird zumeist die Monte-Carlo Abschätzung (\cref{sub:Monte-carlo estimation}) verwendet.
Hierbei ergibt sich (Herleitung: {\color{red} Vorelsung 11 Folie 40}):
\begin{align}
    \begin{bmatrix} \bm\mu\\\bm W \end{bmatrix} &= (\bm Z^T\bm Z)^{-1}\bm Z^T\bm X,\quad
    \text{ mit } \bm Z = \begin{bmatrix} 1 & \bm z_1^T \\ \vdots&\vdots \\ 1 & \bm z_n^T \end{bmatrix} 
    \text{ und } \bm X = \begin{bmatrix} \bm x_1^T \\\vdots\\ \bm x_n^T \end{bmatrix}
    \tag{\ref{eq:EM_for_dimensionality_reduction_maximization-step_1}} \\
    \nomeq{variance} &= \frac{1}{nd}\sum_{i=1}^n\sum_{k=1}^d(y_{ik} - x_{ik})^2,\quad\text{ mit } \bm y_i = \bm W\bm z_i + \bm\mu
    \tag{\ref{eq:EM_for_dimensionality_reduction_maximization-step_2}}
\end{align}

\subsection{Probabilistic \glstopshort{PCA} vs. \glstopshort{PCA}}%
\label{sub:Probabilistic PCA vs. PCA}
In den meisten Fällen ist das normale \gls{PCA} Verfahren (\cref{sub:PCA}) der Anwendung von \say{\gls{EM} for Dimensionality Reduction} (a.k.a. Probabilistic \gls{PCA}) vorzuziehen,
vor allem, da \gls{PCA} ein sehr recheneffizientes Verfahren ist.
Nur wenn eine Dichteabschätzung benötigt wird,
ist die Anwendung der der \say{\gls{EM} for Dimensionality Reduction} sinnvoll.
