Exemple de classification traité avec R en cours d'apprentissage

Mélange de lois de probabilité

Dans cet exemple, nous allons générer un jeu de données comme suit :

1. 10 moyennes $m_k$ issues d'une loi normale bi-variée $\mathcal{N}((1,0),I)$ et étiquetées en BLEU.
2. 10 moyenne supplémentaires issues d'une loi normale bi-variée $\mathcal{N}((0,1),I)$ et étiquetées en ORANGE.
3. Pour chaque classe, on génère 100 observations comme suit:
   1. $m_k$ est tirée aléatoirement avec probabilité 1/10;
   2. l'observation est générée suivant une loi $\mathcal{N}(m_k,I/5)$, conduisant ainsi à un mélange de clusters gaussiens pour chaque classe.

Notre tâche consiste à créer une règle de classification qui détermine si un point est ORANGE ou BLEU à partir de ses coordonnées. Bien entendu, la règle de classification ne doit pas être basée sur la connaissance de la manière dont les données ont été générées.

# quelques déclarations : numpy et pyplot
import numpy as np
from matplotlib import pyplot as plt

%matplotlib inline

# numpy produit beaucoup "warnings", on va les ignorer 
np.warnings.filterwarnings('ignore')

# on peut déclarer des couleurs sous python (j'ai fait un copier-coller)
GRAY1, GRAY4, PURPLE = '#231F20', '#646369', '#A020F0'
BLUE, ORANGE, BLUE1 = '#57B5E8', '#E69E00', '#174A7E'

# on peut déclarer la police d'écriture et l'épaisseur des bords des figures
#plt.rcParams['font.family'] = 'Arial'
plt.rcParams['axes.linewidth'] = 0.5

Générer et préparer les données

Nous avons l'échantillon d'apprentissage qui a été utilisé en cours. De plus, les moyennes de tous les groupes BLUE et ORANGE sont connues. Les données de test sont inconnues et nous les générerons nous-mêmes.

# sklearn est la librairie qu'on va utiliser le plus dans ce cours
# le sous module mixture permet de gérer des mélanges de lois de probabilité
from sklearn.mixture import GaussianMixture
# le sous module metrics fournit un ensemble de mesure de distance : quadratique, absolue et taux de mauvais classement
from sklearn.metrics import accuracy_score
# attention pour générer des données gaussiennes multivariées, on a besoin d'inverser des matrices !!! 
from sklearn.mixture._gaussian_mixture import _compute_precision_cholesky
# on a besoin de pandas pour gérer des structures de données comme les dataframe, la lecture de fichiers etc ...
import pandas as pd
import os

# charger le jeu de données du cours 
df = pd.read_csv("mixture.txt")
print(df.head())
print(df.shape)
X_train = df[['x1', 'x2']].values
y_train = df.y.values

• Les moyennes bleu sont données par [-0.25343316, 1.7414788], [0.26669318, 0.3712341], [2.09646921, 1.2333642], [-0.06127272, -0.2086791], [2.70354085, 0.5968283], [2.37721198, -1.1864147], [1.05690759, -0.6838939], [0.57888354, -0.0683458], [0.62425213, 0.5987384], [1.67335495, -0.2893159]

• Les moyennes orange sont données par [1.19936869, 0.2484086], [-0.30256110, 0.9454190], [0.05727232, 2.4197271], [1.32932203, 0.8192260], [-0.07938424, 1.6138017], [3.50792673, 1.0529863], [1.61392290, 0.6717378], [1.00753570, 1.3683071],[-0.45462141, 1.0860697], [-1.79801805, 1.9297806]

1. Créer une matrice blue_means à l'aide de la fonction array de numpy où chaque linge correspond à une moyenne bleu.

2. Créer une matrice orange_means à l'aide de la fonction array de numpy où chaque linge correspond à une moyenne orange.

3. À l'aide de la fonction d'empilement vertical vstack de numpy, empiler les deux matrice dans une seule matrice all_means.

# création des trois matrices 
blue_means = np.array([[-0.25343316, 1.7414788], [0.26669318, 0.3712341],
                       [2.09646921, 1.2333642], [-0.06127272, -0.2086791],
                       [2.70354085, 0.5968283], [2.37721198, -1.1864147],
                       [1.05690759, -0.6838939], [0.57888354, -0.0683458],
                       [0.62425213, 0.5987384], [1.67335495, -0.2893159]])
orange_means = np.array([[1.19936869, 0.2484086], [-0.30256110, 0.9454190],
                         [0.05727232, 2.4197271], [1.32932203, 0.8192260],
                         [-0.07938424, 1.6138017], [3.50792673, 1.0529863],
                         [1.61392290, 0.6717378], [1.00753570, 1.3683071],
                         [-0.45462141, 1.0860697], [-1.79801805, 1.9297806]])
all_means = np.vstack((blue_means, orange_means))

# Nous allons déclarer la fonction la plus obscure à ce niveau de cours (avant d'avoir fait les lois gaussiennes multivariées)
# même si nous connaissons déjà les moyennes et les covariances, nous devons
# faire un "faux fit", sinon le modèle GaussianMixture ne fonctionnera pas (c'est python !!)
gaussian_mixture_model = GaussianMixture(
    n_components=20,
    covariance_type='spherical',
    means_init=all_means,
    random_state=1
).fit(all_means)
# déclarer les covariances connues 
gaussian_mixture_model.covariances_ = [1/5]*20
# cholesky permet de calculer les inverses des matrices de covariances 
gaussian_mixture_model.precisions_cholesky_ = _compute_precision_cholesky(
    gaussian_mixture_model.covariances_,
    gaussian_mixture_model.covariance_type)

On peut maintenant générer un jeu de données test

1. À l'aide du modèle gaussian_mixture_model, identifier la fonction appropriée pour générer un jeu de données test (X_test et y_test) de taille 10000.

2. Transformer y_test en variable binaire.

X_test, y_test = gaussian_mixture_model.sample(10000)

y_test = 1*(y_test>=10)

Règle de Bayes

Puisque nous connaissons la distribution à partir de laquelle nos données ont été échantillonnées, nous pouvons construire la règle Bayes optimale et calculer le taux d'erreur optimal appelé erreur de Bayes.

1. À l'aide de la fonction predict_proba de la classe gaussian_mixture_model, écrire une fonction optimal_bayes_predict de classification qui prend des échantillons en entrée et produit un vecteur de classes prédites. Calculer l'erreur de Bayes comme suit

def optimal_bayes_predict(X):
    components_proba = gaussian_mixture_model.predict_proba(X)
    # first 10 components are BLEU(0), and others are ORANGE(1)
    blue_proba = np.sum(components_proba[:, :10], axis=1)
    orange_proba = np.sum(components_proba[:, 10:], axis=1)
    y_hat = 1*(blue_proba < orange_proba)
    return y_hat

bayes_error_rate = 1 - accuracy_score(y_test, optimal_bayes_predict(X_test))
print(f'Erreur optimale de Bayes = {bayes_error_rate}')

Fonctions pour tracer des figures

Nous allons tracer de nombreuses figures similaires. Un petit cahier de charges :

• Cercles BLEU et ORANGE représentant les données d'entraînement.
• Petits points de fond BLEU et ORANGE avec la frontière de décision du modèle.
• La frontière de décision optimale de Bayes.
• Informations sur les performances du modèle et les taux d'erreur de test.

# préparer une figure avec un titre et les points du jeu de données d'entrainement
def plot_train_data(title):
    fig, ax = plt.subplots(figsize=(2.8, 2.8), dpi=110)
    ax.set_aspect(1.3)
    ax.scatter(X_train[:, 0], X_train[:, 1], s=18, facecolors='none',
               edgecolors=np.array([BLUE, ORANGE])[y_train])
    ax.tick_params(
        bottom=False, left=False, labelleft=False, labelbottom=False)
    ax.set_xlim(-2.6, 4.2)
    ax.set_ylim(-2.0, 2.9)
    fig.subplots_adjust(left=0, right=1, top=1, bottom=0)
    ax.text(-2.6, 3.2, title, color=GRAY4, fontsize=9)
    for spine in ax.spines.values():
        spine.set_color(GRAY1)
      
    return fig, ax

# test it
_, _ = plot_train_data('Training data')

1. Décortiquer les fonctions suivantes

def fill_prediction_grid(n1, n2, predict):
    x1, x2 = np.linspace(-2.6, 4.2, n1), np.linspace(-2.0, 2.9, n2)
    X = np.transpose([np.tile(x1, n2), np.repeat(x2, n1)])
    y = predict(X)
    return X, y

def fill_prediction_meshgrid(predict):
    n1, n2 = 1000, 1000
    X, y = fill_prediction_grid(n1, n2, predict)
    return X[:, 0].reshape(n1, n2), X[:, 1].reshape(n1, n2), y.reshape(n1, n2)

def plot_model(predict, title):
    fig, ax = plot_train_data(title)
    
    X, y = fill_prediction_grid(69, 99, predict)
    ax.scatter(X[:, 0], X[:, 1], marker='.', lw=0, s=2,
               c=np.array([BLUE, ORANGE])[y])
    
    X0, X1, Y = fill_prediction_meshgrid(predict)
    ax.contour(X0, X1, Y, [0.5], colors=GRAY1, linewidths=[0.7])
    return fig, ax

_, _ = plot_model(optimal_bayes_predict, 'Bayes Optimal Classifier')

# Sauvgarde de la règle de Bayes 
X0_bayes, X1_bayes, Y_bayes = fill_prediction_meshgrid(optimal_bayes_predict)


# Fonction qui ajoute les différentes performances sur une figure 
def plot_model_stat(predict, title):
    fig, ax = plot_model(predict, title)
    ax.contour(X0_bayes, X1_bayes, Y_bayes, [0.5], colors='purple',
               linewidths=[0.5], linestyles='dashed')
    test_error_rate = 1 - accuracy_score(y_test, predict(X_test))
    train_error_rate = 1 - accuracy_score(y_train, predict(X_train))
    parms = {'color': GRAY1, 'fontsize': 7,
             'bbox': {'facecolor': 'white', 'pad': 3, 'edgecolor': 'none'}}
    ax.text(-2.42, -1.35, f'Training Error: {train_error_rate:.3f}', **parms)
    ax.text(-2.42, -1.62, f'Test Error:       {test_error_rate:.3f}', **parms)
    ax.text(-2.42, -1.89, f'Bayes Error:    {bayes_error_rate:.3f}', **parms)
    return fig, ax

1. On a déclarer un ensemble de fonction plus ou moins complexes. Ce choix d'opter pour une programmation fonctionnelle est-il justifié ?

Classification par la règle des plus proches voisins

Nous avons besoin d'appeler un sous module de la librairie sklearn

from sklearn.neighbors import KNeighborsClassifier
from sklearn.model_selection import GridSearchCV

Modèle à 15 plus proches voisins

# ajuster le modèle
neighbors15_classifier = KNeighborsClassifier(n_neighbors=15).fit(X_train, y_train)
# tracer ses performances
_, _ = plot_model_stat(neighbors15_classifier.predict, '15-plus proches voisins')

1. Ajuster et tracer les performances de la règle de classification à 1-plus proche voisins

# ajuster le modèle
neighbors1_classifier = KNeighborsClassifier(n_neighbors=1).fit(X_train, y_train)
# tracer ses performances
_, _ = plot_model_stat(neighbors1_classifier.predict, '1-plus proches voisins')

1. Écrire une boucle for qui permet de remplir les listes train_errors et test_errors pour chaque valeur dans la liste ks. Penser à utiliser les fonctions append associée à un objet list et la fonction score associée à un object KNeighborsClassifier.

2. Identifier k_optim la valeur optimale qui minimise l'erreur de classification test ou qui maximise le taux de bon classement test.

On peut maintenant tracer les deux types d'erreurs sur la même figure comme suit

## déclaration tain et test 
ks = list(range(1, 50))
train_errors = list()
test_errors = list()
for k in ks:
    neighbors_k_classifier = KNeighborsClassifier(n_neighbors=k)
    neighbors_k_classifier.fit(X_train, y_train)
    train_errors.append(neighbors_k_classifier.score(X_train, y_train))
    test_errors.append(neighbors_k_classifier.score(X_test, y_test))
    
    
i_k_optim = np.argmax(test_errors)
k_optim = ks[i_k_optim]
print("Nombre de voisin optimal : %s" % k_optim)

plt.subplot(2, 1, 1)
plt.plot(ks, train_errors, label='Train')
plt.plot(ks, test_errors, label='Test')
plt.vlines(k_optim, plt.ylim()[0], np.max(test_errors), color='k',
           linewidth=3, label='Optimum on test')
plt.legend(loc='lower left')
plt.ylim([0, 1.2])
plt.xlabel('Nombre de voisins ')
plt.ylabel('Performance')

Validation croisée à 5 folds

k_neighbors_grid_search = GridSearchCV(
    KNeighborsClassifier(),
    {'n_neighbors': list(range(1, 50))},
    cv=5
).fit(X_train, y_train)
k_neighbors_grid_search.best_params_