Esempi di cluster

quasiworkshop

________________

Simulazione di 2 cluster veri

Simulo 2 cluster mediamente separati su due indicatori X1 e X2, anche la correlazione tra X1 e X2 differisce in ciascun cluster, così come la varianza

set.seed(20)

rm(list=ls())
# Load packages
library(ggplot2)
library(GGally)
library(mclust)
library(MASS)
library(semTools)
library(lavaan)

N = 1000
g1 = data.frame(mvrnorm(n=N,mu=c(0,0),Sigma=matrix(c(1,.6,
                                                     .6,1),nrow=2)))
g2 = data.frame(mvrnorm(n=N,mu=c(5,5),Sigma=matrix(c(2,-.3,
                                                     -.3,2),nrow=2)))
d = rbind(g1,g2)
d$ClusterVero = rep(c("G1","G2"),each=N)

ggplot(d,aes(x=X1,y=X2,group=ClusterVero,color=ClusterVero,fill=ClusterVero))+
  geom_smooth(method="lm",size=2,alpha=.2)+
  geom_point(size=2.5,alpha=.3)+
  theme(text=element_text(size=25))

Ora usiamo il model-based clustering (Gaussian mixture model) su questi dati, valutiamo soluzioni da 1 fino a un massimo di 5 cluster (componenti)

fit = Mclust(d[,c("X1","X2")],G=1:5)

Quanti cluster/componenti ha rilevato il modello? “2 components”

summary(fit)

## ---------------------------------------------------- 
## Gaussian finite mixture model fitted by EM algorithm 
## ---------------------------------------------------- 
## 
## Mclust VVE (ellipsoidal, equal orientation) model with 2 components: 
## 
##  log-likelihood    n df       BIC      ICL
##       -7475.668 2000 10 -15027.34 -15041.5
## 
## Clustering table:
##    1    2 
## 1002  998

Vediamone i parametri… sembra averli colti molto bene!

fit$parameters

## $pro
## [1] 0.5004552 0.4995448
## 
## $mean
##            [,1]     [,2]
## X1  0.005055991 4.980269
## X2 -0.025455527 5.028909
## 
## $variance
## $variance$modelName
## [1] "VVE"
## 
## $variance$d
## [1] 2
## 
## $variance$G
## [1] 2
## 
## $variance$sigma
## , , 1
## 
##           X1        X2
## X1 1.0010470 0.6496645
## X2 0.6496645 1.0566082
## 
## , , 2
## 
##            X1         X2
## X1  1.9614423 -0.2884035
## X2 -0.2884035  1.9367773
## 
## 
## $variance$scale
## [1] 0.7972769 1.9276152
## 
## $variance$shape
##          [,1]      [,2]
## [1,] 0.474828 1.1609043
## [2,] 2.106026 0.8613974
## 
## $variance$orientation
##            X1         X2
## X1  0.7220534 -0.6918373
## X2 -0.6918373 -0.7220534
## 
## 
## $Vinv
## NULL

________________

________________

Simulazione di 1 solo cluster vero ma in uno scenario “critico”

Simulo i dati, N = 2000, due indicatori, centrati su zero, lieve correlazione positiva, skewness a -0.20 e +0.30 rispettivamente nei due indicatori, cioè molto poco

set.seed(20)

d = data.frame(mvrnonnorm(n=2000,mu=c(0,0),Sigma=matrix(c(1,.2,.2,1),nrow=2),skewness=c(-0.2,0.3)))

Come prima cosa vediamo la distribuzione dei dati

ggpairs(d)+
  theme(text=element_text(size=20))

Ora fittiamo un mixture model e… trova 2 cluster!

fit = Mclust(d)
summary(fit)

## ---------------------------------------------------- 
## Gaussian finite mixture model fitted by EM algorithm 
## ---------------------------------------------------- 
## 
## Mclust EVI (diagonal, equal volume, varying shape) model with 2 components: 
## 
##  log-likelihood    n df       BIC       ICL
##       -5588.801 2000  8 -11238.41 -12425.84
## 
## Clustering table:
##    1    2 
##  918 1082

Ma cosa avrà mai trovato?! Plottiamo

d$class = as.factor(fit$classification)

centers = data.frame(t(fit$parameters$mean))

ggplot(d,aes(x=X1,y=X2))+
  geom_point(aes(group=class,color=class,fill=class),size=2.5,alpha=.5)+
  geom_point(data=centers,color="black",size=6,stroke=2,shape=4)+
  theme(text=element_text(size=25))

Sembra incredibile, ma quel poco di asimmetria ha ingannato il modello. Risimulando senza alcuna asimmetria, infatti, trova correttamente un solo cluster

set.seed(20)

dAlt = data.frame(mvrnonnorm(n=2000,mu=c(0,0),Sigma=matrix(c(1,.2,.2,1),nrow=2),skewness=c(0,0)))
fitAlt = Mclust(dAlt)
summary(fitAlt)

## ---------------------------------------------------- 
## Gaussian finite mixture model fitted by EM algorithm 
## ---------------------------------------------------- 
## 
## Mclust XXX (ellipsoidal multivariate normal) model with 1 component: 
## 
##  log-likelihood    n df      BIC      ICL
##       -5599.596 2000  5 -11237.2 -11237.2
## 
## Clustering table:
##    1 
## 2000