Introduction à l'apprentissage automatique - TP4 exercice 1

SVM sur données synthétiques


Prenez connaissance de la documentation scikit-learn sur les SVM: lisez l'introduction, l'introduction de 1.4, et survolez le reste pour savoir y revenir le jour où vous en avez besoin.

On utilisera la classe SVC.

Notez la valeur par défaut de l'hyperparamètre $C$, et les fonctions noyau disponibles nativement (ainsi que leurs paramètres).


Dans les questions suivantes, nous séparerons les bases de données entre bases d'apprentissage (80% de la base initiale) et base de test (20%) en utilisant model_selection.train_test_split, et nous comparerons les performances des différents classifieurs en calculant le score de classification sur la base de test.


On commence par charger les bibliothèques utiles et définir une fonction de visualisation que nous utiliserons par la suite:

In [ ]:
from sklearn import datasets, model_selection, preprocessing, model_selection, svm
import numpy as np
import matplotlib.pyplot as plt
from matplotlib.colors import ListedColormap
from matplotlib.colors import Normalize

%matplotlib inline

def plot_classif_result_SVM(X,y,clf,title):
    cmap_light = ListedColormap(['#FFAAAA', '#AAFFAA'])
    cmap_bold = ListedColormap(['#FF0000', '#00FF00'])    
    cmap2 = ListedColormap(['#FF8888', '#FFAAAA', '#AAFFAA', '#88FF88'])  
    
    h=0.01 # step size in the mesh
    
    x_min, x_max = X[:, 0].min() - 1, X[:, 0].max() + 1
    y_min, y_max = X[:, 1].min() - 1, X[:, 1].max() + 1
    xx, yy = np.meshgrid(np.arange(x_min, x_max, h),
                         np.arange(y_min, y_max, h))
    Zdf = clf.decision_function(np.c_[xx.ravel(), yy.ravel()])
    Z=Zdf>0
    Zdfbin = (np.abs(Zdf)<=1)   # 0 if inside margin,  1 if outside 
    Color=np.zeros(Z.shape)  # colors for each region 
    for i in range(len(Z)):
        if (Z[i]):
            if Zdfbin[i]: Color[i]=2
            else: Color[i]=3
        else:
            if Zdfbin[i]: Color[i]=1 
            else: Color[i]=0
                
    # Put the result into a color plot
    Z = Z.reshape(xx.shape)
    Color = Color.reshape(xx.shape)    
    plt.figure(figsize=[10,8])
    plt.pcolormesh(xx, yy, Color, cmap=cmap2)
    
    # Plot also the training points
    plt.scatter(X[:, 0], X[:, 1], c=y, cmap=cmap_bold,
                edgecolor='k', s=20)
    
    # Plot the support vectors (stars):
    plt.scatter(X[clf.support_, 0], X[clf.support_, 1], c=y[clf.support_], cmap=cmap_bold,edgecolor='k',s=80, marker='*')    
    
    plt.xlim(xx.min(), xx.max())
    plt.ylim(yy.min(), yy.max())
    plt.title(title);
    plt.axis("equal")
    plt.show();

1. Jeu de données "blobs"

Commençons par un jeu de données constitué de 1000 points du plan dans deux classes linéairement séparables, obtenu avec make_blobs (comme dans un des TP précédents).

In [ ]:
# génération dataset
# on précise random_state pour travailler tous sur le même jeu de données
# avec cluster_std=1.5, les classes sont linéairement séparables
X_dataset, y_dataset = datasets.make_blobs(n_features=2, centers=2, cluster_std=1.5, n_samples=1000, random_state=1)

X_train, X_test, y_train, y_test = model_selection.train_test_split(X_dataset,y_dataset,test_size=.2)

# affichage dataset train (gros points) + test (petits points)
plt.figure(figsize=[10,8])
plt.scatter(X_train[:, 0], X_train[:, 1], marker='o', c=y_train, s=25, edgecolor='k')
plt.scatter(X_test[:, 0], X_test[:, 1], marker='o', c=y_test, s=10)
plt.title('dataset')
plt.axis('equal')
plt.grid()
plt.show();

Ensuite, nous entraînons une SVM à noyau linéaire, et nous visualisons les vecteurs supports à l'aide de la fonction plot_classif_result_SVM.

On affiche également le score de classification sur la base test, ainsi que le nombre de vecteurs supports pour chaque classe.

Question 1. Les résultats dans le cas où les classes sont linéairement séparables vous semblent-ils cohérents avec le cours?

Que peut-on dire lorsque les deux classes ne sont plus linéairement séparables? (augmentez la valeur de cluster_std)

En particulier, comparez la marge, les vecteurs supports, leur nombre, le score sur la base de test, et le temps de calcul.

(pour obtenir deux classes non linéairement séparables, augmentez la valeur de cluster_std dans make_blobs dans la cellule précédente, et relancez la cellule suivante)

In [ ]:
SVM=svm.SVC(kernel='linear')
print('Apprentissage...\n')
%time SVM.fit(X_train,y_train)
print('\nAffichage...\n')
plot_classif_result_SVM(X_train,y_train,SVM,"SVM")
%time print("score test SVM %.3f" % SVM.score(X_test, y_test) )
print("nombre de vecteurs supports: %d pour classe 0 et %d pour classe 1" % (SVM.n_support_[0],SVM.n_support_[1]))

Réponse:

</font>

2. Jeu de données "moons":

Reprenons le jeu de données synthétique make_moons:

In [ ]:
# génération dataset
X_dataset, y_dataset = datasets.make_moons(noise=0.2, n_samples=1000)

X_train, X_test, y_train, y_test = model_selection.train_test_split(X_dataset,y_dataset,test_size=.2)

# affichage dataset train+test
plt.figure(figsize=[10,8])
plt.scatter(X_train[:, 0], X_train[:, 1], marker='o', c=y_train, s=25, edgecolor='k')
plt.scatter(X_test[:, 0], X_test[:, 1], marker='o', c=y_test, s=10)
plt.title('dataset')
plt.axis('equal')
plt.grid()
plt.show();

Noyau linéaire

La celulle suivante permet l'apprentissage par SVM linéaire, avec différentes valeurs de l'hyperparamètre $C$.

Question 2. Vérifiez que diminuer la valeur de l'hyperparamètre $C$ augmente le nombre de vecteurs supports (et change la surface de séparation). Est-ce cohérent avec la discussion du cours?

In [ ]:
SVM=svm.SVC(kernel='linear',C=100)
SVM.fit(X_train,y_train)
plot_classif_result_SVM(X_train,y_train,SVM,"SVM, C=100")
print("score test SVM %.3f" % SVM.score(X_test, y_test) )
print("nombre de vecteurs supports: %d pour classe 0 et %d pour classe 1" % (SVM.n_support_[0],SVM.n_support_[1]))

SVM=svm.SVC(kernel='linear')  # par défaut, C=1
SVM.fit(X_train,y_train)
plot_classif_result_SVM(X_train,y_train,SVM,"SVM, C=1")
print("score test SVM %.3f" % SVM.score(X_test, y_test) )
print("nombre de vecteurs supports: %d pour classe 0 et %d pour classe 1" % (SVM.n_support_[0],SVM.n_support_[1]))

SVM=svm.SVC(kernel='linear',C=0.01)
SVM.fit(X_train,y_train)
plot_classif_result_SVM(X_train,y_train,SVM,"SVM, C=0.01")
print("score test SVM %.3f" % SVM.score(X_test, y_test) )
print("nombre de vecteurs supports: %d pour classe 0 et %d pour classe 1" % (SVM.n_support_[0],SVM.n_support_[1]))

SVM=svm.SVC(kernel='linear',C=1e-4)
SVM.fit(X_train,y_train)
plot_classif_result_SVM(X_train,y_train,SVM,"SVM, C=0.0001")
print("score test SVM %.3f" % SVM.score(X_test, y_test) )
print("nombre de vecteurs supports: %d pour classe 0 et %d pour classe 1" % (SVM.n_support_[0],SVM.n_support_[1]))

Réponse:

</font>

Noyau RBF

Question 3 (observation plutôt que question). Observez le résultat de la classification avec un noyau gaussien (RBF), et l'évolution selon différentes valeurs du paramètre gamma (on gardera $C=1$, valeur par défaut de SVC).

In [ ]:
SVM=svm.SVC(kernel='rbf',gamma=0.01)
SVM.fit(X_train,y_train)
%time plot_classif_result_SVM(X_train,y_train,SVM,"SVM, gamma=0.01")
print("gamma=0.01, score test SVM %.3f" % SVM.score(X_test, y_test) )

SVM=svm.SVC(kernel='rbf',gamma=1)
SVM.fit(X_train,y_train)
%time plot_classif_result_SVM(X_train,y_train,SVM,"SVM, gamma=1")
print("gamma=1, score test SVM %.3f" % SVM.score(X_test, y_test) )

SVM=svm.SVC(kernel='rbf',gamma=100)
SVM.fit(X_train,y_train)
%time plot_classif_result_SVM(X_train,y_train,SVM,"SVM, gamma=100")
print("gamma=100, score test SVM %.3f" % SVM.score(X_test, y_test) )

Avec un noyau gaussien, $k(x,y)=\exp(-\gamma||x-y||^2)$ est différent de 0 seulement lorsque $x$ et $y$ sont proches relativement à $\gamma$ (disons $||x-y||<1/\sqrt{\gamma}$).

Le coefficient $\sqrt{\gamma}$ peut être vu comme l'inverse de la distance d'influence des vecteurs supports. En effet, le classifieur s'écrit: $$f(x) = \sum_i \lambda_i y_i k(x_i,x) + b$$

Lorsque $\gamma$ est grand, $k(x_i,x)$ a une contribution significative seulement pour $x$ proche d'un des vecteurs supports $x_i$. Ainsi, le signe de $f$ sera différent de celui de $b$ seulement pour $x$ proche d'un vecteur support $x_i$ associé à la classe de signe contraire à celui de $b$. On a intérêt à ce que toutes les observations soient des vecteurs supports de manière à minimiser le nombre de mauvaises classifications, donc la somme des variables d'écart. La surface de séparation est donc la superposition de disques autours des points d'une des classes, ce que l'on observe bien ici.

Lorsque $\gamma$ est petit, tous les $k(x_i,x)$ dans l'expression de $f$ ont une contribution. Le modèle est alors trop "moyenné" et on a une surface entre classes très régulière (presque une droite ici).

Justification au passage (en complément): si $\gamma$ est petit, on identifie l'exponentielle et son développement limité à l'ordre 1: $\exp(-\gamma||x-y||^2)=1-\gamma||x-y||^2$. Ensuite: $$f(x) = \sum_i \lambda_i y_i k(x_i,x) + b = \sum_i \lambda_i y_i (1-\gamma||x_i-x||^2) + b = \sum_i \lambda_i y_i (1-\gamma||x_i||^2-\gamma||x||^2-2\gamma x_i \cdot x) + b $$ Donc: $$f(x) = \sum_i \lambda_i y_i (1-\gamma ||x_i||^2) + b - \gamma ||x||^2 \sum_i \lambda_i y_i - 2\gamma \left(\sum_i \lambda_i y_i x_i\right) \cdot x$$ Comme $\sum_i \lambda_i y_i = 0$ (contrainte primale), $f(x)$ vérifie bien une relation du type $f(x) = B+ w\cdot x$, où $B=b-\sum_i \lambda_i y_i \gamma ||x_i||^2$ est une constante ne dépendant pas de $x$, et $w=- 2\gamma \left(\sum_i \lambda_i y_i x_i\right)$. Cela justifie que la frontière de séparation est linéaire!


Il y a sous-apprentissage avec $\gamma$ petit, et sur-apprentissage avec $\gamma$ grand.

Lorsqu'on ne spécifie pas $\gamma$ dans svm.SVC, une valeur est calculée à partir des observations. Voyez dans la documentation comment la valeur par défaut de $\gamma$ est adaptée en fonction des observations.

</font>

Question 4. Pour le noyau RBF et la valeur par défaut de $\gamma$, la cellule suivante présente différentes classifications selon les valeurs de $C$. Retrouvez les situations identifiées dans le cas linéaire.

Rappel : d'après la documentation : "The C parameter trades off misclassification of training examples against simplicity of the decision surface. A low C makes the decision surface smooth, while a high C aims at classifying all training examples correctly by giving the model freedom to select more samples as support vectors."

In [ ]:
SVM=svm.SVC(kernel='rbf',C=0.01)
SVM.fit(X_train,y_train)
plot_classif_result_SVM(X_train,y_train,SVM,"SVM, C=0.01")
print("C=0.01, score test SVM %.3f" % SVM.score(X_test, y_test) )

SVM=svm.SVC(kernel='rbf',C=1)  # valeur par défaut
SVM.fit(X_train,y_train)
plot_classif_result_SVM(X_train,y_train,SVM,"SVM, C=1")
print("C=1, score test SVM %.3f" % SVM.score(X_test, y_test) )

SVM=svm.SVC(kernel='rbf',C=100)  
SVM.fit(X_train,y_train)
plot_classif_result_SVM(X_train,y_train,SVM,"SVM, C=100")
print("C=100, score test SVM %.3f" % SVM.score(X_test, y_test) )

SVM=svm.SVC(kernel='rbf',C=10000)  
SVM.fit(X_train,y_train)
plot_classif_result_SVM(X_train,y_train,SVM,"SVM, C=10000")
print("C=10000, score test SVM %.3f" % SVM.score(X_test, y_test) )

Réponse:

</font>

La cellule suivante permet de trouver la valeur optimale de l'hyperparamètre $C$ par validation croisée à 5 plis (5-fold cross validation) sur la base d'apprentissage.

In [ ]:
# estimation de C par 5-fold cross validation 
accuracy=[]
arrayC=10**(np.arange(-1,5.5,.5))
for C in arrayC:
    SVM=svm.SVC(kernel='rbf',C=C)      
    scores = model_selection.cross_val_score(SVM, X_train, y_train, cv=5)
    accuracy.append(scores.mean())
    print('C=%.4f  score de validation croisée = %.4f +/- %.4f'%(C,scores.mean(),scores.std()))
plt.figure()
plt.semilogx(arrayC,accuracy,'-*')
plt.title('score de validation croisée contre C')
plt.grid()
plt.show();

Pour trouver une valeur optimale aux hyperparamètres $\gamma$ et $C$ par validation croisée, on dispose de la fonction GridSearchCV.

Cette fonction va calculer le score de validation croisée pour différentes valeurs des paramètres: ici, $C$ et $\gamma$ peuvent prendre des valeurs dans $\{10^{-3}, 10^{-2}, 10^{-1}, 1, 10, 100, 1000\}$. Pour gagner un peu sur le temps de calcul, on va se contenter de 3 plis (par défaut, cv=5 dans la dernière version de scikit-learn).

Visualisez les résultats de la cellule suivante:

In [ ]:
gamma_range=10**(np.arange(-3.,3.5,1))
C_range=10**(np.arange(-3.,3.5,1)) 
parameters = { 'gamma': gamma_range, 'C': C_range }
SVM = svm.SVC(kernel='rbf')
gridsearch=model_selection.GridSearchCV(SVM, parameters,cv=3)
gridsearch.fit(X_train,y_train)
print("Meilleur estimateur trouvé:")
print(gridsearch.best_estimator_)
print("Meilleurs paramètres:")
print(gridsearch.best_params_)

Question 5. Quel est le score de classification sur la base de test du meilleur classifieur identifié?

In [ ]:
# votre code ici

Question 6. (révision): Comparez au résultat de la classification de la base de test par les algorithmes de classification aux 1,5, 10 plus proches voisins, au classifieur de la régression logistique, et au classifieur naïf gaussien.

In [ ]:
from sklearn import  neighbors, naive_bayes, linear_model, neural_network

# votre code ici

réponse:

</font>

In [ ]: