import sys
import os
cwd = os.getcwd()
sys.path.insert(0, cwd+os.sep+".."+os.sep)

import histogrammes as hist

from corriges.gaussien import *

from corriges.melanges import *

import scipy

import pandas

input_file = ".." + os.sep + "PearsonLee.txt"

x = pandas.read_csv(input_file, skiprows=0, sep=" ", header=25)

noms_col = list(x.keys()) # Noms de colonnes

stature = scipy.array(x["parent"])

n = len(stature) # taille d'échantillon

hist.histogramme(stature, "histogramme des statures")
# Loi de melange gaussien pour estimer la variance

# theta = [[65., 5.], [66., 5.]]  # Valeur initiale des parametres
# Donne essentiellement une seule composante

theta = [[60., 5.], [70., 5.]]  # Valeur initiale des parametres

theta = [[60., 5.], [65., 5.], [70., 5.]]  # Valeur initiale des parametres

# Hyperparametres (loi a priori normale inverse-chi2)
# On prend ici les memes hyperparametres pour les deux gaussiennes
hyper_param = [0., 1., 2, 1.]

nit = 4000  # nombre d'iterations

K = len(theta)  # nombre d'etats du melange

## 2 etats
# Echantillonneur de Gibbs
mcmc2 = gibbs(nit, stature, theta, densite_norm, x_sim_norm,
              rpost_norm_invchi2, hyper_param)

# En sortie : labels des individus, loi predictive,
# liste des proportions estimees au cours des iterations
# liste des parametres d'emission estimes au cours des iterations
# tri permettant de compenser les permutations de labels

z, xpred, tri_prop, lambdap, tri = mcmc2

from matplotlib.pyplot import figure
import matplotlib.pyplot as plt

burn = 200

# Loi predictive a posteriori

fig = figure(figsize=(13, 7))
fig.add_subplot(1,2,1)
plt.hist(stature)
plt.title("histogramme des donnees")
fig.add_subplot(1,2,2)
plt.hist(xpred[burn:,0])
plt.title("histogramme de la loi predictive a posteriori")

# Proportions du melange

tri_prop[burn:,1].mean()

# Parametres d'emission pour l'etat 0

lambdap = scipy.array(lambdap)

fig = figure(figsize=(13, 7))

fig.add_subplot(1,2,1)

plt.hist(lambdap[burn:,0,0])

plt.title("Esperance a posteriori etat 0")

fig.add_subplot(1,2,2)

plt.hist(lambdap[burn:,0,1])

plt.title("Ecart-type a posteriori etat 0")


# Parametres d'emission pour l'etat 1

lambdap = scipy.array(lambdap)

fig = figure(figsize=(13, 7))

fig.add_subplot(1,2,1)

plt.hist(lambdap[burn:,1,0])

plt.title("Esperance a posteriori etat 1")

fig.add_subplot(1,2,2)

plt.hist(lambdap[burn:,1,1])

plt.title("Ecart-type a posteriori etat 1")

## 3 etats
theta = [[60., 5.], [65., 5.], [70., 5.]]  # Valeur initiale des parametres

K = len(theta)  # nombre d'etats du melange

nit = 30
burn = 10
# Echantillonneur de Gibbs
mcmc3 = gibbs(nit, stature, theta, densite_norm, x_sim_norm,
              rpost_norm_invchi2, hyper_param)

z, xpred, tri_prop, lambdap, tri = mcmc3

# Loi predictive a posteriori

fig = figure(figsize=(13, 7))
fig.add_subplot(1,2,1)
plt.hist(stature)
plt.title("histogramme des donnees")
fig.add_subplot(1,2,2)
plt.hist(xpred[burn:,0])
plt.title("histogramme de la loi predictive a posteriori")
