COMPITO SU CLUSTERING E ACP
Erik De Luca
30/1/2022
COMPITO SU CLUSTERING E ACP
1. Analisi dei gruppi sul dataset feeling.csv
a) scegliere il numero di gruppi
Come primo passaggio importo i dati dal file .csv e osservo com’è strutturato.
feeling <- read.csv("feeling.csv", header=T, sep=";")
colnames(feeling)=stringr::str_remove(colnames(feeling),"ft_")
colnames(feeling)=stringr::str_remove(colnames(feeling),"_2016")
show(feeling[c(1:5),])## presvote16post black white hisp asian muslim jew christ fem immig blm wallst
## 1 2 95 94 90 90 89 94 50 100 89 79 21
## 2 1 97 99 100 100 75 97 100 52 97 29 75
## 3 2 77 83 77 90 51 75 50 69 75 89 25
## 4 1 91 75 90 93 65 95 91 36 90 2 30
## 5 2 100 97 98 98 99 100 26 100 100 100 NA
## gays unions police altright
## 1 99 78 89 2
## 2 84 48 97 51
## 3 50 99 94 0
## 4 70 30 96 1
## 5 100 100 52 NASi nota subito che effettuare un’analisi dei gruppi su questo dataframe risulterebbe di difficile interpretazione. Quello che ci interessa sapere è: “Chi ha votato cosa?”, quindi nelle righe del nuovo dataframe inserirò i vari gruppi mentre nelle colonne cosa hanno votato.
Per analizzare i gruppi calcolo le distane tra ogni unità statistica. Prima di fare ciò standardizzo i valori per evitare che qualche modalità risulti più rilevante di altre solo per l’ordine di grandezza. In questo caso l’effetto non sarebbe comunque significativo.
Provo a utilizzare due metodi diversi per la ricerca del numero ottimale di gruppi.
kMax <- 8
wss <- sapply(1:kMax,
function(k,data) kmeans(data, centers = k)$tot.withinss,
data=feelScale)
ggp <- ggplot(data = data.frame(x=1:kMax, y=wss), mapping = aes(x=x,y=y)) +
geom_point(colour = "red") +
geom_line(colour = "blue") +
xlab("Numero di gruppi") +
ylab("Somma dei quadrati entro gruppi (WGSS)")
print(ggp)
# Ricerca del livello ottimale di gruppi con un altro metodo
minC <- 2
maxC <- 10
res <- sapply(minC:maxC,
function(nc, data) index.G1(data, kmeans(data,centers = nc)$cluster),
data = feelScale)
ggp <- ggplot(data=data.frame(x=2:(length(res)+1), y= res), mapping = aes(x=x,y=y)) +
geom_point() +
geom_line() +
xlab("Numero di cluster") +
ylab("Statistica pseudo-F di Calinski-Harabasz")
print(ggp)
Il numero ottimale di gruppi è 3. Calcolo le matrici di dissimiliarità con diveri metodi di linkage. Visualizzo i dendogrammi dei diversi metodi.
d <- dist(feelScale)
# matrici di dissimiliarità
feel_sl <- hclust(d, method = "single")
feel_cl <- hclust(d, method = "complete")
feel_al <- hclust(d, method = "average")
feel_ward <- hclust(d, method = "ward.D2")
op <- par(mfrow = c(2, 2))
plot(feel_sl, labels = feeling$comp_short, cex = .7,
main = "feeling data (single linkage)", xlab = "feeling")
plot(feel_cl, labels = feeling$comp_short, cex = .7,
main = "feeling data (complete linkage)", xlab = "feeling")
plot(feel_al, labels = feeling$comp_short, cex = .7,
main = "feeling data (average linkage)", xlab = "feeling")
plot(feel_ward, labels = feeling$comp_short, cex = .7,
main = "feeling data (Ward)", xlab = "feeling")
rect.hclust(feel_ward, 3)
c3 <- cutree(feel_ward, 3)b) commentarli in base alla % di votanti repubblicani e democratici nei gruppi
Aggiungo al dataframe le colonne con le rispettive percentuali di votanti per i 2 partiti.
newf=cbind(newf,newf[,1]/(newf[,1]+newf[,2]),newf[,2]/(newf[,1]+newf[,2]))
colnames(newf)=c("tot rep","tot dem","rep%","dem%")Per rappresentare la suddivisione dei vari gruppi e per comprendere al meglio le loro caratteristiche si è scelto un grafico a dispersione anziché un grafico a torta o un grafico a barre che potevano risultare a primo impatto più significativi. Infatti con il grafico a torta e quello a barre sarebbe stato difficile comprendere le vicinanze tra ogni gruppo ma sarebbe stata più facile l’interpretazione su cosa votassero. Inoltre con il grafico a torta avremmo perso il valore assoluto dei voti.
Nel grafico infatti sotto la retta \(y=x\) si trovano le unità statistiche che hanno un percentuale maggiore di votanti reppublicani. Al di sopra, invece, si trovano le unità statistiche che hanno un percentuale maggiore di votanti reppublicani. Le unità statistiche più vicine all’origine saranno i gruppi rappresentati da meno persone, mentre quelli con una distanza maggiore saranno quelli con di più.
Gruppo=paste(rownames(newf),"\n",
"Percentuale repubblicani: ",round(newf$`rep%`*100,digits = 2),"%",
"\n", "Percentuale democratici: ",round(newf$`dem%`*100,digits = 2),"%")
ggp <- ggplot(data= newf,
mapping = aes(x = `tot rep`, y = `tot dem`,label=Gruppo,colour=as.factor(c3))) +
geom_point(size=1) +
xlab("Totale voti repubblicani") +
ylab("Totale voti democratici")+
geom_abline(slope = 1,col="blue",size=0.1)
ggplotly(ggp)2. Analisi delle componenti principali sul dataset autovetture
a) che % di varianza nei dati spiegano le prime 2 componenti principali?
Prima di iniziare a fare qualsiasi lavoro è importante controllare la qualità dei dati. Quindi va fatta una analisi esplorativa dei dati e fare le opportune modifiche per avere dei dati significativi e rappresentativi.
#lettura dati
dataAuto = read.csv("autovetture.csv",
sep = ";",
dec = ",",
row.names = 1)
summary(dataAuto)## marca modello versione cilindri
## Length:76 Length:76 Length:76 Min. : 3.000
## Class :character Class :character Class :character 1st Qu.: 4.000
## Mode :character Mode :character Mode :character Median : 4.000
## Mean : 4.868
## 3rd Qu.: 5.000
## Max. :12.000
##
## cilindrata compressione potenza giri_pot_max coppia
## Min. : 5 Min. : 6.00 Min. : 7.0 Min. : 8 Min. : 9.0
## 1st Qu.:1595 1st Qu.:10.50 1st Qu.: 74.0 1st Qu.:4000 1st Qu.:150.0
## Median :1991 Median :11.10 Median :101.5 Median :5450 Median :225.0
## Mean :2380 Mean :13.26 Mean :130.3 Mean :5096 Mean :263.7
## 3rd Qu.:2916 3rd Qu.:17.25 3rd Qu.:147.8 3rd Qu.:6000 3rd Qu.:320.0
## Max. :6749 Max. :22.00 Max. :640.0 Max. :8500 Max. :720.0
## NA's :25
## giri_coppia_max marce diam_sterzata diam_freni_ant
## Min. : 10 Min. : 4.000 Min. : 4.40 Min. : 13.0
## 1st Qu.:2000 1st Qu.: 5.000 1st Qu.:10.47 1st Qu.:257.0
## Median :3000 Median : 5.000 Median :10.90 Median :284.0
## Mean :3018 Mean : 5.464 Mean :10.76 Mean :280.3
## 3rd Qu.:3925 3rd Qu.: 6.000 3rd Qu.:11.43 3rd Qu.:301.0
## Max. :6000 Max. :11.000 Max. :12.55 Max. :398.0
## NA's :7 NA's :20 NA's :44
## diam_freni_post lung larg alt passo
## Min. : 14.0 Min. : 15 Min. : 16 Min. : 17 Min. : 18
## 1st Qu.:228.0 1st Qu.:4052 1st Qu.:1735 1st Qu.:1430 1st Qu.:2476
## Median :251.0 Median :4360 Median :1804 Median :1510 Median :2600
## Mean :250.2 Mean :4289 Mean :1785 Mean :1545 Mean :2588
## 3rd Qu.:290.0 3rd Qu.:4674 3rd Qu.:1895 3rd Qu.:1704 3rd Qu.:2735
## Max. :355.0 Max. :5830 Max. :2200 Max. :1960 Max. :3302
## NA's :44 NA's :8
## peso rim_frenato rim_non_frenato bagag
## Min. : 19 Min. : 20 Min. : 21.0 Min. : 22.0
## 1st Qu.:1291 1st Qu.:1050 1st Qu.:650.0 1st Qu.: 215.0
## Median :1440 Median :1300 Median :750.0 Median : 400.0
## Mean :1532 Mean :1478 Mean :673.5 Mean : 506.2
## 3rd Qu.:1849 3rd Qu.:2000 3rd Qu.:750.0 3rd Qu.: 473.0
## Max. :2845 Max. :3500 Max. :800.0 Max. :3200.0
## NA's :38 NA's :63 NA's :19
## serbatoio vel acc con_urbano
## Min. : 23.00 Min. : 19.4 Min. : 3.40 Min. : 5.26
## 1st Qu.: 53.50 1st Qu.:167.2 1st Qu.: 8.70 1st Qu.: 7.85
## Median : 62.00 Median :185.0 Median :10.50 Median : 9.90
## Mean : 64.28 Mean :192.5 Mean :10.99 Mean :12.33
## 3rd Qu.: 70.00 3rd Qu.:215.2 3rd Qu.:12.60 3rd Qu.:12.25
## Max. :195.00 Max. :340.0 Max. :25.00 Max. :82.00
## NA's :5 NA's :3 NA's :1
## con_extraurb con_misto emissioni benzina
## Min. : 3.400 Min. : 4.090 Min. : 27.3 Min. : 0.000
## 1st Qu.: 5.150 1st Qu.: 6.200 1st Qu.:149.0 1st Qu.: 0.000
## Median : 6.200 Median : 7.600 Median :184.0 Median : 1.000
## Mean : 7.171 Mean : 8.733 Mean :194.4 Mean : 1.118
## 3rd Qu.: 7.750 3rd Qu.: 9.425 3rd Qu.:221.0 3rd Qu.: 1.000
## Max. :27.000 Max. :28.000 Max. :495.0 Max. :30.000
## NA's :1 NA's :7
## diesel traz_ant traz_post traz_int
## Min. : 0.0000 Min. : 0.000 Min. : 0.0000 Min. : 0.0000
## 1st Qu.: 0.0000 1st Qu.: 0.000 1st Qu.: 0.0000 1st Qu.: 0.0000
## Median : 0.0000 Median : 1.000 Median : 0.0000 Median : 0.0000
## Mean : 0.6711 Mean : 1.092 Mean : 0.5132 Mean : 0.6842
## 3rd Qu.: 1.0000 3rd Qu.: 1.000 3rd Qu.: 0.0000 3rd Qu.: 0.2500
## Max. :31.0000 Max. :32.000 Max. :33.0000 Max. :34.0000
##
## alim trazione
## Length:76 Length:76
## Class :character Class :character
## Mode :character Mode :character
##
##
##
## dataAuto=dataAuto[-1,]
# elimino le colonne con molti na
colnames(dataAuto)## [1] "marca" "modello" "versione" "cilindri"
## [5] "cilindrata" "compressione" "potenza" "giri_pot_max"
## [9] "coppia" "giri_coppia_max" "marce" "diam_sterzata"
## [13] "diam_freni_ant" "diam_freni_post" "lung" "larg"
## [17] "alt" "passo" "peso" "rim_frenato"
## [21] "rim_non_frenato" "bagag" "serbatoio" "vel"
## [25] "acc" "con_urbano" "con_extraurb" "con_misto"
## [29] "emissioni" "benzina" "diesel" "traz_ant"
## [33] "traz_post" "traz_int" "alim" "trazione"colRimosse = c(6,12:14,20:22)
dataAuto = dataAuto[,-colRimosse]
# Seleziono solo colonne numeriche
colNumeriche = unlist(lapply(dataAuto, is.numeric))
dataAuto = na.omit(dataAuto[, colNumeriche])
#rimuovo altre variabili che non mi interessa analizzare
colRimosse = c(4,6,7,11,13,22:24)
dataAuto = dataAuto[,-colRimosse]autoPCA <- princomp(dataAuto, cor = TRUE)
#summary(autoPCA, loadings = TRUE)
par(mfrow=c(1,1))
plot(autoPCA, type = "l")
abline(h = 1, lty = 2)
Dovrò utilizzare almeno 4 componenti secondo il criterio della “proporzione della varianza spiegata”. Nel grafico si vede quanta percentuale di varianza nei dati spiegano per ogni componente principale. oppure nella tabella qua sotto:
autoPCA$sdev[1:7]## Comp.1 Comp.2 Comp.3 Comp.4 Comp.5 Comp.6 Comp.7
## 2.9195338 1.8446704 1.0357538 0.9525344 0.7991014 0.6128047 0.4887481b) Commentare la collocazione delle variabili rispetto alle componenti principali
La prima componente principale sembra contrapporre misure di “meccanica” (traz_post e ant, marce, acc, …) con misure della dimensione (larghezza, lunghezza, …), mentre la seconda componente principale sembra contrapporre la potenza e le dimensioni del motore con le dimensioni dell’auto.La prima componente principale sembra contrapporre misure di “meccanica” (traz_post e ant, marce, acc, …) con misure della dimensione (larghezza, lunghezza, …), mentre la seconda componente principale sembra contrapporre la potenza e le dimensioni del motore con le dimensioni dell’auto.
ggcorr(cbind(dataAuto, autoPCA$scores[,1:5]), label = TRUE, cex = 3,label_size = 2)
- Commentare la collocazione degli individui nello spazio delle variabili
A causa della grandezza di osservazioni conviene usare un grafico interattivo.
Nel gruppo 3 si trovano le auto potenti che infatti sono più lontane sull’asse delle x, si nota che l’Hummer H1 ha un basso valore della seconda componente per le sue dimensioni che infatti la posizionano alla stessa altezza di altri suv che si trovano nel gruppo 2.
Nel gruppo 1 si trovano le auto leggere, di piccole dimensioni e poco potenti: le city car. Mentre nel gruppo 2 si trovano auto più potenti.
km4 <- kmeans(scale(dataAuto), centers = 4)
Cluster <- as.character(km4$cluster)
autoScores <- cbind(data.frame(autoPCA$scores[, c("Comp.1", "Comp.2")]),
marca=rownames(autoPCA$scores),
Cluster=Cluster)
ggp <- ggplot(data= autoScores,
mapping = aes(x = Comp.1, y=Comp.2, label=marca,colour=Cluster)) +
geom_point(size=1) +
xlab("Prima dimensione PCA") +
ylab("Seconda dimensione PCA")
#geom_text(hjust=0.5, vjust=-0.5, size=2)
ggplotly(ggp) 3. Realizzare un’analisi dei gruppi commentando i diversi raggruppamenti ottenuti mediante il metodo del legame singolo, completo, medio, centroidi e Ward
dataCrimini <- read.table("UScrime.txt", header = T)
# rimuoviamo eventuali valori mancanti e standardizziamo i valori
dataCrimini = scale(na.omit(dataCrimini))Prima di analizzare i dati: ometto eventuali dati mancanti (NA) e standardizzo i valori, in modo che alcune variabili non influiscano maggiormente su altre.
# uso la clustering gerarchica:
clust_single <- hclust(dist(dataCrimini),method="single") #metodo del legame singolo
plot(clust_single)
rect.hclust(clust_single ,5, border = 3)
Con il metodo del legame singolo risulta difficile stabilire il numero di gruppi, poiché se ne forma uno molto popoloso e gli altri composti da singoli stati. Il numero scelto di gruppi è stato 5.
clust_compl <- hclust(dist(dataCrimini),method="complete") #metodo del legame completo
plot(clust_compl)
rect.hclust(clust_compl ,3, border = 3)
Utilizzando il metodo del legame completo la situazione migliora e il numero ottimale di gruppi è 3, così da averli sufficientemente popolati. Se aggiungessimo un ulteriore gruppo, codesto avrebbe solo uno stato all’interno.
clust_medio <- hclust(dist(dataCrimini),method="average") #metodo del legame medio
plot(clust_medio)
rect.hclust(clust_medio ,5, border = 3)
Con il metodo del legame medio il raggruppamento migliore si ottiene con 5.
clust_centr <- hclust(dist(dataCrimini),method="centroid") #metodo del centroide
plot(clust_centr)
rect.hclust(clust_centr ,4, border = 3)
Il metodo dei centroidi non è adatto all’analisi di questo dataset, come nel caso del metodo singolo, il raggruppamento migliore lo si ha con 4 ma non è comunque soddisfacente e utile nella pratica.
clust_ward <- hclust(dist(dataCrimini),method="ward.D") #metodo di Ward
plot(clust_ward)
rect.hclust(clust_ward ,6, border = 3)
Nel metodo di Ward si potrebbero utilizzare 6 gruppi e avremmo tutti e 6 gruppi sufficientemente popolati.
Il metodo del legame medio e il metodo di Ward suddividono in maniera più proporzionale il grafico, confrantandoli con tutti gli altri metodi sono i più sensati.