Latent profilanalys (LPA) försöker identifiera kluster av individer (dvs. latenta profiler) baserat på svar på en serie kontinuerliga variabler (dvs. indikatorer). LPA antar att det finns obemärkta latenta profiler som genererar svarsmönster på indikatorposter.
Här kommer jag att gå igenom ett snabbt exempel på LPA för att identifiera grupper av människor baserat på deras intressen/hobbyer. Uppgifterna kommer från Ungdomsundersökningen, tillgänglig fritt på Kaggle.com.
Här är en smygtitt på vad vi går för:
terminologi Obs: människor använder termerna kluster, profiler, klasser och grupper omväxlande, men det finns subtila skillnader. Jag håller mig mest till profilen för att hänvisa till en gruppering av fall, i enlighet med LPA-terminologi. Vi bör notera att LPA är en gren av Gaussisk ändlig Blandningsmodellering, som inkluderar Latent klassanalys (LCA). Skillnaden mellan LPA och LCA är konceptuell, inte beräknings: LPA använder kontinuerliga indikatorer och LCA använder binära indikatorer. LPA är en probabilistisk modell, vilket innebär att den modellerar sannolikheten för fall som tillhör en profil. Detta är överlägset ett tillvägagångssätt som K-medel som använder avståndsalgoritmer.
med det åt sidan, låt oss ladda in data.
library(tidyverse)
survey <- read_csv("https://raw.githubusercontent.com/whipson/tidytuesday/master/young_people.csv") %>% select(History:Pets)
data är på 32 intressen/hobbies. Varje objekt rankas 1 (inte intresserad) till 5 (mycket intresserad).
beskrivningen på Kaggle föreslår att det kan vara slarvigt att svara (t.ex. deltagare som valt samma värde om och om igen). Vi kan använda det slarviga paketet för att identifiera ”strängsvar”. Låt oss också leta efter multivariata outliers med Mahalanobis avstånd (Se mitt tidigare inlägg på Mahalanobis för att identifiera outliers).
library(careless)library(psych)interests <- survey %>% mutate(string = longstring(.)) %>% mutate(md = outlier(., plot = FALSE))
vi täcker strängen som svarar på högst 10 och använder en Mahalanobis D cutoff av alpha = .001.
cutoff <- (qchisq(p = 1 - .001, df = ncol(interests)))interests_clean <- interests %>% filter(string <= 10, md < cutoff) %>% select(-string, -md)
paketet mclust utför olika typer av modellbaserad kluster och dimensionsreduktion. Dessutom är det verkligen intuitivt att använda. Det kräver fullständig data (ingen saknas), så i det här exemplet tar vi bort fall med NAs. Detta är inte det föredragna tillvägagångssättet; vi skulle vara bättre att tillskriva. Men för illustrativa ändamål fungerar det bra. Jag kommer också att standardisera alla indikatorer så när vi plottar profilerna är det tydligare att se skillnaderna mellan kluster. Att köra den här koden tar några minuter.
library(mclust)interests_clustering <- interests_clean %>% na.omit() %>% mutate_all(list(scale))BIC <- mclustBIC(interests_clustering)
vi börjar med att plotta Bayesian informationskriterier för alla modeller med profiler från 1 till 9.
plot(BIC)
det är inte omedelbart klart vilken modell som är bäst eftersom y-axeln är så stor och många av modellerna poäng nära varandra. sammanfattning (BIC) visar de tre bästa modellerna baserade på 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.461
den högsta BIC kommer från VVE, 3. Detta säger att det finns 3 kluster med variabel volym, variabel form, lika orientering och ellipsodial fördelning (se Figur 2 från detta papper för en visuell). VEE, 3 är dock inte långt efter och kan faktiskt vara en mer teoretiskt användbar modell eftersom den begränsar fördelningens form att vara lika. Av den anledningen går vi med VEE, 3.
om vi vill titta närmare på den här modellen sparar vi den som ett objekt och inspekterar det med sammanfattning().
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 210
utgången beskriver profilernas geometriska egenskaper och antalet fall som klassificeras i var och en av de tre klusterna.
BIC är ett av de bästa passformsindexen, men det rekommenderas alltid att leta efter mer bevis på att lösningen vi har valt är den rätta. Vi kan också jämföra värden för det integrerade färdiga Likelikood-kriteriet (ICL). Se detta dokument för mer information. ICL skiljer sig inte mycket från BIC, förutom att det lägger till ett straff på lösningar med större entropi eller klassificeringsosäkerhet.
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.203
Vi ser liknande resultat. ICL föreslår att modellen VEE, 3 passar ganska bra. Slutligen utför vi Bootstrap Likelihood Ratio Test (BLRT) som jämför modellpassning mellan k-1 och k klustermodeller. Med andra ord, det ser ut att se om en ökning av profiler ökar passformen. Baserat på simuleringar av Nylund, Asparouhov och Muth Usbi (2007) är BIC och BLRT de bästa indikatorerna för hur många profiler det finns. Denna kodrad tar lång tid att springa, så om du bara följer med föreslår jag att du hoppar över den om du inte vill gå ut för en kaffepaus.
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.000
BLRT föreslår också att en 3-profillösning är idealisk.
visualisera LPA
Nu när vi är säkra på vårt val av en 3-profillösning, låt oss plotta resultaten. Specifikt vill vi se hur profilerna skiljer sig åt på indikatorerna, det vill säga de objekt som utgör profilerna. Om lösningen är teoretiskt meningsfull bör vi se skillnader som är vettiga.
först extraherar vi medel för varje profil (kom ihåg att vi valde att dessa skulle standardiseras). Sedan använder vi pivot_longer för att krossa den i lång form. Observera att jag trimmar värden som överstiger +1 SD, annars stöter vi på plottningsproblem.
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))
här är koden för tomten. Jag ordnar om indikatorerna så att liknande aktiviteter ligger nära varandra.
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")
Vi har många indikatorer (mer än typiska för LPA), men vi ser några intressanta skillnader. Det är uppenbart att den röda gruppen är intresserad av vetenskap och den blå gruppen visar större intresse för konst och humaniora. Den gröna gruppen verkar ointresserad av både vetenskap och konst, men måttligt intresserad av andra saker.
Vi kan göra denna plot mer informativ genom att ansluta profilnamn och proportioner. Jag kommer också att spara denna plot som ett objekt så att vi kan göra något riktigt coolt med det!
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ågot riktigt coolt som jag vill göra är att göra en interaktiv plot. Varför skulle jag vilja göra det här? Tja, ett av problemen med den statiska tomten är att med så många indikatorer är det svårt att läsa värdena för varje indikator. En interaktiv plot låter läsaren smala in på specifika indikatorer eller profiler av intresse. Vi använder plotly för att göra vår statiska plot till en interaktiv.
library(plotly)ggplotly(p, tooltip = c("Interest", "Mean")) %>% layout(legend = list(orientation = "h", y = 1.2))
det finns ett snabbt exempel på LPA. Sammantaget tycker jag att LPA är ett bra verktyg för undersökande analys, även om jag ifrågasätter dess Reproducerbarhet. Det som är viktigt är att statistikern tar hänsyn till både passande index och teori när man beslutar om antalet profiler.
referenser & resurser
Bertoletti, M., Friel, N., & Rastelli, R. (2015). Välja antalet kluster i en ändlig blandningsmodell med hjälp av ett exakt integrerat slutfört Sannolikhetskriterium. https://arxiv.org/pdf/1411.4257.pdf.
Nylund, K. L., Asparouhov, T., & Muth Ukrainian, B. O. (2007). Beslut om antalet klasser i Latent klassanalys och Tillväxtblandningsmodellering: en Monte Carlo-Simuleringsstudie. Strukturell Ekvationsmodellering, 14, 535-569.
Scrucca, L., Fop, M., Murphy, T. B., & Raftery, A. E. (2016). mclust5: klustring, klassificering och Densitetsuppskattning med gaussiska finita Blandningsmodeller. R Journal, 8, 289-317.