Interessons-nous aujourd'hui aux séquences d’ADN. Nous utiliserons le dataset téléchargeable ici : https://www.kaggle.com/datasets/nageshsingh/dna-sequence-dataset
L'ensemble des fichiers nécessaire à cet article sont disponibles ici.
Vous trouverez dans ce lot de données un ensemble de séquences d’ADN issues de 3 espèces : l’homme, le chien et le chimpanzé. Chacune de ces séquences appartient à une des 7 familles de séquences proposées. Ces 7 classes correspondent respectivement à :
- 0 : Récepteurs couplés aux protéines G
- 1 : Recepteurs-transmembranaires à activité tyrosine kinase
- 2 : Tyrosine phosphatase
- 3 : Asparagine synthetase
- 4 : ADN mitochondrial
- 5 : Récepteur vaniloïde
- 6 : Homéo-domaine TF1
L’idée sera donc ici de nous placer dans un contexte de classification à classes multiples, qui diffère de la classification binaire que nous avions réalisée ici : https://bioinfo-fr.net/bioinformatique-et-ia-un-premier-pas.
Nous aurons recours à un encodage des séquences afin qu’elles puissent être analysées par notre réseau. Je vous propose aujourd’hui d’utiliser la théorie du chaos (ou jeux du chaos) que nous développerons un peu plus loin dans cet article. Pour le moment, étudions nos données.
Étude préalable des données
Nous nous attarderons sur les séquences humaines qui sont au nombre de 4020, et que nous traiterons selon une représentation en 5‑mers (k=5), soit en matrices de 32x32 pixels.
Comme pour toute approche d’apprentissage nous divisons notre jeu de données en 3 sets : entraînement, test et validation selon le découpage 80%, 10% et 10% (pour plus d'information sur le découpage des données c'est ici : https://fr.wikipedia.org/wiki/Jeux_d%27entrainement,_de_validation_et_de_test). Nous obtiendrons après redimensionnement des tenseurs de taille 3216x32x32x1 pour le jeu d’entraînement, 402x32x32x1 pour les jeux de test et de validation avec une répartition stratifiée et respectueuse des proportions des différentes classes.
En nous intéressant à la répartition des classes de séquences à l’échelle du jeu de données global, on remarque une surreprésentation des séquences de classe 0 et une sous représentation des classes 3 et 6.
Le code nécessaire au chargement des données et à l'étude préalable des données est le suivant :
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 |
[crayon-672342a67a35b959593324 ]from collections import defaultdict import numpy as np import matplotlib.pyplot as plt import pandas as pd from sklearn.model_selection import train_test_split import tensorflow as tf import keras from tensorflow.keras.layers import Conv2D, Dense, Dropout, BatchNormalization, MaxPool2D from tensorflow.keras.models import load_model from keras.utils import to_categorical # SVM part from sklearn.datasets import make_classification from sklearn.svm import SVC from sklearn.metrics import accuracy_score from sklearn.metrics import ConfusionMatrixDisplay """# Mount drive""" from google.colab import drive drive.mount('/content/drive') """# Load data""" working_folder = "/content/drive/MyDrive/Sequences/" human_sequences = "/content/drive/MyDrive/Sequences/human.txt" dog_sequences = "/content/drive/MyDrive/Sequences/dog.txt" chimpanze_sequences = "/content/drive/MyDrive/Sequences/chimpanzee.txt" df_human = pd.read_csv(human_sequences, sep = "\t") df_dog = pd.read_csv(dog_sequences, sep = "\t") df_chimpanze = pd.read_csv(chimpanze_sequences, sep = "\t") """# Exploratory Data Analysis""" labels = [1, 2, 3, 4, 5, 6] plt.figure(figsize = (15, 10)) plt.subplot(131) plt.pie(df_human["class"].value_counts(sort = False).sort_index() , autopct='%1.1f%%') plt.title("Proportion of sequences classes in human") plt.legend(title="Class", labels=labels, loc="lower right",) plt.subplot(132) labels_dog = df_dog["class"].unique() plt.pie(df_dog["class"].value_counts(sort = False).sort_index() , autopct='%1.1f%%') plt.title("Proportion of sequences classes in dog") plt.legend(title="Class", labels=labels, loc="lower right",) plt.subplot(133) labels_chimpanze = df_chimpanze["class"].unique() plt.pie(df_chimpanze["class"].value_counts(sort = False).sort_index() , autopct='%1.1f%%') plt.title("Proportion of sequences classes in chimpanze") plt.legend(title="Class", labels=labels, loc="lower right",) plt.savefig(working_folder + "EDA.png") plt.show() |
La théorie du chaos
Très étudiée en physique depuis un paquet d’années, la théorie du chaos est étroitement liée à la création de fractales (https://fr.wikipedia.org/wiki/Fractale).
Dans notre cas nous allons créer des représentations matricielles carrées (graphiques) des séquences. Ces représentations à l'allure de fractales sont également étroitement liées à la composition en k‑mers de nos séquences. Voyons cela plus en détail.
Pour cela nous allons définir un carré dont chacun des sommets correspond à l’un des 4 fantastiques nucléotides qui font la richesse de notre métier, j’ai nommé A, T, C et G. Le centre du carré sera notre origine de repère soit (0, 0), ainsi nous disposerons les 4 nucléotides mentionnées aux coordonnées A (-1, 1), G (1, 1), T (1, ‑1) et enfin C (-1, ‑1).
L’idée est de se déplacer au sein de ce carré, à mesure que nous parcourons la séquence. Oui mais comment ? Et bien c’est très simple. En partant du centre du carré, nous allons nous déplacer jusqu’au point situé à mi-parcours sur le segment qui relie l’origine de notre repère (O(0, 0)) au premier nucléotide de la séquence. Prenons l’exemple où T serait le premier nucléotide. Alors nous nous déplacerons jusqu’au point médian entre (0, 0) et (1, ‑1), soit (0.5, ‑0.5). Appelons ce point Q puis dirigeons nous dans la direction du nucléotide suivant. Ce sera un G, alors nous suivrons le segment reliant Q (0.5, ‑0.5) à G(1, 1), pour atteindre le point médian du segment, soit le point R (0.75, 0.25). Le nucléotide d’après est un A. Nous parcourons alors le segment RA en direction de A (-1, 1) jusqu’au point médian S(-0.125, 0.625), et ainsi de suite…
Ainsi nous placerons un ensemble de points (autant que de nucléotides dans la séquence) dans le carré de départ. La représentation obtenue sera alors relative à l'abondance de l’ensemble des k‑mers ( à k préalablement défini).
Les représentations numériques des espaces n’étant pas continues mais discrètes, il convient de définir notre espace numérique, autrement dit la taille de la matrice de sortie. Pour cela, il convient de définir la taille k des oligomères que nous considérerons. Ainsi la matrice (ou espace) de sortie aura une cardinalité (ou aire) de 4^k éléments, soit un carré avec un côté de longueur (sqrt(4^k)).
Le parcours que nous avons vu précédemment permet ainsi de localiser les k‑mers au sein de cette matrice. En effet, si nous reprenons le chemin préalablement suivi nous nous dirigeons d’abord vers le cadran inférieur droit du carré pour le T, puis dans le cadran supérieur droit du cadran précédent pour le G, puis dans le cadran supérieur gauche de ce dernier cadran et ainsi de suite. À chaque nucléotide parcouru le cadran courant est alors re-divisé en 4 en conservant les orientations relatives aux sommets des nucléotides.
On observera alors que l’ensemble des séquences, et a fortiori des k‑mers commençant par TG se situeront dans le cadran supérieur droit du cadran inférieur droit. Selon cette découpe successive en 4 cadrans, il est alors possible de placer l’ensemble des k‑mers dans espace préalablement défini (ou éléments de la matrice i.e. pixels de l’image).
Pour un espace représentatif des 2‑mers et des 3‑mers, nous aurions les découpages suivants :
Très bien, mais quelles valeurs affecter à ces cadrans ou éléments de la matrice ? Et bien maintenant que l’ensemble des espaces dédiées aux k‑mers sont définis, nous affecterons aux éléments correspondants dans la matrice, la probabilité d’apparition de ces k‑mers dans la matrice. Pour plus d’informations sur les k‑mers et le calcul des probabilité d’occurence, je vous renvoie ici (https://fr.wikipedia.org/wiki/K‑m%C3%A8re).
━━
Ainsi, dans le cadre d’une représentation en 3‑mers, nous obtiendrons une matrice de dimension 8x8 dont les éléments auront une valeur comprise entre 0 et 1. Pour notre tutoriel ce sera des 5‑mers avec une représentation 32x32.
Pour en savoir plus sur les représentations de la théorie du chaos c’est ici : https://pubmed.ncbi.nlm.nih.gov/2336393/.
Le code nécessaire à la production de telles représentation est le suivant :
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 |
[crayon-672342a67a360288234853 ]"""## Sanitize data ### Human """ ## Delete sequences that contains N nb_sequences_to_drop = 0 rows_indexes_to_drop = list() for idx, seq in enumerate(df_human['sequence']): if 'N' in seq : nb_sequences_to_drop +=1 # display(df.loc[df['sequence'] == seq]) rows_indexes_to_drop.append(idx) print(f"Number of sequence to be dropped : {nb_sequences_to_drop}") print(f"Rows to be dropped : {rows_indexes_to_drop}") print(f"Shape of the data before dropping sequences containing Ns : {df_human.shape}") df_human.drop(rows_indexes_to_drop, inplace = True) print(f"Shape of the data after dropping sequences containing Ns : {df_human.shape}") """### Dog""" nb_sequences_to_drop = 0 rows_indexes_to_drop = list() for idx, seq in enumerate(df_dog['sequence']): if 'N' in seq : nb_sequences_to_drop +=1 rows_indexes_to_drop.append(idx) print(f"Number of sequence to be dropped : {nb_sequences_to_drop}") print(f"Rows to be dropped : {rows_indexes_to_drop}") print(f"Shape of the data before dropping sequences containing Ns : {df_dog.shape}") df_dog = df_dog.drop(rows_indexes_to_drop, inplace = True) """### Chimpanze""" nb_sequences_to_drop = 0 rows_indexes_to_drop = list() for idx, seq in enumerate(df_chimpanze['sequence']): if 'N' in seq : nb_sequences_to_drop +=1 rows_indexes_to_drop.append(idx) print(f"Number of sequence to be dropped : {nb_sequences_to_drop}") print(f"Rows to be dropped : {rows_indexes_to_drop}") print(f"Shape of the data before dropping sequences containing Ns : {df_human.shape}") df_chimpanze.drop(rows_indexes_to_drop, inplace = True) """## CGR transformation""" def cgr_encoding(sequence, k): """ sequence : Sequence to encode as CGR k = k‑mers size """ kmer_count = defaultdict(int) for i in range(len(sequence)-(k-1)): kmer_count[str(sequence[i :i+k])] +=1 for key in kmer_count.keys(): if "N" in key : del kmer_count[key] probabilities = defaultdict(float) N = len(sequence) for key, value in kmer_count.items(): probabilities[key] = float(value) / (N - k + 1) array_size = int(np.sqrt(4**k)) chaos = [] for i in range(array_size): chaos.append([0]*array_size) maxx = array_size maxy = array_size posx = 1 posy = 1 for key, value in probabilities.items(): for char in key : if char == "T": posx += maxx / 2 elif char == "C": posy += maxy / 2 elif char == "G": posx += maxx / 2 posy += maxy / 2 maxx = maxx / 2 maxy /= 2 chaos[int(posy)-1][int(posx)-1] = value maxx = array_size maxy = array_size posx = 1 posy = 1 return tf.constant(chaos) |
Voici ce à quoi ressemblent nos séquences (non, non il ne s’agit pas de la TV qui ne capte pas).
Le code nécessaire à la production de ces images est le suivant :
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 |
[crayon-672342a67a364215519049 ]"""## Plot data""" # Randomly selected one example per class random_idx_class_0 = np.random.choice(df_human[df_human["class"] == 0].index) random_idx_class_1 = np.random.choice(df_human[df_human["class"] == 1].index) random_idx_class_2 = np.random.choice(df_human[df_human["class"] == 2].index) random_idx_class_3 = np.random.choice(df_human[df_human["class"] == 3].index) random_idx_class_4 = np.random.choice(df_human[df_human["class"] == 4].index) random_idx_class_5 = np.random.choice(df_human[df_human["class"] == 5].index) # Plot some data plt.figure(figsize = (10, 10)) plt.subplot(231) plt.imshow(cgr_encoding(sequence = df_human["sequence"][random_idx_class_0], k = 5), cmap = "gray") plt.title("CGR of class 0") plt.subplot(232) plt.imshow(cgr_encoding(sequence = df_human["sequence"][random_idx_class_1], k = 5), cmap = "gray") plt.title("CGR of class 1") plt.subplot(233) plt.imshow(cgr_encoding(sequence = df_human["sequence"][random_idx_class_2], k = 5), cmap = "gray") plt.title("CGR of class 2") plt.subplot(234) plt.imshow(cgr_encoding(sequence = df_human["sequence"][random_idx_class_3], k = 5), cmap = "gray") plt.title("CGR of class 3") plt.subplot(235) plt.imshow(cgr_encoding(sequence = df_human["sequence"][random_idx_class_4], k = 5), cmap = "gray") plt.title("CGR of class 4") plt.subplot(236) plt.imshow(cgr_encoding(sequence = df_human["sequence"][random_idx_class_5], k = 5), cmap = "gray") plt.title("CGR of class 5") plt.subplots_adjust(wspace=0.12, hspace=0) plt.show() |
Il apparaît que les représentations des séquences d’ADN sont caractéristiques des espèces, des types de séquences, etc. Voyons alors comment utiliser ces représentations dans notre objectif de classification.
Pour ce faire, nous allons utiliser un petit réseau à 2 blocs convolutifs suivis d’une couche dense et d’une seconde couche dense de classification.
Afin d'entraîner notre réseau, il convient de préparer les données afin qu’elles puissent être fournies au réseau. Il convient aussi de diviser le jeu de donnée de départ en 3 sets : entrainement, test et validation.
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 |
[crayon-672342a67a368084546649 ]"""## Learning part ### Generate data """ encoded_list = [] for seq in df_human["sequence"]: encoded_list.append(cgr_encoding(seq, k =5)) X, y = np.array(encoded_list), to_categorical(df_human["class"]) """### Data spliting""" X = X.reshape(4020, 32, 32, 1) X_train, X_test, y_train, y_test = train_test_split(X, y, random_state = 42, stratify = y, shuffle = True, test_size = 0.2) # Split data into train and test with respect of class proportions X_test, X_validation, y_test, y_validation = train_test_split(X_test, y_test, random_state = 42, stratify = y_test, shuffle = True, test_size = 0.5) # Split data into test and validation with respect of class proportions |
Paramètres du réseau
Comme mentionné précédemment , nous alors utiliser un petit réseau de type convolutionnel. Contrairement au réseau utilisé la fois dernière, celui ci sera basé sur des bloc convolutionnels 2D, nénamoins le fonctionnement est similaire. Pour quelques rappels sur la convolution c’est ici : https://medium.com/@CharlesCrouspeyre/comment-les‑r%C3%A9seaux-de-neurones-%C3%A0-convolution-fonctionnent-b288519dbcf8.
Notre réseau comprendra donc deux blocs convolutionnels (une couche de convolution suivie d’une couche de pooling elle même suivie par une couche de dropout), lesquels seront suivis par une première couche dense contenant et aplatissant les filtres précédemment obtenus par convolution, et une couche dense à 7 neurones pour la classification.
Contrairement à la fois précédente, la couche de classification n’utilisera pas la fonction sigmoïde comme activation, mais la fonction softmax afin de rendre une probabilité d'appartenance à une classe donnée plutôt qu’une réponse de type “oui/non”.
Pour ce qui concerne les paramètres d’optimisation, nous utiliserons une descente de gradient (pour l’ajustement des poids du réseau) de type ADAM avec un taux d’apprentissage de 0.001. La fonction de perte sera la “categorical cross entropy” afin de tenir compte de l’erreur dans un contexte de classification à plusieurs classes.
Le code nécessaire à la construction du réseau et à son entrainement est le suivant :
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 |
[crayon-672342a67a36c295805156 ]model = tf.keras.Sequential() model.add(Conv2D(filters = 32, kernel_size = 3, activation = 'relu', padding = "valid", input_shape = (32, 32, 1), name = "Conv_1")) # Input shape : Batch_size, width, height, channels model.add(Dropout(0.3, name = "Dropout_1")) model.add(MaxPool2D(pool_size=4, strides=None, padding='valid', name = "MaxPool_1")) model.add(Conv2D(filters = 64, kernel_size = 2, activation = 'relu', padding = "valid", name = "Conv_2")) # Input shape : Batch_size, width, height, channels model.add(Dropout(0.3, name = "Dropout_2")) model.add(tf.keras.layers.Flatten(name = "Flatten")) model.add(tf.keras.layers.Dense(units = 7, activation = 'softmax', name = "Output")) model.summary() # Compile model model.compile(loss='categorical_crossentropy', optimizer=tf.keras.optimizers.Adam(learning_rate=0.001), metrics=['accuracy']) # Define early stopping to avoid overfitting early_stop = keras.callbacks.EarlyStopping(monitor = 'val_accuracy', min_delta = 0.0005, patience=8, restore_best_weights=True ) # Learning from data history = model.fit(X_train, y_train, batch_size = 8, validation_data= (X_test, y_test),epochs=200) |
Analyse des résultats
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 |
[crayon-672342a67a36f442850282 ]"""## Plot results""" plt.figure(figsize = (10,4)) plt.subplot(121) plt.plot(history.epoch, history.history["loss"], 'g', label='Training loss') plt.plot(history.epoch, history.history["val_loss"], 'b', label='Validation loss') plt.title('Training loss') plt.xlabel('Epochs') plt.ylabel('Loss') plt.legend() plt.subplot(122) plt.plot(history.epoch, history.history["accuracy"], 'r', label='Training accuracy') plt.plot(history.epoch, history.history["val_accuracy"], 'orange', label='Validation accuracy') plt.title('Training accuracy') plt.xlabel('Epochs') plt.ylabel('Accuracy') plt.legend() plt.savefig(working_folder + "base_modelpng") plt.show() """### Save model""" working_folder = "/content/drive/MyDrive/Sequences/" model.save(working_folder) """### Model evaluation""" loss, acc = model.evaluate(X_validation, y_validation, verbose=2) print("Basic model, accuracy : {:5.2f}%".format(100 * acc)) """## Plot results""" plt.figure(figsize = (10,4)) plt.subplot(121) plt.plot(history.epoch, history.history["loss"], 'g', label='Training loss') plt.plot(history.epoch, history.history["val_loss"], 'b', label='Validation loss') plt.title('Training loss') plt.xlabel('Epochs') plt.ylabel('Loss') plt.legend() plt.subplot(122) plt.plot(history.epoch, history.history["accuracy"], 'r', label='Training accuracy') plt.plot(history.epoch, history.history["val_accuracy"], 'orange', label='Validation accuracy') plt.title('Training accuracy') plt.xlabel('Epochs') plt.ylabel('Accuracy') plt.legend() plt.savefig(working_folder + "base_modelpng") plt.show() """### Save model""" working_folder = "/content/drive/MyDrive/Sequences/" model.save(working_folder) """### Model evaluation""" loss, acc = model.evaluate(X_validation, y_validation, verbose=2) print("Basic model, accuracy : {:5.2f}%".format(100 * acc)) |
Notre entraînement est terminé, que pouvons-nous en déduire ?
La qualité globale de l'entraînement est plutôt satisfaisante. Nous atteignons un début de plateau au niveau de la justesse, avec 0.9 au bout de 200 epochs (phases d’apprentissages) en entraînement. La justesse en validation est moindre, mais convenable. Notre réseau semble plutôt bien avoir appris. En évaluant le modèle avec le set de donnée de validation, nous obtenons une performance globale de 74% ce qui est relativement satisfaisant, mais pourrait-on faire mieux ? Voyons si le transfert learning (ou apprentissage transféré) pourrait nous aider.
Tranfert d'apprentissage
Le transfert learning représente l'utilisation d’un apprentissage préalablement réalisé afin de répondre à notre question. Pour autant nous pourrions être les premiers à nous poser cette question de classification de séquences, et encore plus si l’on considère notre encodage selon la théorie du chaos.
En pratique, le transfert learning est plutôt mis en place afin de répondre à des contraintes de temps et/ou de performances. Il s’agit de réutiliser des modèles semblables au nôtre et ayant déjà été entraînés et donnant des résultats satisfaisants. Vous me direz que tout cela est bien joli, mais comment trouver un modèle strictement identique au mien ? Et bien il n’est pas nécessaire que ce modèle soit strictement identique. Il faut simplement que les dimensions d’entrée du réseau soient compatibles avec nos données, ou alors, adapter nos données au réseau pré-entraîné.
Se pose alors la question de la question auquel répond notre réseau versus celle à laquelle répond le réseau pré-entrainé. Pour répondre à ces interrogations, je vous propose de plonger directement dans notre exemple.
Nous allons utiliser un réseau convolutionnel déjà entraîné sur le set de données MNIST et l’adapter à nos besoins. MNIST ? Vous l’avez sans doute déjà croisé il s’agit de ce célèbre jeu de données contenant des images de chiffres écrits à la main, célèbre dans tous les tutoriaux d’apprentissage profond. En principe, ces images sont de taille 28x28x1. Je me suis permis d’adapter ce réseau afin qu’il accepte nos représentations fractales de séquences de dimensions 32x32x1. Le réseau que nous allons utiliser est une adaptation du réseau LeNet5 proposé par Yann Lecun (plus d'informations).
Cependant, une différence majeure se présente. En effet, la couche de classification de ce réseau comporte 10 neurones afin de rendre compte de l’appartenance de l’une des images de chiffres à une classe entre 0 et 9. Pas de panique nous allons adapter cela à nos besoins.
Voyons d’abord dans les grandes lignes ce que nous allons faire. Tout d’abord il va s’agir de recharger le modèle déjà entraîné. (pour plus d’informations sur le modèle et sa conception, je vous renvoie vers ce notebook : https://github.com/bioinfo-fr/data_ia/blob/main/article_ia_2/LeNet5.ipynb).
Nous souhaitons conserver l’ensemble des couches du réseau, mais aussi profiter des poids déjà pré-établis par l'entraînement précédent. Pour autant le réseau n’est pas prêt à l’emploi pour notre problématique. Il va s’agir de “fine-tuner” ce modèle afin de tirer partie des descripteurs (filtres) établis lors de l'entraînement sur la base de données MNIST, en n'entraînant (ajustant les poids) que d’un sous-ensemble de couches de notre réseau. Dans notre exemple, le LeNet étant un modèle plutôt petit, nous ne conserverons que les filtres de la première couche. Il va alors falloir “geler” les poids (coefficients des filtres convolutifs) de la première couche et laisser le réseau n’ajuster que les autres poids.
Par ailleurs, nous avons évoqué que le LeNet5 permet une classification dans un contexte de 10 classes. Nous allons alors ôter la dernière couche du réseau dont le rôle est de réaliser cette classification, afin de la remplacer avec une couche dense de classification à 7 classes adaptée à notre problème, qui sera à entraîner. Nous allons garder les mêmes paramètres d’optimisation que ceux de notre premier réseau : ADAM avec un taux d’apprentissage de 0.001, fonction de perte de type “categorical cross entropy”.
Le code nécessaire à la réalisation de ces opérations est le suivant :
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 |
[crayon-672342a67a376443557576 ]"""## Load LeNet5""" lenet_folder = "/content/drive/MyDrive/LeNet/" lenet_model = load_model(lenet_folder) """### Check LeNet model""" lenet_model.summary() """### Freeze bottom layers for transfer learning""" # Freeze the first layers for idx, layer in enumerate(lenet_model.layers): if idx <= 1 : layer.trainable = False """### Pop last dense layer and replace it for relevent goal""" # Remove last layer lenet_model.pop() # Replace it with a more appropriate lenet_model.add(tf.keras.layers.Dense(units = 7, activation = 'softmax', name = "Output")) """### Check modified model""" lenet_model.summary() """### Compile model and train""" lenet_model.compile(loss='categorical_crossentropy', optimizer=tf.keras.optimizers.Adam(learning_rate=0.001), metrics=['accuracy']) # Define early stopping to avoid overfitting early_stop = keras.callbacks.EarlyStopping(monitor = 'val_accuracy', min_delta = 0.0005, patience=8, restore_best_weights=True ) history_lenet = lenet_model.fit(X_train, y_train, batch_size = 8, validation_data= (X_test, y_test), epochs=200) transfer_learning_directory = "/content/drive/MyDrive/Sequences/transfert_learning/" lenet_model.save(transfer_learning_directory) |
Analyse des résultats
Voyons nos résultats :
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 |
[crayon-672342a67a37a382207508 ]"""# Plot final results""" plt.figure(figsize = (14,4)) plt.subplot(121) plt.plot(history.epoch, history.history["loss"], 'g', linestyle = 'dotted', label='Training — LeNet') plt.plot(history_lenet.epoch, history_lenet.history["loss"], 'g', label='Training — LeNet') plt.plot(history.epoch, history.history["val_loss"], 'b', label='Validation loss') plt.plot(history_lenet.epoch, history_lenet.history["val_loss"], 'b' , linestyle = 'dotted', label='Validation — LeNet') plt.title('Training loss') plt.xlabel('Epochs') plt.ylabel('Loss') plt.legend() plt.subplot(122) plt.plot(history.epoch, history.history["accuracy"], 'r', label='Training accuracy') plt.plot(history_lenet.epoch, history_lenet.history["accuracy"], 'r',linestyle = 'dotted', label='Training — LeNet') plt.plot(history.epoch, history.history["val_accuracy"], 'orange', label='Validation accuracy') plt.plot(history_lenet.epoch, history_lenet.history["val_accuracy"], 'orange', linestyle = 'dotted', label='Validation — LeNet') plt.title('Training accuracy') plt.xlabel('Epochs') plt.ylabel('Accuracy') plt.legend() plt.savefig(transfer_learning_directory + "dna_lenet_model.png") plt.show() """### Evaluate model""" loss, acc = lenet_model.evaluate(X_validation, y_validation, verbose=2) print("DNA LeNet model, accuracy : {:5.2f}%".format(100 * acc)) |
On note ici deux différences avec le modèle précédent. La première étant que la justesse à l’issue de nos 200 epochs est légèrement supérieure en utilisant notre DNA-LeNet. La justesse en validation est semblable. Cependant la valeur de la fonction de perte sur la fin de l'entraînement semble s’affoler pour ce même réseau.
Par ailleurs, on observe un petit plateau au début de l'entraînement de nos deux modèles. Toutefois, le plateau dans le cadre de notre DNA-LeNet est plus court, le réseau se met plus rapidement à apprendre, mais apprend plus lentement que notre première version de réseau.
En évaluant ce nouveau modèle avec le set de données de validation, nous obtenons une performance globale de 72%. Il semble que changer de modèle ne nous ait que peu aidé.
Peut être que cela est imputable à notre couche de classification qui peine à déduire la bonne classe de l’ensemble des descripteurs fournis par les couches convolutives ? Voyons si nous pouvons classer nos séquences autrement que via cette couche.
Couplage apprentissage profond et apprentissage machine
À la différence de l’apprentissage machine, l’apprentissage profond présente l’avantage de nous affranchir partiellement de l’extraction de descripteurs de nos données. En effet, dans le cadre d’un réseau convolutif, ces derniers sont appris. Il s’agit des filtres de convolution dont les coefficients sont ajustés tout au long de l'entraînement.
Pourrait-on alors ne garder que ces descripteurs, et les fournir à un algorithme de classification plus conventionnel ? C’est ce que nous allons voir.
Nous utiliserons la SVM (Machine à Support de Vecteur). Il s’agit d’un algorithme de classification binaire dont le principe général est assez simple. Il va s’agir de placer nos données dans un espace vectoriel de grande dimension afin d'établir l’équation d’un hyperplan, lequel représentera la frontière séparant au mieux les différentes classes liées au différents points dans cet espace. Autrement dit dans notre cas, séparant l’ensemble des représentations des différentes classes de séquence, en maximisant la marge entre ces représentations à la marge.
Puisqu’il s’agit d’un hyperplan séparateur, nous ne pourrons que séparer deux classes à la fois. Néanmoins les SVMs offrent une stratégie appelée “one vs one” qui permet de construire un système de classification à plusieurs classes. La séparation de classe va donc s’effectuer deux à deux dans notre espace cherchant à d’abord séparer la classe 0 de la classe 1, puis la classe 0 de la classe 2, puis la classe 0 de la classe 3 et ainsi de suite. Nous allons alors entraîner autant de classifieurs binaires qu’il existe de combinaisons de classes.
Une fois ces classifieurs entraînés, l’étape de classification se fera ainsi. Une des représentations sera placée dans l’espace de classification. Une classe sera attribuée selon chacun des hyperplans calculés. Pour une représentation donnée, nous aurons donc un ensemble de classes plausibles. La classe prédite sera alors déterminée à la majorité. En d’autres termes, si la classe 0 est prédite 9 fois, la classe 1 elle 3 fois et la classe 5 à son tour 5 fois, alors la classe 0 sera affectée à notre représentation.
Trêve de théorie, passons à la pratique.
En premier lieu il s’agit de transformer notre réseau convolutionnel initial en extracteur de descripteurs. Pour cela nous allons retirer la couche de classification. Ainsi la sortie du réseau sera la sortie de la couche “Flatten” soit un vecteur de dimension 2304 par séquence.
1 2 3 4 5 6 7 8 9 10 |
[crayon-672342a67a37e581874438 ]"""## Change classification layer for SVM ### Adapt intial CNN for SVM classification """ # Chop off head of CNN modified_model = model modified_model.pop() modified_model.summary() |
L’étape suivante consiste en le ré-encodage des données de sorties. On va alors fournir l’ensemble de nos représentations CGR à notre réseau afin d’obtenir les vecteurs correspondants, et ce pour les sets d'entraînement, de test et de validation :
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 |
[crayon-672342a67a381771062857 ]"""### Prepare data for svm""" # Go back on labels from one hot encoding X_train_svm = modified_model.predict(X_train) y_train_svm = np.argmax(y_train, axis=1) X_test_svm = modified_model.predict(X_test) y_test_svm = np.argmax(y_test, axis=1) X_validation_svm = modified_model.predict(X_validation) y_validation_svm = np.argmax(y_validation, axis=1) """### Checking dimensions""" print(f"Shape of x_train_svm : {X_train_svm.shape}") print(f"Shape of y_train_svm : {y_train_svm.shape}") print(f"Shape of x_test_svm : {X_test_svm.shape}") print(f"Shape of y_test_svm : {y_test_svm.shape}") print(f"Shape of x_validation_svm : {X_validation_svm.shape}") print(f"Shape of y_validation_svm : {y_validation_svm.shape}") |
Enfin nous allons définir notre classifieur, soit une SVM en mode “one vs one” et entraîner ce modèle :
1 2 3 |
[crayon-672342a67a385311406034 ]"""### Instanciate SVM classifier and train it""" SVM_classifier = SVC(decision_function_shape = 'ovo', class_weight = "balanced", C = 125.0) |
Notre modèle étant à présent entraîné en utilisant les descripteurs extraits par le le réseau convolutionnel initial modifié, examinons les résultats :
1 |
[crayon-672342a67a388467152598 ]SVM_classifier.fit(X_train_svm, y_train_svm,) |
Analyse des résultats
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 |
[crayon-672342a67a38b818450324 ]"""### Evaluate classification""" # Make predictions on test set y_pred_svm = SVM_classifier.predict(X_test_svm) accuracy = accuracy_score(y_test_svm, y_pred_svm) print(f"Acuracy using CNN coupled to SVM : {accuracy}") """## Confusion matrix""" class_names = np.arange(7) # Plot non-normalized confusion matrix titles_options = [ ("Confusion matrix, without normalization", None), ("Normalized confusion matrix", "true"), ] for title, normalize in titles_options : disp = ConfusionMatrixDisplay.from_estimator( SVM_classifier, X_test_svm, y_test_svm, display_labels=class_names, cmap=plt.cm.Blues, normalize=normalize, ) disp.ax_.set_title(title) print(title) print(disp.confusion_matrix) plt.savefig(working_folder + "svm_confusion_matrix.png") plt.show() |
Comme toute procédure d’apprentissage machine, il aura fallu paramétrer les hyper-paramètres du modèle. Ici, je ne me suis intéressé qu’au coefficient “C” qui est le coefficient de marge d’un classifieur type SVM.
Les résultats sont ici très intéressants. Nous obtenons une justesse de 83.8% contre environ 70–75% avec nos deux précédentes démarches.
En observant la matrice de confusion, on peut s'apercevoir que les classes 0, 1 et 6 sont bien séparées des autres, tandis que les classes 2 et 3 semblent être plus difficilement distinguables. Néanmoins notre système de classification couplant CNN et SVM semble efficace.
Conclusion
Au travers de ce tutoriel, nous aurons parcouru la chaîne standard d’une (mini) analyse basée sur l’apprentissage automatique.
La première étape aura consisté en l’import des données, une analyse exploratoire des données et leur encodage selon la théorie des jeux du chaos. Dans une deuxième étape nous aurons entraîné un petit réseau de convolution avec nos séquences transformées en images.
Dans une troisième étape, nous aurons tenté d’améliorer les performances de ce modèle, en utilisant un autre modèle du même acabit entraîné lui, sur des images de chiffres, afin d’aborder la notion de transfert d’apprentissage.
Il semblerait que les étapes de classification de ces deux réseaux précédents aient du mal à remplir leur rôle. Cela peut éventuellement être lié à l'entraînement pour lequel 200 epochs ne seraient pas suffisants. En effet, on peut observer sur les courbes d’apprentissage que le plateau n’est pas tout à fait atteint.
Enfin dans une dernière étape, nous aurons cherché à améliorer notre classifieur, en utilisant notre réseau convolutionnel initial en tant qu’extracteur de descripteurs et en réalisant la classification avec une machine à vecteur de support, abordant ainsi le couplage apprentissage profond et apprentissage machine.
L'ensemble des données, des scripts, notebooks et sauvegardes de modèles sont disponibles sur ce repo git.
N.B. : Si vous reproduisez ce tutoriel il se peut que vos résultats diffèrent de ceux présentés ici en raison de l'aspect stochastique des procédures d'apprentissage et découpe du jeu de donnée initial.
Un grand merci à nos relecteurs pour cet article : Léopold, Samuel, Azerin & ZaZo0o !
Laisser un commentaire