%pylab inline

import numpy as np
from matplotlib import pyplot as plt

def pdf(x, mu, sigma_sqr):
    """
    Calculates the normal distribution's probability density 
    function (PDF).  
        
    """
    term1 = 1.0 / ( math.sqrt(2*np.pi * sigma_sqr))
    term2 = np.exp( (-1/2.0) * ( ((x-mu)**2 / sigma_sqr) ))
    return term1 * term2

# generating some sample data
x = np.arange(0, 50, 0.05)

def posterior(likelihood, prior):
    """
    Calculates the posterior probability (after Bayes Rule) without
    the scale factor p(x) (=evidence).  
        
    """
    return likelihood * prior

# probability density functions
posterior1 = posterior(pdf(x, mu=2, sigma_sqr=0.25), prior=2/3.0)
posterior2 = posterior(pdf(x, mu=1.5, sigma_sqr=0.04), prior=1/3.0)

# Class conditional densities (likelihoods)
plt.plot(x, posterior1)
plt.plot(x, posterior2)
plt.title('Posterior Probabilities w. Decision Boundary')
plt.ylabel('P(w)')
plt.xlabel('random variable x')
plt.legend(['P(w_1|x)', 'p(w_2|x)'], loc='upper right')
plt.ylim([0,1])
plt.xlim([0,4])
plt.axvline(1.08625, color='r', alpha=0.8, linestyle=':', linewidth=2)
plt.axvline(1.72328, color='r', alpha=0.8, linestyle=':', linewidth=2)
plt.annotate('R1', xy=(0.5, 0.8), xytext=(0.5, 0.8))
plt.annotate('R2', xy=(1.3, 0.8), xytext=(1.3, 0.8))
plt.annotate('R1', xy=(2.3, 0.8), xytext=(2.3, 0.8))
plt.show()

# Parameters
mu_1 = 2
mu_2 = 1.5
sigma_1_sqr = 0.25
sigma_2_sqr = 0.04

# Generating 10 random samples drawn from a Normal Distribution for class 1 & 2
x1_samples = sigma_1_sqr**0.5 * np.random.randn(10) + mu_1
x2_samples = sigma_1_sqr**0.5 * np.random.randn(10) + mu_2
y = [0 for i in range(10)]

# Plotting sample data with a decision boundary

plt.scatter(x1_samples, y, marker='o', color='green', s=40, alpha=0.5)
plt.scatter(x2_samples, y, marker='^', color='blue', s=40, alpha=0.5)
plt.title('Classifying random example data from 2 classes')
plt.ylabel('P(x)')
plt.xlabel('random variable x')
plt.legend(['w_1', 'w_2'], loc='upper right')
plt.ylim([-1,1])
plt.xlim([0,4])
plt.axvline(1.08625, color='r', alpha=0.8, linestyle=':', linewidth=2)
plt.axvline(1.72328, color='r', alpha=0.8, linestyle=':', linewidth=2)
plt.annotate('R1', xy=(0.5, 0.8), xytext=(0.5, 0.8))
plt.annotate('R2', xy=(1.3, 0.8), xytext=(1.3, 0.8))
plt.annotate('R1', xy=(2.3, 0.8), xytext=(2.3, 0.8))
plt.show()

w1_as_w2, w2_as_w1 = 0, 0
for x1,x2 in zip(x1_samples, x2_samples):
    if x1 > 1.08625 and x1 < 1.72328:
        w1_as_w2 += 1
    if x2 <= 1.08625 and x2 >= 1.72328:
        w2_as_w1 += 1
  
emp_err =  (w1_as_w2 + w2_as_w1) / float(len(x1_samples) + len(x2_samples))
    
print('Empirical Error: {}%'.format(emp_err * 100))