'''
Programmation du réseau de neurones proposé dans le sujet CCINP PSI Modélisation 2024
'''

## Application du sujet

# Import

import numpy as np
from matplotlib import pyplot as plt
plt.close('all')

# Fonctions

def f(x):
    return 1/(1+np.exp(-x))

def f_prime(x):
    return np.exp(-x)/(1+np.exp(-x))**2
# = f(x)*(1-f(x))

# Fonction d'affichage

def Affiche(fig,X,Y,type):
    plt.figure(fig)
    plt.plot(X,Y,type)
    plt.grid(True)
    plt.show()
    plt.pause(0.01)

# Tracé de la fonction f sur [-20,20] et -[1,1]

LX = np.linspace(-20,20,1000)
LY = [f(t) for t in LX]
Affiche(1,LX,LY,'-')

LX = np.linspace(-1,1,1000)
LY = [f(t) for t in LX]
Affiche(2,LX,LY,'-')

# Création des matrices W1 W2 W3 B1 B2 B3

X = [[0],
     [1],
     [0],
     [0]]
X = np.array(X)

W1 = [[ 0.3, 0.2 , 0.25,-0.04],
      [ 0.4,-0.3 , 0.6 ,-0.3 ],
      [-0.4, 0.3 ,-0.6 ,-0.4 ],
      [-1  ,-0.75,-0.1 ,-0.5]]
W1 = np.array(W1)

B1 = [[0.0],
     [0],
     [0],
     [0]]
B1 = np.array(B1)

W2 = W1.copy()
B2 = B1.copy()
W3 = W1[0,:].copy().reshape(1,4)
B3 = np.array([[0.]])

'''
Veiller à créer des copies car les matrices seront modifiées en rétropropagation indépendemment les unes des autres
'''

# Inférence d'une couche

# Fonction inference_couche

def inference_couche(X,W,B):
    return f(np.dot(W,X)+B)

# Inférence de la première couche sur l'exemple proposé dans le sujet

LX = np.linspace(-1,1,1000)
LY = [f(t) for t in LX]
Affiche(3,LX,LY,'-')
LX = np.dot(W1,X)+B1
LY = [f(t) for t in LX]
Affiche(3,LX,LY,'or')

A1 = inference_couche(X,W1,B1)
print('A1:')
print(A1)

# Inférence du réseau de neurones

'''
Cette fonction évalue la valeur de delta (sortie 30 couche) pour une entrée de première couche (x1,x2,x3,x4)
'''

def inference(X):
    A1 = inference_couche(X,W1,B1)
    A2 = inference_couche(A1,W2,B2)
    A3 = inference_couche(A2,W3,B3)
    Delta = A3[0,0]
    return Delta

Delta = round(inference(X),2)
print(Delta)

# Données du sujet

Donnees = [
    [0, 0, 0, 1, -1],
    [0, 0, 1, 1, -1],
    [0, 0, 2, 0, -1],
    [0, 0, 0, 2, -1],
    [0, 0, 1, 2, -1],
    [0, 1, 0, 0, 0.18],
    [2, 2, 0, 0, 0.32],
    [2, 2, 2, 0, 0.49],
    [2, 2, 2, 2, -1],
    [0, 0, 1, 0, -1],
    [0, 0, 2, 1, -1],
    [0, 0, 2, 2, -1],
    [1, 0, 0, 0, 0.22],
    [2, 2, 1, 0, 0.42]
]

Donnees = np.array(Donnees)
Donnees_A = Donnees[:9]
Donnees_T = Donnees[9:]

# Mise en place de la phase d'entrainement

def retropropagation_couche(X,W,B,EA,alpha):
    Y = np.dot(W,X) + B
    fp = f_prime(Y)
    EW = np.dot(EA*fp,X.T)
    W -= alpha*EW # Modification en place
    EB = EA*fp
    B -= alpha*EB # Modification en place
    EX = np.dot(W.T,EA*fp)
    return EX

def epoch_ligne(X,Deltaopt,alpha):
    A1 = inference_couche(X,W1,B1)
    A2 = inference_couche(A1,W2,B2)
    A3 = inference_couche(A2,W3,B3) # = Delta
    EX3 = 2*(A3-Deltaopt) # Dérivée de l'erreur
    EX2 = retropropagation_couche(A2,W3,B3,EX3,alpha) # 3° couche
    EX1 = retropropagation_couche(A1,W2,B2,EX2,alpha) # 2° couche
    EX0 = retropropagation_couche(X,W1,B1,EX1,alpha) # 1° couche

def epoch(LD,alpha):
    N = len(LD)
    for i in range(N):
        X = LD[i,:4].reshape(-1,1) # Vertical
        Deltaopt = LD[i,4]
        epoch_ligne(X,Deltaopt,alpha)

def erreur_totale(LD):
    N = len(LD)
    S = 0
    for i in range(N):
        X = LD[i,:4].reshape(-1,1) # Vertical
        Deltaopt = LD[i,4]
        Delta = inference(X)
        Err = Deltaopt - Delta
        S += Err**2
    Err = S/N
    return Err

Err = erreur_totale(Donnees_A)
print("Erreur tot départ",Err)

def entrainement(LD,Nit,alpha):
    Lit = []
    Lerr = []
    for i in range(Nit):
        Lit.append(i)
        epoch(Donnees_A,alpha)
        err = erreur_totale(Donnees_A)
        Lerr.append(err)
    return Lit,Lerr

alpha = 0.01
Nit = 100
Lit,Lerr = entrainement(Donnees_A,Nit,alpha)
Affiche(4,Lit,Lerr,'-')

Err = erreur_totale(Donnees_A)
print("Ltot finale entraînement:",Err)

# Reconstruction de delta opt

def reconstruction(LD):
    '''Renvoie la table des données avec une colonne supplémentaire donnant la valeur de delta obtenue avec le réseau'''
    Nl,Nc = LD.shape
    Res = np.zeros([Nl,Nc+1])
    Res[:,:Nc] = LD
    for i in range(Nl):
        D = LD[i]
        X = D[:4].reshape(-1,1)
        Delta = inference(X)
        Res[i,Nc] = Delta
    return Res

Deltaopt_Rec = reconstruction(Donnees)
print(Deltaopt_Rec)

# Erreur sur les données de test

err = erreur_totale(Donnees_T)
print("Ltot finale test:",err)

## Application à un autre set de données avec des Wi et Bi aléatoires

'''
Les données du sujets ne mènent pas à la convergence du réseau
Testons un grand dataset différent dont les données sont en plus normalisées
'''

from random import random as rd

def cree_dataset(N):
    dataset = []
    for i in range(N):
        x1,x2,x3,x4 = rd(),rd(),rd(),rd()
        d = (x1+x2+x3+x4)/4
        data = [x1,x2,x3,x4,d]
        dataset.append(data)
    return np.array(dataset)

N = 1000
Ind = int((75/100)*N)
Donnees = cree_dataset(N)
Donnees_A = Donnees[:Ind]
Donnees_T = Donnees[Ind:]

W1 = np.random.rand(4,4)
W2 = np.random.rand(4,4)
W3 = np.random.rand(1,4)
B1 = np.random.rand(4,1)
B2 = np.random.rand(4,1)
B3 = np.random.rand(1,1)

Err = erreur_totale(Donnees_A)
print("Erreur tot départ",Err)

alpha = 0.01
Nit = 1000
Lit,Lerr = entrainement(Donnees_A,Nit,alpha)
Affiche(10,Lit,Lerr,'-')
# Première erreur de Lerr après une itération

Err = erreur_totale(Donnees_A)
print("Ltot finale entraînement:",Err)

Deltaopt_Rec = reconstruction(Donnees)
print(Deltaopt_Rec)

'''
En théorie, on trouve la droite y=x
Je trace à la fois les données d'apprentissage et de test
'''

Lx = np.linspace(0,1,2)
Affiche(11,Lx,Lx,'--')
Lx = Deltaopt_Rec[:,4]
Ly = Deltaopt_Rec[:,5]
Affiche(11,Lx,Ly,'o')

## Utilisation de sklearn

from sklearn.neural_network import MLPRegressor
from joblib import dump

rn = MLPRegressor(solver='sgd')
# rn = MLPRegressor(solver='sgd',hidden_layer_sizes=(4),activation='logistic')

'''
Il faut augmenter le nombre de données pour que le réseau converge
'''

N = 100000
Ind = int((75/100)*N)
Donnees = cree_dataset(N)
Donnees_A = Donnees[:Ind]
Donnees_T = Donnees[Ind:]

trn_x = [L[:4] for L in Donnees_A]
trn_y = [L[4] for L in Donnees_A]
tst_x = [L[:4] for L in Donnees_T]
tst_y = [L[4] for L in Donnees_T]

rn.fit(trn_x,trn_y)
dump(rn,'rn.joblib')

'''
from joblib import load
rn = load('rn.joblib')
'''

Score = rn.score(tst_x,tst_y)
Score = int(Score*10000)/100
print("Score:",Score,"%")

Lx = np.linspace(0,1,2)
Affiche(21,Lx,Lx,'--')
Lx = tst_y
Ly = [rn.predict([Lx])[0] for Lx in tst_x]
Affiche(21,Lx,Ly,'o')

# Analyse

def analyse():
    print("Dimensions du réseau:")
    print("Couche d'entrée:",rn.n_features_in_)
    nb_cachees = len(rn.hidden_layer_sizes)
    for i in range(nb_cachees):
        nb = rn.hidden_layer_sizes[i]
        print("Couches cachée",i+1,":",nb)
    print("Couche de sortie:",rn.n_outputs_)

analyse()