import time
import numpy as np

import sklearn
from sklearn.pipeline import Pipeline
from sklearn.metrics import mean_squared_error
from sklearn.model_selection import train_test_split

import matplotlib.pyplot as plt

###############################################################################
# 05-regression_polynomiale.py
# @title: Fondamentaux - Apprentissage par régression polynomiale
# @project: Mes scripts de ML
# @lang: fr
# @authors: Philippe Roy <philippe.roy@ac-grenoble.fr>
# @copyright: Copyright (C) 2023 Philippe Roy
# @license: GNU GPL
###############################################################################

###
# Commandes NumPy :
# - np.array : créer un tableau à partir d'une liste de listes
# - np.linspace : créer un tableau 1D de la valeur de début à la valeur de fin avec n valeurs
# - .reshape : réagencer le tableau avec le nombre de lignes et le nombre de colonnes
###

###
# Commandes Scikit-Learn :
# - sklearn.linear_model.LinearRegression : créer un modèle de régression linéaire (méthode des moindres carrés)
# - sklearn.preprocessing.PolynomialFeatures : ajouter les carrées, les cubes, ..., x**n
# - .fit : entrainement du modèle
# - .predict : prédiction du modèle
###

###############################################################################
# Initialisation
###############################################################################

# Init du temps
t_debut = time.time()

# Init des plots
fig = plt.figure(figsize=(15, 5))
fig.suptitle("Régression polynomiale")
donnees_ax = fig.add_subplot(121) # Observations : x1 et cibles : y
apts_ax = fig.add_subplot(122) # Courbes d'apprentissage

###############################################################################
# Observations
###############################################################################

# Observations d'apprentisage
m = 100 # Nombre d'observations
bg = 1 # Quantité du bruit gaussien
x1 = 6*np.random.rand(m, 1) -3 # Liste des observations x1
y = 0.5 * x1**2 + x1 + 2 + bg * np.random.rand(m, 1) # Liste des cibles y
donnees_ax.plot(x1, y, 'b.', label="Observations")

# Séparation du jeu d'entrainement et du jeu de validation
x1_train, x1_val, y_train, y_val = train_test_split(x1, y, test_size=0.2)

# Nouvelles observations
x1_new=np.linspace(-3, 3, 100).reshape(100, 1)

###############################################################################
# Phase d'apprentissage
###############################################################################

# Polynome de degré 2
poly_features_2 = sklearn.preprocessing.PolynomialFeatures(degree=2, include_bias=False) # Ajout des carrés
std_scaler_2 = sklearn.preprocessing.StandardScaler()
lin_reg_2 = sklearn.linear_model.LinearRegression() # Modèle régression linéaire
model_2 = Pipeline([("poly_features", poly_features_2),("std_scaler", std_scaler_2),("lin_reg", lin_reg_2)])
model_2.fit(x1, y) # Entrainement

# Polynome de degré 1
poly_features_1 = sklearn.preprocessing.PolynomialFeatures(degree=1, include_bias=False) # Ajout des carrés
std_scaler_1 = sklearn.preprocessing.StandardScaler()
lin_reg_1 = sklearn.linear_model.LinearRegression() # Modèle régression linéaire
model_1 = Pipeline([("poly_features", poly_features_1),("std_scaler", std_scaler_1),("lin_reg", lin_reg_1)])
model_1.fit(x1, y) # Entrainement

# Polynome de degré 30
poly_features_30 = sklearn.preprocessing.PolynomialFeatures(degree=30, include_bias=False) # Ajout des carrés
std_scaler_30 = sklearn.preprocessing.StandardScaler()
lin_reg_30 = sklearn.linear_model.LinearRegression() # Modèle régression linéaire
model_30 = Pipeline([("poly_features", poly_features_30),("std_scaler", std_scaler_30),("lin_reg", lin_reg_30)])
model_30.fit(x1, y) # Entrainement

###############################################################################
# Phase d'inférence
###############################################################################

y_predict_2=model_2.predict(x1_new) # Prédiction polynome de degrée 2
y_predict_1=model_1.predict(x1_new) # Prédiction polynome de degrée 1
y_predict_30=model_30.predict(x1_new) # Prédiction polynome de degrée 300

###############################################################################
# Courbes d'apprentissage (fonction de la taille de lot et du degrée du modèle)
###############################################################################

# Polynome de degré 2

train_errors_2, val_errors_2 = [], []
for m in range (1, len (x1_train)):
    model_2.fit(x1_train[:m], y_train[:m]) # Entrainement sur jeu d'entrainement de taille variable
    y_train_predict_2=model_2.predict(x1_train[:m]) # Prédiction sur jeu d'entrainement de taille variable
    y_val_predict_2=model_2.predict(x1_val) # Prédiction sur jeu de validation
    train_errors_2.append(mean_squared_error(y_train[:m], y_train_predict_2))
    val_errors_2.append(mean_squared_error(y_val, y_val_predict_2))

# Polynome de degré 1

train_errors_1, val_errors_1 = [], []
for m in range (1, len (x1_train)):
    model_1.fit(x1_train[:m], y_train[:m]) # Entrainement sur jeu d'entrainement de taille variable
    y_train_predict_1=model_1.predict(x1_train[:m]) # Prédiction sur jeu d'entrainement de taille variable
    y_val_predict_1=model_1.predict(x1_val) # Prédiction sur jeu de validation
    train_errors_1.append(mean_squared_error(y_train[:m], y_train_predict_1))
    val_errors_1.append(mean_squared_error(y_val, y_val_predict_1))

# Polynome de degré 30

train_errors_30, val_errors_30 = [], []
for m in range (1, len (x1_train)):
    model_30.fit(x1_train[:m], y_train[:m]) # Entrainement sur jeu d'entrainement de taille variable
    y_train_predict_30=model_30.predict(x1_train[:m]) # Prédiction sur jeu d'entrainement de taille variable
    y_val_predict_30=model_30.predict(x1_val) # Prédiction sur jeu de validation
    train_errors_30.append(mean_squared_error(y_train[:m], y_train_predict_30))
    val_errors_30.append(mean_squared_error(y_val, y_val_predict_30))

###############################################################################
# Résultats
###############################################################################

# Plot des donnéesx
donnees_ax.set_title("Données")
donnees_ax.plot(x1_new, y_predict_1, 'g-', label="Prédictions - degré 1")
donnees_ax.plot(x1_new, y_predict_2, 'r-', label="Prédictions - degré 2")
donnees_ax.plot(x1_new, y_predict_30, 'k-', label="Prédictions - degré 30")
donnees_ax.set_xlabel(r'$x_1$')
donnees_ax.set_ylabel(r'$y$', rotation=0)
donnees_ax.legend()
donnees_ax.set(xlim=(-3, 3), ylim=(0, 10))

# Plot des courbes d'apprentissage
apts_ax.set_title("Courbes d'apprentissage")

apts_ax.plot(np.sqrt(train_errors_1), 'g:', label="Entrainement - degré 1")
apts_ax.plot(np.sqrt(val_errors_1), 'g-', label="Validation - degré 1")

apts_ax.plot(np.sqrt(train_errors_2), 'r:', label="Entrainement - degré 2")
apts_ax.plot(np.sqrt(val_errors_2), 'r-', label="Validation - degré 2")

apts_ax.plot(np.sqrt(train_errors_30), 'k:', label="Entrainement - degré 30")
apts_ax.plot(np.sqrt(val_errors_30), 'k-', label="Validation - degré 30")

apts_ax.set_xlabel("Taille du lot d'entrainement")
apts_ax.set_ylabel("RMSE")
apts_ax.set(xlim=(0, 80), ylim=(0, 10))
# apts_ax.set(xlim=(0, 80), ylim=(0, 3))
apts_ax.legend()

plt.show()

# Performances
print ("Temps : "+str(time.time()-t_debut))