Travaux pratiques - Réseaux de neurones multicouches (2)

(correspond à 1 séance de TP)

Cahier Jupyter

Références externes utiles :

• Documentation NumPy

• Documentation SciPy

• Documentation MatPlotLib

• Site scikit-learn

• Site langage python

L’objectif de cette séance de TP est d’examiner les algorithmes d’optimisation employés et les méthodes de régularisation disponibles.

Préparation des données

Les données employées sont les mêmes que dans la séance précédents (ensemble de 150 observations de fleurs d’iris). Vous pouvez les lire avec

import numpy as np
# valeurs séparées par virgule, colonnes 0 à 3, y compris 1ère ligne
data = np.loadtxt('iris.data', delimiter=',', usecols=[0,1,2,3], skiprows=0)

Les étiquettes étant sous forme textuelle dans le fichier iris.data (c’est à dire Iris-setosa, Iris-versicolor et respectivement Iris-virginica) nous tirons profit du fait que ces trois classes se succèdent dans le fichier de données et qu’il y a exactement 50 observations par classe pour générer directement des étiquettes de classe. Nous construirons des modèles permettant de séparer les iris en trois classes, Iris-setosa, Iris-versicolor et Iris-virginica. Nous pouvons obtenir les étiquettes de classe avec :

# générer étiquettes de classe
l1c = np.ones(50, dtype=int)
l2c = np.zeros(50, dtype=int)
target = np.concatenate((l1c, l2c, 2*l1c))

Dans un cas général il faudra utiliser la fonction LabelEncoder.

Nous pouvons maintenant séparer les ensebles d’apprentissage et de test, et visualiser les données :

np.random.seed(10)
from sklearn.model_selection import train_test_split
X_train1, X_test1, y_train1, y_test1 = train_test_split(data, target, test_size=0.4)
import matplotlib.pyplot as plt
plt.figure()
cmp = np.array(['r','g','b'])
plt.scatter(X_train1[:,2],X_train1[:,3],c=cmp[y_train1],s=50,edgecolors='none')
plt.scatter(X_test1[:,2],X_test1[:,3],c='none',s=50,edgecolors=cmp[y_test1])
plt.show()

Etude des algorithmes d’optimisation et de certains de leurs paramètres

Comparaison des algorithmes L-BFGS, SGD et Adam

Nous commençons par comparer les algorithmes d’optimisation proposés, L-BFGS, SGD et Adam, avec les paramètres par défaut.

Il faut noter que pour SGD (descente de gradient stochastique) la taille par défaut du mini-batch dans l’implémentation Scikit-learn (paramètre batch_size) est le minimum entre 200 et le nombre d’exemples. En conséquence, dans le cas des données Iris (150 exemples au total), avec les valeurs par défaut des paramètres, SGD est une descente de gradient par lot global (batch) et non stochastique.

from sklearn.neural_network import MLPClassifier
# L-BFGS :
clf = MLPClassifier(solver='lbfgs', max_iter=2000, random_state=100)
clf.fit(X_train1[:,2:4], y_train1)
print("Nombre itérations avec L-BFGS : {}".format(clf.n_iter_))
plt.figure()
plt.xlabel('itérations')
plt.ylabel('erreur')
#avec lbfgs, MLPClassifier n'a pas d'attribut loss_curve_
# SGD :
clf = MLPClassifier(solver='sgd', max_iter=2000, random_state=100)
clf.fit(X_train1[:,2:4], y_train1)
print("Nombre itérations avec SGD : {}".format(clf.n_iter_))
plt.plot(clf.loss_curve_, label='SGD')
# Adam :
clf = MLPClassifier(solver='adam', max_iter=2000, random_state=100)
clf.fit(X_train1[:,2:4], y_train1)
print("Nombre itérations avec Adam : {}".format(clf.n_iter_))
plt.plot(clf.loss_curve_, label='Adam')
plt.legend(loc='upper right')
plt.show()

Nous constatons que le nombre d’itérations est bien plus faible pour L-BFGS que pour les deux autres algorithmes, mais chaque itération est plus coûteuse. Avec L-BFGS l’évolution de l’erreur d’apprentissage n’est pas accessible à travers clf.loss_curve_. Ensuite, le nombre d’itérations est plus faible pour Adam que pour SGD. Adam est recommandé pour des ensembles de données volumineux, pour des cas plus simples L-BFGS est en général plus rapide (comme ici).

Effet de la vitesse d’apprentissage

En utilisant l’algorithme SGD, examinons l’effet de la variation de la vitesse d’apprentissage (initiale) en modifiant le paramètre learning_rate_init :

train_scores = []
test_scores = []
for lr_init in [0.1, 0.01, 0.001, 0.0001]:
    clf = MLPClassifier(solver='sgd', learning_rate_init=lr_init, max_iter=2000, random_state=100)
    clf.fit(X_train1[:,2:4], y_train1)
    train_scores.append(clf.score(X_train1[:,2:4], y_train1))
    test_scores.append(clf.score(X_test1[:,2:4], y_test1))

plt.figure()
plt.xlabel('log10(learning_rate_init)')
plt.ylabel('scores')
plt.plot(np.log10([0.1, 0.01, 0.001, 0.0001]), train_scores, 'go')
plt.plot(np.log10([0.1, 0.01, 0.001, 0.0001]), test_scores, 'r+')
plt.show()

Question :

Expliquez ce que vous constatez. Pour cela il peut être utile d’afficher les courbes d’apprentissage.

Correction :

Les résultats se détériorent avec une vitesse d’apprentissage trop élevée ou trop faible. Si la vitesse d’apprentissage est trop faible il faut augmenter le nombre d’itérations pour permettre la convergence. Si la vitesse d’apprentissage est trop élevée les modifications lors de certaines itérations sont trop fortes et l’algorithme diverge.

# avec une vitesse trop élevée
clf = MLPClassifier(solver='sgd', learning_rate_init=0.1, max_iter=2000, random_state=100)
clf.fit(X_train1[:,2:4], y_train1)
plt.figure()
plt.plot(clf.loss_curve_)
plt.show()

Question :

La vitesse d’apprentissage évolue-t-elle lors des itérations successives ? Pouvons-nous améliorer la convergence pour un cas où la vitesse d’apprentissage est initialisée à une valeur élevée (par ex. learning_rate_init=0.1) ?

Correction :

Non, la vitesse d’apprentissage reste constante, la valeur par défaut du paramètre learning_rate est 'constant'. En donnant à ce paramètre la valeur 'adaptive' on obtient de meilleurs résultats.

clf = MLPClassifier(solver='sgd', learning_rate_init=0.1, learning_rate='adaptive', max_iter=2000, random_state=100)
clf.fit(X_train1[:,2:4], y_train1)
plt.figure()
plt.plot(clf.loss_curve_)
plt.show()

Effet de l’inertie ou momentum

Question :

Avec une vitesse d’apprentissage fixée (par ex. à la valeur par défaut), examinez l’impact de différentes valeurs pour le paramètre momentum (ou inertie ; par ex. considérer les valeurs 0.1, 0.3, 0.7, 0.9) et tirez une conclusion.

Correction :

L’augmentation de l’inertie a un impact considérable sur l’amélioration des résultats et sur la réduction du nombre d’itérations avant convergence.

# momentum avec SGD
train_scores = []
test_scores = []
for mmt in [0.1, 0.3, 0.7, 0.9]:
    clf = MLPClassifier(solver='sgd', momentum=mmt, max_iter=6000, tol=0.00005, random_state=100)
    clf.fit(X_train1[:,2:4], y_train1)
    print("Nombre itérations avec momentum = {} : {}".format(mmt,clf.n_iter_))
    train_scores.append(clf.score(X_train1[:,2:4], y_train1))
    test_scores.append(clf.score(X_test1[:,2:4], y_test1))

plt.figure()
plt.xlabel('momentum')
plt.ylabel('scores')
plt.plot([0.1, 0.3, 0.7, 0.9], train_scores, 'go')
plt.plot([0.1, 0.3, 0.7, 0.9], test_scores, 'r+')
plt.show()

Etude des méthodes de régularisation disponibles

L’implémentation des perceptrons multi-couches dans Scikit-learn (MLPClassifier pour la classification) rend accessibles deux méthodes de régression:

• L’emploi d’un terme d’oubli (ou régularisation L2) dans la fonction minimisée, contrôlé par le paramètre alpha,

• L’arrêt précoce (early stopping) qui extrait une partie des données d’apprentissage pour les employer comme données de validation, contrôlé par le booléen early_stopping et par validation_fraction.

Nous examinons dans la suite l’impact de la régularisation.

Effet de la régularisation L2 (oubli)

Nous examinons l’effet de la régularisation L2 sur les mêmes données Iris déjà employées, avec l’algorithme d’optimisation SGD (descente de gradient stochastique). Nous affichons aussi des histogrammes des poids des connexions (les seuils sont laissés de côté mais vous pouvez les ajouter dans les histogrammes) et pour avoir accès à leurs valeurs nous devons réaliser d’abord une itération d’apprentissage. Les erreurs d’apprentissage et de test sont également affichées.

Pour commencer, nous faisons l’apprentissage sans régularisation, c’est à dire avec ``alpha``=0.0 (la valeur par défaut est 0.0001, pas assez faible pour ne pas avoir d’effet) :

from sklearn.neural_network import MLPClassifier
clf = MLPClassifier(hidden_layer_sizes=(200,), solver='sgd', max_iter=10000, learning_rate_init=0.01, alpha=0.0, random_state=100, tol=1e-6)
# afficher l'histogramme des poids après 1 itération
clf.partial_fit(X_train1[:,0:2], y_train1,np.unique(y_train1))
train_score = clf.score(X_train1[:,0:2], y_train1)
print("Le score en train est {}".format(train_score))
test_score = clf.score(X_test1[:,0:2], y_test1)
print("Le score en test est {}".format(test_score))# on ignore les biais, peuvent être obtenus avec clf.intercepts_
# forme des matrices affichées pour faciliter la compréhension de l'affichage de l'histogramme
print("Forme matrice poids entrée -> couche cachée : {}".format(clf.coefs_[0].shape))
print("Forme matrice poids couche cachée -> sorties: {}".format(clf.coefs_[1].shape))
plt.figure()
plt.hist(np.concatenate((clf.coefs_[0].reshape(400,),clf.coefs_[1].reshape(600,))),histtype='step')
plt.show()

L’histogramme montre que les poids sont initialisés suivant une loi uniforme sur un intervalle restreint autour de 0 (et une itération d’apprentissage ne les modifie pas beaucoup).

L’apprentissage est relancé ensuite jusqu’à un arrêt sur un critère d’absence de changement supérieur à tol (nous avons modifié la valeur qui était par défaut de 0.0001) :

clf.fit(X_train1[:,0:2], y_train1)
print("Nombre itérations avec SGD : {}".format(clf.n_iter_))
train_score = clf.score(X_train1[:,0:2], y_train1)
print("Le score en train est {}".format(train_score))
test_score = clf.score(X_test1[:,0:2], y_test1)
print("Le score en test est {}".format(test_score))
plt.figure()
plt.hist(np.concatenate((clf.coefs_[0].reshape(400,),clf.coefs_[1].reshape(600,))),histtype='step')
plt.show()

On constate un écart assez important entre l’erreur d’apprentissage et l’erreur de test, cette dernière étant sensiblement inférieure à la première, ce qui indique un sur-apprentissage. Il est possible d’afficher les courbes de probabilité = 0.5 pour chacun des trois neurones de sortie :

plt.figure()
# afficher les nuages de points d'apprentissage (pleins) et de test (creux)
plt.scatter(X_train1[:,0],X_train1[:,1],c=cmp[y_train1],s=50,edgecolors='none')
plt.scatter(X_test1[:,0],X_test1[:,1],c='none',s=50,edgecolors=cmp[y_test1])
# calculer la probabilité de sortie du MLP pour tous les points du plan
nx, ny = 200, 200
x_min, x_max = plt.xlim()
y_min, y_max = plt.ylim()
xx, yy = np.meshgrid(np.linspace(x_min, x_max, nx),np.linspace(y_min, y_max, ny))
Z = clf.predict_proba(np.c_[xx.ravel(), yy.ravel()])
Z0 = Z[:, 0].reshape(xx.shape)
Z1 = Z[:, 1].reshape(xx.shape)
Z2 = Z[:, 2].reshape(xx.shape)
# dessiner le contour correspondant à la frontière proba = 0,5
plt.contour(xx, yy, Z0, [0.5])
plt.contour(xx, yy, Z1, [0.5])
plt.contour(xx, yy, Z2, [0.5])
plt.show()

ainsi que l’évolution de l’erreur d’apprentissage :

plt.plot(clf.loss_curve_)
plt.show()

Question :

Refaites l’expérience en augmentant progressivement la valeur de alpha, par ex. à 0.0001, 0.01, 0.1, 2.0, 5.0, et examinez l’évolution des erreurs de test et d’apprentissage, l’histogramme des valeurs des poids des connexions ainsi que la forme des courbes de probabilités = 0.5. Que constatez-vous ?

Correction :

L’écart entre erreur d’apprentissage et erreur de test diminue avec l’augmentation de la régularisation (augmentation de alpha), au-delà d’une certaine valeur l’erreur d’apprentissage et l’erreur de test augmenteront car la régularisation sera excessive (les frontières deviendront presque linéaires). L’histogramme des valeurs des poids devient de plus en plus étroit car l’oubli a comme effet de pousser les valeurs des poids vers 0. Pour ces données les courbes de probabilités à 0.5 n’évoluent pas de façon significative avant la valeur de 5.0 pour alpha, qui produit une régularisation forte avec comme effet une linéarisation des frontières.

Effet de la régularisation par arrêt précoce

Pour examiner l’effet de l’arrêt précoce il est nécessaire de séparer un échantillon des données d’apprentissage à utiliser comme ensemble de validation. Les données des Iris étant déjà peu nombreuses, l’erreur sur un échantillon de validation nécessairement trop faible serait un mauvais indicateur du sur-apprentissage. Pour cette raison, nous nous servons d’un ensemble de données synthétique, généré suivant :

rot = np.array([[0.94, -0.34], [0.34, 0.94]])
sca = np.array([[4, 0], [0, 3]])
sca2 = np.array([[2, 0], [0, 1]])
sca3 = np.array([[1, 0], [0, 2]])
np.random.seed(10)
# générer données classe 1
c1d = (np.random.randn(200,2)).dot(sca).dot(rot)
# générer données classes 2 et 3
c2d1 = np.random.randn(50,2).dot(sca2)+[-10, 2]
c2d2 = np.random.randn(50,2).dot(sca2)+[-7, -2]
c2d3 = np.random.randn(50,2).dot(sca3)+[-2, -6]
c2d4 = np.random.randn(50,2).dot(sca3)+[5, -7]
data = np.concatenate((c1d, c2d1, c2d2, c2d3, c2d4))
# générer étiquettes de classe
l1c = np.zeros(200, dtype=int)
l2c = np.ones(50, dtype=int)
target = np.concatenate((l1c, l2c, 2*l2c, l2c, 2*l2c))
# séparation en ensemble d'apprentissage et ensemble de test
np.random.seed(100)
from sklearn.model_selection import train_test_split
X_train1, X_test1, y_train1, y_test1 = train_test_split(data, target, test_size=0.4)
# visualisation des données
cmp = np.array(['r','g','b'])
plt.figure()
plt.scatter(X_train1[:,0],X_train1[:,1],c=cmp[y_train1],s=20,edgecolors='none')
plt.scatter(X_test1[:,0],X_test1[:,1],c='none',s=20,edgecolors=cmp[y_test1])
plt.show()

Nous nous limiterons à des perceptrons multi-couches (PMC ou MLP) sans régularisation par oubli, donc avec alpha=0.0. Nous considérns d’abord, pour référence, un PMC ayant appris sans régularisation, donc avec early_stopping=False (valeur pas défaut, donc non précisé lors de l'appel) :

clfes = MLPClassifier(hidden_layer_sizes=(1000,), solver='sgd', max_iter=1000, learning_rate_init=0.05, alpha=0.0, random_state=100, tol=1e-6)
clfes.fit(X_train1[:,0:2], y_train1)
print("Forme matrice poids entrée -> couche cachée : {}".format(clfes.coefs_[0].shape))
print("Forme matrice poids couche cachée -> sorties: {}".format(clfes.coefs_[1].shape))
print("Nombre itérations avec SGD : {}".format(clfes.n_iter_))
train_score = clfes.score(X_train1[:,0:2], y_train1)
print("Le score en train est {}".format(train_score))
test_score = clfes.score(X_test1[:,0:2], y_test1)
print("Le score en test est {}".format(test_score))
plt.figure()
plt.hist(np.concatenate((clfes.coefs_[0].reshape(2000,),clfes.coefs_[1].reshape(3000,))),histtype='step')
plt.show()

L’erreur de test obtenue est plus faible que l’erreur d’apprentissage et l’écart est relaivement important. Un écart bien plus important peut être obtenu, en général, pour des données plus complexes (dimension supérieure, forme plus complexe de la frontière).

Nous pouvons afficher les courbes de probabilité = 0.5 pour chacun des trois neurones de sortie :

plt.figure()
# afficher les nuages de points d'apprentissage (pleins) et de test (creux)
plt.scatter(X_train1[:,0],X_train1[:,1],c=cmp[y_train1],s=20,edgecolors='none')
plt.scatter(X_test1[:,0],X_test1[:,1],c='none',s=20,edgecolors=cmp[y_test1])
# calculer la probabilité de sortie du MLP pour tous les points du plan
nx, ny = 200, 200
x_min, x_max = plt.xlim()
y_min, y_max = plt.ylim()
xx, yy = np.meshgrid(np.linspace(x_min, x_max, nx),np.linspace(y_min, y_max, ny))
Z = clfes.predict_proba(np.c_[xx.ravel(), yy.ravel()])
Z0 = Z[:, 0].reshape(xx.shape)
Z1 = Z[:, 1].reshape(xx.shape)
Z2 = Z[:, 2].reshape(xx.shape)
# dessiner le contour correspondant à la frontière proba = 0,5
plt.contour(xx, yy, Z0, [0.5])
plt.contour(xx, yy, Z1, [0.5])
plt.contour(xx, yy, Z2, [0.5])
plt.show()

Question :

Refaites l’expérience en introduisant l’utilisation de l’arrêt précoce avec early_stopping=True. Le paramètre validation_fraction devrait être mis à une valeur de 0.2 ou 0.3 car la valeur par défaut de 0.1 est trop faible. Que constatez-vous pour l’évolution des erreurs de test et d’apprentissage, l’histogramme des valeurs des poids des connexions ainsi que la forme des courbes de probabilités = 0.5 ?

Correction :

L’écart entre les erreurs d’apprentissage et de test diminue avec l’arrêt précoce (l’erreur d’apprentissage finale augmente, l’erreur de test diminue ; le score d’apprentissage final diminue, le score de test augmente). Le nombre d’itérations d’apprentissage diminue également, parfois de façon significative. L’histogramme des valeurs des poids des connexions devient plus étroit avec apprentissage précoce (en absence de terme d’oubli) car plus l’apprentissage de poursuit, plus les poids ont tendance à s’éloigner de 0 (en l’absence de terme d’oubli).

.. code-block:: python

   clfes = MLPClassifier(hidden_layer_sizes=(1000,), solver='sgd', max_iter=1000, learning_rate_init=0.05, early_stopping=True, validation_fraction=0.2, alpha=0.0, random_state=100, tol=1e-6)
   clfes.fit(X_train1[:,0:2], y_train1)
   print("Forme matrice poids entrée -> couche cachée : {}".format(clfes.coefs_[0].shape))
   print("Forme matrice poids couche cachée -> sorties: {}".format(clfes.coefs_[1].shape))
   print("Nombre itérations avec SGD : {}".format(clfes.n_iter_))
   train_score = clfes.score(X_train1[:,0:2], y_train1)
   print("Le score en train est {}".format(train_score))
   test_score = clfes.score(X_test1[:,0:2], y_test1)
   print("Le score en test est {}".format(test_score))
   plt.figure()
   plt.hist(np.concatenate((clfes.coefs_[0].reshape(2000,),clfes.coefs_[1].reshape(3000,))),histtype='step')
   plt.show()