Latent Profile Analysis (LPA) se snaží identifikovat shluky jednotlivců (tj. LPA předpokládá, že existují nepozorované latentní profily, které generují vzorce odpovědí na položky indikátoru.
zde projdu rychlý příklad LPA k identifikaci skupin lidí na základě jejich zájmů / koníčků. Data pocházejí z průzkumu mladých lidí, volně k dispozici na Kaggle.com.
zde je náhled na to, co hledáme:
terminologie poznámka: lidé používají pojmy klastry, profily, třídy a skupiny zaměnitelně, ale existují jemné rozdíly. Většinou se budu držet profilu, abych odkazoval na seskupení případů, v souladu s terminologií LPA. Měli bychom poznamenat, že LPA je větev Gaussova modelování konečných směsí, která zahrnuje analýzu latentní třídy (LCA). Rozdíl mezi LPA a LCA je koncepční, nikoli výpočetní: LPA používá spojité indikátory a LCA používá binární indikátory. LPA je pravděpodobnostní model, což znamená, že modeluje pravděpodobnost, že případ patří do profilu. To je lepší než přístup jako K-znamená, že používá algoritmy vzdálenosti.
s tím stranou, pojďme načíst data.
library(tidyverse)survey <- read_csv("https://raw.githubusercontent.com/whipson/tidytuesday/master/young_people.csv") %>% select(History:Pets)údaje jsou o 32 zájmech/koníčcích. Každá položka je hodnocena 1 (nezajímá) na 5 (velký zájem).
popis na Kaggle naznačuje, že může být neopatrný reagovat (např účastníci, kteří si vybrali stejnou hodnotu znovu a znovu). Neopatrný balíček můžeme použít k identifikaci „odpovědi na řetězec“. Podívejme se také na vícerozměrné odlehlé hodnoty s vzdáleností Mahalanobis (viz můj předchozí příspěvek na Mahalanobis pro identifikaci odlehlých hodnot).
library(careless)library(psych)interests <- survey %>% mutate(string = longstring(.)) %>% mutate(md = outlier(., plot = FALSE))Budeme cap řetězec reagovat na maximálně 10 a použít Mahalanobisova D cutoff alfa = .001.
cutoff <- (qchisq(p = 1 - .001, df = ncol(interests)))interests_clean <- interests %>% filter(string <= 10, md < cutoff) %>% select(-string, -md)balíček mclust provádí různé typy shlukování založených na modelu a redukce rozměrů. Navíc je to opravdu intuitivní použití. Vyžaduje kompletní data (chybí), takže v tomto příkladu odstraníme případy s NAs. To není preferovaný přístup; bylo by lepší, kdybychom to přičítali. Ale pro ilustrativní účely to funguje dobře. Budu také standardizovat všechny ukazatele, takže když vykreslíme profily, je jasnější vidět rozdíly mezi klastry. Spuštění tohoto kódu bude trvat několik minut.
library(mclust)interests_clustering <- interests_clean %>% na.omit() %>% mutate_all(list(scale))BIC <- mclustBIC(interests_clustering)začneme tím, že vykreslování Bayesovského Informačního Kritéria pro všechny modely s profily v rozmezí od 1 do 9.
plot(BIC)
To není hned jasné, který model je nejlepší od osy y je tak velký, a mnoho modelů skóre blízko u sebe. shrnutí (BIC) ukazuje tři nejlepší modely založené na 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.461nejvyšší BIC pochází z VVE, 3. To říká, že tam jsou 3 shluky s proměnlivým objemem, variabilní tvar, stejné orientace, a ellipsodial distribuce (viz Obrázek 2 z tohoto dokumentu pro vizuální). VEE, 3 však není pozadu a ve skutečnosti může být teoreticky užitečnějším modelem, protože omezuje tvar rozdělení tak, aby byl stejný. Z tohoto důvodu půjdeme s VEE, 3.
Pokud se chceme na tento model podívat blíže, uložíme jej jako objekt a zkontrolujeme jej pomocí shrnutí ().
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 210výstup popisuje geometrické vlastnosti profilů a počet případů klasifikovaných do každé ze tří clusterů.
BIC je jedním z nejlepších indexů, ale vždy se doporučuje hledat další důkazy, že řešení, které jsme zvolili, je správné. Můžeme také porovnat hodnoty kritéria integrovaného dokončeného Likelikood (ICL). Viz tento dokument pro více informací. ICL se příliš neliší od BIC, kromě toho, že přidává trest na řešení s větší entropií nebo klasifikační nejistotou.
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.203Můžeme vidět podobné výsledky. ICL naznačuje, že model VEE, 3 se hodí docela dobře. Nakonec provedeme test pravděpodobnosti Bootstrap (BLRT), který porovnává přizpůsobení modelu mezi modely klastrů k-1 a k. Jinými slovy, vypadá to, zda zvýšení profilů zvyšuje fit. Na základě simulací Nylund, Asparouhov a Muthén (2007) jsou BIC a BLRT nejlepšími ukazateli toho, kolik profilů existuje. Spuštění tohoto řádku kódu bude trvat dlouho, takže pokud jen sledujete, navrhuji jej přeskočit, pokud nechcete vystoupit na přestávku na kávu.
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 také naznačuje, že 3-profil řešení je ideální.
vizualizace LPA
Nyní, když jsme přesvědčeni o výběru 3-profilového řešení, vykreslíme výsledky. Konkrétně chceme vidět, jak se profily liší na indikátorech,tj. Pokud je řešení teoreticky smysluplné, měli bychom vidět rozdíly, které dávají smysl.
nejprve extrahujeme prostředky pro každý profil(nezapomeňte, že jsme je vybrali jako standardizované). Poté použijeme pivot_longer, abychom ho hádali do dlouhé formy. Všimněte si, že ořezávám hodnoty přesahující + 1 SD, jinak narazíme na problémy s vykreslováním.
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))zde je Kód pro spiknutí. Přeskupuji ukazatele tak, aby podobné aktivity byly blízko sebe.
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")
Máme spoustu ukazatele (více než typické pro LPA), ale můžeme vidět některé zajímavé rozdíly. Je zřejmé, že červená skupina se zajímá o vědu a modrá Skupina projevuje větší zájem o umění a humanitní vědy. Zdá se, že zelená skupina nemá zájem o vědu i umění, ale mírně se zajímá o jiné věci.
Tento graf můžeme učinit informativnějším zapojením názvů profilů a proporcí. Chystám se také uložit tento spiknutí jako objekt, abychom s ním mohli udělat něco opravdu 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
něco opravdu cool, že chci udělat je interaktivní pozemku. Proč bych to měl dělat? Jedním z problémů se statickým grafem je, že s tolika indikátory je těžké číst hodnoty pro každý indikátor. Interaktivní graf umožňuje čtenáři zúžit konkrétní ukazatele nebo profily zájmu. Použijeme plotly, abychom změnili náš statický pozemek na interaktivní.
library(plotly)ggplotly(p, tooltip = c("Interest", "Mean")) %>% layout(legend = list(orientation = "h", y = 1.2))existuje rychlý příklad LPA. Celkově si myslím, že LPA je skvělý nástroj pro průzkumnou analýzu, i když zpochybňuji její Reprodukovatelnost. Důležité je, že statistik při rozhodování o počtu profilů bere v úvahu jak fit indexy, tak teorii.
Odkazy & Zdroje
Bertoletti, M., Friel, N., & Rastelli, R. (2015). Výběr počtu shluků v modelu konečné směsi pomocí přesného integrovaného kritéria vyplněné pravděpodobnosti. https://arxiv.org/pdf/1411.4257.pdf.
Nylund, K. L., Asparouhov, T., & Muthén, B. O (2007). Rozhodování o počtu tříd v analýze latentní třídy a modelování růstových směsí: simulační studie Monte Carlo. Modelování Strukturní Rovnice, 14, 535-569.
Scrucca, L., Fop, m., Murphy, T. b., & Raftery, A. E. (2016). mclust5: shlukování, klasifikace a odhad hustoty pomocí Gaussovských modelů konečných směsí. Deník R, 8, 289-317.