TP++ - Analyse factorielle discriminante - Étape descriptive

(les TP++ sont des compléments à l’enseignement de base et ne font pas partie des séances de travaux pratiques avec encadrement)

(une variante de ce TP utilisant le langage Scala)

Références externes utiles :

L’objectif de cette séance de TP est de montrer une mise en œuvre avec Spark de l’étape descriptive de l’analyse factorielle discriminante (AFD) ou Linear Discriminant Analysis (LDA). Cette méthode d’analyse n’est pas (encore) proposée par Spark dans MLLib.

L’analyse factorielle discriminante traite des observations caractérisées par des variables quantitatives plus une variable nominale, de « classe », et cherche à trouver dans l’espace initial constitué par les variables quantitatives un sous-espace linéaire dans lequel la séparation entre classes est optimale. Une présentation rapide de l’AFD se trouve dans cette section de RCP216 et une présentation plus détaillée dans cette section de RCP208.

L’approche employée ici pour réaliser l’AFD consiste à

  • effectuer une ACP normée des observations,

  • projeter les observations et les centres des classes sur les axes principaux correspondant à des valeurs propres non nulles,

  • appliquer aux observations une transformation qui égalise les variances des projections sur les axes principaux (whitening),

  • calculer les centres de gravité des classes,

  • appliquer une nouvelle ACP aux seuls centres de gravité, les axes principaux résultants étant ainsi les axes discriminants de l’analyse discriminante des observations,

  • projeter les observations sur les trois premiers axes discriminants.

La mise en œuvre de l’AFD est illustrée ici sur les données « textures », voir par exemple cette présentation. Ces données consistent en 5500 observations caractérisées par 40 variables numériques de départ et appartenant à 11 classes.

L’approche proposée ici peut être adaptée à la visualisation d’autres données vectorielles qui sont caractérisées aussi par des étiquettes de classes (adaptées aux modèles décisionnels de classement) ou des étiquettes de groupes (résultant de l’application de modèles descriptifs issus de la classification automatique, par ex. k-means).

Pour préparer les répertoires et télécharger les données, ouvrez une fenêtre terminal et entrez les commandes suivantes :

%%bash
mkdir -p tpafd/data
wget -nc http://cedric.cnam.fr/vertigo/Cours/RCP216/docs/texture.csv -P tpafd/data/
pyspark

Sous Windows, créez le répertoire tpafd et le sous-répertoire data, téléchargez le fichier de données dans le sous-répertoire data, ouvrez une fenêtre invite de commandes et lancez pyspark.

Ensuite il suffit d’exécuter les instructions ci-dessous (dans pyspark ou dans le cahier Jupyter du TP++, suivant votre choix initial). Noter que la séquence d’opérations indiquée ci-dessous peut être optimisée.

# Etape 1 :
# - charger les données textures
# - trouver le centre de gravité de chacune des 11 classes

from pyspark.sql.types import DoubleType, StructType, StructField

texturesFields = [StructField("val"+str(i), DoubleType(), True) for i in range(0, 40)]
texturesSchema = StructType(texturesFields).add("label", DoubleType(), True)

# lire les données textures dans un DataFrame
texturesDF = spark.read.format("csv").option("sep"," ").schema(texturesSchema) \
                  .load("tpafd/data/texture.csv").cache()

# vérifier la lecture : il y a 500 données (lignes) pour chacune des 11 classes
texturesDF.groupBy("label").count().show()

# calculer le centre de gravité de chacune des 11 classes
from pyspark.sql.functions import col, avg
cgclassesDF = texturesDF.groupBy("label").avg()
print(cgclassesDF.columns)

# vérification
cgclassesDF.groupBy("label").count().show()

# Etape 2 :
# - assembler les colonnes 1D pour préparer la normalisation et l'ACP,
# - appliquer la normalisation et l'ACP
# - projeter les données et les centres des classes sur les composantes principales

from pyspark.ml.linalg import Vectors
from pyspark.ml.feature import VectorAssembler

# définir l'assembler pour les données
colsEntree = ["val"+str(i) for i in range(0, 40)]
assembleur = VectorAssembler().setInputCols(colsEntree).setOutputCol("features")

# "assembler" les données
texturesDFA = assembleur.transform(texturesDF)
# vérification
texturesDFA.printSchema()
texturesDFA.select("features").show()

# définir l'assembler pour les centres des classes
colsEntreeCG = ["avg(val"+str(i)+")" for i in range(0, 40)]
assembleurCG = VectorAssembler().setInputCols(colsEntreeCG).setOutputCol("features")

# "assembler" les centres des classes
cgclassesDFA = assembleurCG.transform(cgclassesDF).select("label","features")

from pyspark.ml.feature import StandardScaler
from pyspark.ml.feature import PCA

# définir la normalisation
scaler = StandardScaler().setInputCol("features") \
                         .setOutputCol("scaledFeatures") \
                         .setWithStd(True).setWithMean(True)
# définir l'ACP sur les données (ACP1)
pca = PCA().setInputCol("scaledFeatures").setOutputCol("pcaFeatures").setK(40)

# calculer la normalisation et inspecter
scalerModel = scaler.fit(texturesDFA)
scalerModel.mean

# appliquer la normalisation aux données
scaledTexturesDFA = scalerModel.transform(texturesDFA).select("label","scaledFeatures")
scaledTexturesDFA.printSchema()

# calculer l'ACP (ACP1) et inspecter
pcaModel = pca.fit(scaledTexturesDFA)
pcaModel.explainedVariance

# Etape 3 :
# - "blanchir" les centres des classes,
# - appliquer une ACP aux seuls centres
#   (les axes principaux résultants sont les axes discriminants des données initiales)

from pyspark.ml.feature import ElementwiseProduct
from pyspark.ml.linalg import Vectors
import numpy as np

# préparer le vecteur de "blanchiment" (whitening)
wv = Vectors.dense(1/(np.sqrt(pcaModel.explainedVariance)))

# définir le blanchiment
whitener = ElementwiseProduct().setScalingVec(wv) \
                               .setInputCol("pcaFeatures") \
                               .setOutputCol("whitenedFeatures")

# appliquer normalisation, projection ACP1 et blanchiment aux centres des classes
cgclassesDFAwhitened = whitener.transform(pcaModel.transform(scalerModel \
                                                  .transform(cgclassesDFA)))

# définir nouvelle ACP (ACP2) sur centers "blanchis"
#  the resulting PCs are the discriminant components of the initial data
pcaCG = PCA().setInputCol("whitenedFeatures").setOutputCol("ldaFeatures").setK(3)

# calculer nouvelle ACP (ACP2) et inspecter
ldaModel = pcaCG.fit(cgclassesDFAwhitened)
ldaModel.explainedVariance

# Etape 4:
# - appliquer aux données projection ACP1, puis blanchiment et enfin projection ACP2
# - enregistrer éventuellement les résultats avant visualization

# appliquer aux données projection ACP1, puis blanchiment et enfin projection ACP2
ldaTexturesDFA = ldaModel.transform(whitener.transform(pcaModel \
                         .transform(scaledTexturesDFA))) \
                         .select("label","ldaFeatures")
ldaTexturesDFA.write.parquet("tpafd/data/ldaTexturesDFA")

Nous pouvons ensuite réaliser une visualisation 2D des résultats :

from pyspark.ml.functions import vector_to_array
from pyspark.sql.functions import col
import pandas as pd

# Extraction de deux colonnes à partir de l'unique colonne de vecteurs
#  et création de dataframe pandas à partir du DataFrame Spark
donneesAafficher = ldaTexturesDFA.withColumn("x", vector_to_array("ldaFeatures")) \
                     .select(["label"] + [col("x")[i] for i in range(2)]).toPandas()
donneesAafficher.head(5)

# Visualisation de type scatter plot avec une couleur par étiquette
import matplotlib.pyplot as plt
donneesAafficher.plot.scatter('x[0]', 'x[1]', c='label', cmap='coolwarm')

Les données appartenant à une même classe sont représentées avec une même couleur. Vu le nombre élevé de classes (11), certaines nuances sont proches entre classes différentes. Vous pouvez trouver une projection des mêmes données sur les 2 premiers axes discriminants dans la Fig. 13 de cette section de RCP216.