L’Analyse de profil latent (LPA) tente d’identifier des groupes d’individus (c’est-à-dire des profils latents) en fonction des réponses à une série de variables continues (c’est-à-dire des indicateurs). LPA suppose qu’il existe des profils latents non observés qui génèrent des modèles de réponses sur les éléments de l’indicateur.
Ici, je vais passer par un exemple rapide de LPA pour identifier des groupes de personnes en fonction de leurs intérêts / passe-temps. Les données proviennent de l’Enquête sur les jeunes, disponible gratuitement sur Kaggle.com .
Voici un aperçu de ce que nous recherchons:
Note terminologique: Les gens utilisent les termes clusters, profils, classes et groupes de manière interchangeable, mais il existe des différences subtiles. Je vais surtout m’en tenir au profil pour faire référence à un regroupement de cas, conformément à la terminologie LPA. Il convient de noter que LPA est une branche de la Modélisation des Mélanges Finis Gaussiens, qui inclut l’Analyse de Classes Latentes (ACV). La différence entre LPA et LCA est conceptuelle et non informatique: LPA utilise des indicateurs continus et LCA utilise des indicateurs binaires. LPA est un modèle probabiliste, ce qui signifie qu’il modélise la probabilité que le cas appartienne à un profil. C’est supérieur à une approche comme K-means qui utilise des algorithmes de distance.
Avec cela de côté, chargeons les données.
library(tidyverse)survey <- read_csv("https://raw.githubusercontent.com/whipson/tidytuesday/master/young_people.csv") %>% select(History:Pets)Les données portent sur 32 intérêts / passe-temps. Chaque article est classé de 1 (pas intéressé) à 5 (très intéressé).
La description sur Kaggle suggère qu’il peut y avoir une réponse négligente (par exemple, les participants qui ont sélectionné la même valeur encore et encore). Nous pouvons utiliser le package careless pour identifier « chaîne répondant ». Cherchons également des valeurs aberrantes multivariées avec la distance de Mahalanobis (voir mon post précédent sur Mahalanobis pour identifier les valeurs aberrantes).
library(careless)library(psych)interests <- survey %>% mutate(string = longstring(.)) %>% mutate(md = outlier(., plot = FALSE))Nous allons plafonner la chaîne répondant à un maximum de 10 et utiliser une coupure Mahalanobis D de alpha =.001.
cutoff <- (qchisq(p = 1 - .001, df = ncol(interests)))interests_clean <- interests %>% filter(string <= 10, md < cutoff) %>% select(-string, -md)Le paquet mclust effectue différents types de clustering basé sur le modèle et de réduction de dimension. De plus, il est vraiment intuitif à utiliser. Il nécessite des données complètes (pas manquantes), donc pour cet exemple, nous supprimerons les cas avec NAs. Ce n’est pas l’approche privilégiée; nous ferions mieux d’imputer. Mais à des fins d’illustration, cela fonctionne bien. Je vais également standardiser tous les indicateurs afin que lorsque nous traçons les profils, il soit plus clair de voir les différences entre les clusters. L’exécution de ce code prendra quelques minutes.
library(mclust)interests_clustering <- interests_clean %>% na.omit() %>% mutate_all(list(scale))BIC <- mclustBIC(interests_clustering)Nous allons commencer par tracer des critères d’information bayésiens pour tous les modèles avec des profils allant de 1 à 9.
plot(BIC)
On ne sait pas immédiatement quel modèle est le meilleur puisque l’axe des ordonnées est si grand et que de nombreux modèles se rapprochent les uns des autres. résumé (BIC) montre les trois premiers modèles basés sur BIC.
summary(BIC)## Best BIC values:## VVE,3 VEE,3 EVE,3## BIC -75042.7 -75165.1484 -75179.165## BIC diff 0.0 -122.4442 -136.461Le BIC le plus élevé provient de VVE, 3. Cela dit, il y a 3 grappes avec un volume variable, une forme variable, une orientation égale et une distribution ellipsodiale (voir la figure 2 de cet article pour un visuel). Cependant, VEE, 3 n’est pas loin derrière et peut en fait être un modèle plus théoriquement utile car il contraint la forme de la distribution à être égale. Pour cette raison, nous irons avec VEE, 3.
Si nous voulons examiner ce modèle de plus près, nous l’enregistrons en tant qu’objet et l’inspectons avec summary().
mod1 <- Mclust(interests_clustering, modelNames = "VEE", G = 3, x = BIC)summary(mod1)## ---------------------------------------------------- ## Gaussian finite mixture model fitted by EM algorithm ## ---------------------------------------------------- ## ## Mclust VEE (ellipsoidal, equal shape and orientation) model with 3 components: ## ## log-likelihood n df BIC ICL## -35455.83 874 628 -75165.15 -75216.14## ## Clustering table:## 1 2 3 ## 137 527 210La sortie décrit les caractéristiques géométriques des profils et le nombre de cas classés dans chacun des trois groupes.
BIC est l’un des meilleurs indices d’ajustement, mais il est toujours recommandé de rechercher plus de preuves que la solution que nous avons choisie est la bonne. Nous pouvons également comparer les valeurs du critère intégré Completed Likelikood (ICL). Voir cet article pour plus de détails. ICL n’est pas très différent de BIC, sauf qu’il ajoute une pénalité sur les solutions avec une plus grande incertitude d’entropie ou de classification.
ICL <- mclustICL(interests_clustering)plot(ICL)
summary(ICL)## Best ICL values:## VVE,3 VEE,3 EVE,3## ICL -75134.69 -75216.13551 -75272.891## ICL diff 0.00 -81.44795 -138.203Nous voyons des résultats similaires. ICL suggère que le modèle VEE, 3 convient assez bien. Enfin, nous allons effectuer le test de rapport de vraisemblance Bootstrap (BLRT) qui compare l’ajustement du modèle entre les modèles de cluster k-1 et k. En d’autres termes, il cherche à voir si une augmentation des profils augmente l’ajustement. Sur la base de simulations de Nylund, Asparouhov et Muthén (2007), BIC et BLRT sont les meilleurs indicateurs du nombre de profils. Cette ligne de code prendra beaucoup de temps à s’exécuter, donc si vous ne faites que suivre, je vous suggère de la sauter, sauf si vous voulez sortir pour une pause-café.
mclustBootstrapLRT(interests_clustering, modelName = "VEE")## Warning in mclustBootstrapLRT(interests_clustering, modelName = "VEE"): some## model(s) could not be fitted!## ------------------------------------------------------------- ## Bootstrap sequential LRT for the number of mixture components ## ------------------------------------------------------------- ## Model = VEE ## Replications = 999 ## LRTS bootstrap p-value## 1 vs 2 197.0384 0.001## 2 vs 3 684.8743 0.001## 3 vs 4 -124.1935 1.000BLRT suggère également qu’une solution à 3 profils est idéale.
Visualisation de LPA
Maintenant que nous sommes confiants dans notre choix d’une solution à 3 profils, tracons les résultats. Plus précisément, nous voulons voir en quoi les profils diffèrent sur les indicateurs, c’est-à-dire les éléments qui composent les profils. Si la solution est théoriquement significative, nous devrions voir des différences qui ont du sens.
Tout d’abord, nous allons extraire les moyennes pour chaque profil (rappelez-vous, nous les avons choisies pour être standardisées). Ensuite, nous utilisons pivot_longer pour le transformer en forme longue. Notez que je coupe les valeurs supérieures à + 1 SD, sinon nous rencontrons des problèmes de traçage.
means <- data.frame(mod1$parameters$mean) %>% rownames_to_column() %>% rename(Interest = rowname) %>% pivot_longer(cols = c(X1, X2, X3), names_to = "Profile", values_to = "Mean") %>% mutate(Mean = round(Mean, 2), Mean = ifelse(Mean > 1, 1, Mean))Voici le code du tracé. Je réordonne les indicateurs afin que les activités similaires soient rapprochées.
means %>% ggplot(aes(Interest, Mean, group = Profile, color = Profile)) + geom_point(size = 2.25) + geom_line(size = 1.25) + scale_x_discrete(limits = c("Active sport", "Adrenaline sports", "Passive sport", "Countryside, outdoors", "Gardening", "Cars", "Art exhibitions", "Dancing", "Musical instruments", "Theatre", "Writing", "Reading", "Geography", "History", "Law", "Politics", "Psychology", "Religion", "Foreign languages", "Biology", "Chemistry", "Mathematics", "Medicine", "Physics", "Science and technology", "Internet", "PC", "Celebrities", "Economy Management", "Fun with friends", "Shopping", "Pets")) + labs(x = NULL, y = "Standardized mean interest") + theme_bw(base_size = 14) + theme(axis.text.x = element_text(angle = 45, hjust = 1), legend.position = "top")
Nous avons beaucoup d’indicateurs (plus que typiques pour LPA), mais nous voyons des différences intéressantes. De toute évidence, le groupe rouge s’intéresse aux sciences et le groupe bleu montre un plus grand intérêt pour les arts et les sciences humaines. Le groupe vert semble désintéressé à la fois de la science et de l’art, mais modérément intéressé par d’autres choses.
Nous pouvons rendre ce tracé plus informatif en branchant les noms et les proportions des profils. Je vais aussi enregistrer cette intrigue en tant qu’objet afin que nous puissions en faire quelque chose de vraiment cool!
p <- means %>% mutate(Profile = recode(Profile, X1 = "Science: 16%", X2 = "Disinterest: 60%", X3 = "Arts & Humanities: 24%")) %>% ggplot(aes(Interest, Mean, group = Profile, color = Profile)) + geom_point(size = 2.25) + geom_line(size = 1.25) + scale_x_discrete(limits = c("Active sport", "Adrenaline sports", "Passive sport", "Countryside, outdoors", "Gardening", "Cars", "Art exhibitions", "Dancing", "Musical instruments", "Theatre", "Writing", "Reading", "Geography", "History", "Law", "Politics", "Psychology", "Religion", "Foreign languages", "Biology", "Chemistry", "Mathematics", "Medicine", "Physics", "Science and technology", "Internet", "PC", "Celebrities", "Economy Management", "Fun with friends", "Shopping", "Pets")) + labs(x = NULL, y = "Standardized mean interest") + theme_bw(base_size = 14) + theme(axis.text.x = element_text(angle = 45, hjust = 1), legend.position = "top")p
Le quelque chose de vraiment cool que je veux faire est de créer un tracé interactif. Pourquoi voudrais-je faire ça? Eh bien, l’un des problèmes avec le tracé statique est qu’avec autant d’indicateurs, il est difficile de lire les valeurs de chaque indicateur. Une intrigue interactive permet au lecteur de se concentrer sur des indicateurs ou des profils d’intérêt spécifiques. Nous utiliserons plotly pour transformer notre intrigue statique en une intrigue interactive.
library(plotly)ggplotly(p, tooltip = c("Interest", "Mean")) %>% layout(legend = list(orientation = "h", y = 1.2))Il y a un exemple rapide de LPA. Dans l’ensemble, je pense que l’APL est un excellent outil d’analyse exploratoire, bien que je remette en question sa reproductibilité. Ce qui est important, c’est que le statisticien considère à la fois les indices d’ajustement et la théorie lorsqu’il décide du nombre de profils.
Références &Ressources
Bertoletti, M., Friel, N., &Rastelli, R. (2015). Choisir le nombre de grappes dans un modèle de mélange fini en utilisant un critère de vraisemblance complet intégré exact. https://arxiv.org/pdf/1411.4257.pdf.
Nylund, K. L., Asparouhov, T., &Muthén, B. O. (2007). Décider du Nombre de Classes dans l’Analyse des Classes Latentes et la Modélisation des Mélanges de Croissance: Une Étude de Simulation de Monte Carlo. Modélisation d’équations structurelles, 14, 535-569.
Scrucca, L., Fop, M., Murphy, T. B., &Raftery, A. E. (2016). mclust5: Clustering, Classification et Estimation de la Densité À l’Aide de Modèles de Mélanges Finis Gaussiens. Le Journal R, 8, 289-317.