# Exercise - 1D Gaussian Mixture Model and Expectation Maximization

## Introduction

In this notebook you will implement an EM algorithm to estimate the mean and variance of a gaussian mixture of two one-dimensional gaussians.

Remark: In order to detect errors in your own code, execute the notebook cells containing assert or assert_almost_equal. These statements raise exceptions, as long as the calculated result is not yet correct.

## Requirements

### Knowledge

To complete this exercise notebook, you should possess knowledge about the following topics:

• Maximum-Likelihood
• Notes by Christian Herta [HER18a]
• Expectation Maximization:
• Notes by Christian Herta [HER18b]
• Coursera Video [NOV18]

### Python Modules

import numpy as np
from matplotlib import pyplot as plt
from scipy.stats import norm

%matplotlib inline
np.random.seed(42)

## Exercises

### Data Generation

Given is a mixture of two one dimensional Gaussian distributions with mean$\mu_1$,$\mu_2$ and variance$\sigma^2_1$ and$\sigma^2_2$:

mu_1 = -2.5
sigma_square_1 = 4.
sigma_1 = np.sqrt(sigma_square_1)

mu_2 = 3.5
sigma_square_2 = 4.
sigma_2 = np.sqrt(sigma_square_2)
prob_1 = 0.4
size=40
def get_data(mu_1 = mu_1, sigma_1 = sigma_1,
mu_2 = mu_2, sigma_2 = sigma_2,
prob_1 = prob_1, size=size):

# latent variale - assignments
z = np.random.binomial(1, prob_1, size)

x_1 = np.random.normal(loc=mu_1, scale=sigma_1, size=size)
x_2 = np.random.normal(loc=mu_2, scale=sigma_2, size=size)
x = z * x_1 + (1-z) * x_2

return z, x
z, x = get_data()
x_1 = x[z==1]
x_2 = x[1-z==1]

plt.plot(x_1, np.zeros_like(x_1), 'rx')
plt.plot(x_2, np.zeros_like(x_2), 'bx')
plt.figure(figsize=(10,5))
plt.plot(x_1, np.zeros_like(x_1), 'rx')
plt.plot(x_2, np.zeros_like(x_2), 'bx')

rv_1 = norm(loc = mu_1, scale = sigma_1)
rv_2 = norm(loc = mu_2, scale = sigma_2)
x_ = np.arange(-6, 14, .1)
#plot the pdfs of the normal distributions
plt.plot(x_, prob_1 * rv_1.pdf(x_), "r-")
plt.plot(x_, (1-prob_1) * rv_2.pdf(x_), "b-")
plt.plot(x_, prob_1 * rv_1.pdf(x_) + (1-prob_1) * rv_2.pdf(x_), "g-")
plt.figure(figsize=(10,5))
plt.plot(x, np.zeros_like(x), 'gx')
plt.plot(x_, prob_1 * rv_1.pdf(x_) + (1-prob_1) * rv_2.pdf(x_), "g-")

### Exercise - Expectation Maximization

#### Recap

• E-Step: Compute the responsability of the data points$q(z_i=j) = p(z_i=j \mid x^{(i)}, \theta^{(t)})$
• M-Step: Maximize the Likelihood (MLE) w.r.t. the parameters, i.e. compute new$\theta^{(t+1)}$

The parameters to estimate here are mean and variance of one of the participating Gaussians.

for E-Step: $q(z_i=j) = p(z_i=j \mid x^{(i)}, \theta) = \frac{p(x^{(i)} \mid z_i=j, \theta) p (z_i=j \mid \theta)}{p(x^{(i)} \mid \theta)} = \frac{p(x^{(i)} \mid z_i=j, \theta) p(z_i=j \mid \theta)}{\sum_{j'} p(x^{(i)}, z_i=j' \mid \theta)} = \frac{p(x^{(i)} \mid z_i=j, \theta) p(z_i=j \mid \theta)}{\sum_{j'} p(x^{(i)} \mid z_i=j', \theta) p(z_i=j' \mid \theta)}$ $q(z_i=j)$ reflects the responsibility the$j$'th class (Gaussian) has for the$i$'th data point.

for M-Step: ${\mu_j}^{(t+1)} = \frac{ \sum_{i} x^{(i)} q(z_i=j)}{\sum_{i} q(z_i=j) }$

${\sigma_j}^{(t+1)} = \left( \frac{ \sum_{i} ({\mu_k}^{(t+1)}-x^{(i)})^2 q(z_i=j)}{\sum_{i} q(z_i=j) }\right)^{1/2}$

$p_j^{(t+1)} = \frac{\sum_i q(z_i=j)}{m}$

-$m$: number of training points

• Implement the EM-Algorithm for a mixture of two 1D-Gaussians. Apply the algorithm on the train set. Plot the result.

• Monitor the log-likelihood during the training. When your implementation is correct, it should increase after each step.

mu1=-0.5
mu2=0.5
sigma1=1
sigma2=1
prob1 = px1 = 0.5
px2 = 1-px1
px1, px2 = e_step_(x, mu1, sigma1, mu2, sigma2, prob1)
print (px1, px2)
mu1, sigma1, mu2, sigma2, prob1 = m_step(x, px1, px2)
print (mu1, sigma1, mu2, sigma2, prob1 )

## Literature

### Notebook License (CC-BY-SA 4.0)

The following license applies to the complete notebook, including code cells. It does however not apply to any referenced external media (e.g., images).

Exercise - 1D Gaussian Mixture Model and Expectation Maximization
by Christian Herta