### Cluster Analysis: K-means & PAM

### 30.04.25
## K-Means
library(cluster)
dat<-read.csv("food_spending.csv", header = TRUE)
head(dat)
summary(dat)
x.scaled <- scale(dat) 
?kmeans
## K=2 (fissato)
km.out <- kmeans(x.scaled, centers = 2)
km.out$cluster
km.out$centers
km.out$withinss
km.out$tot.withinss
plot(dat, col=km.out$cluster, pch=19)
## K=3
km.out <- kmeans(x.scaled, centers = 3)
km.out$cluster
km.out$centers
km.out$withinss
km.out$tot.withinss
plot(dat, col=km.out$cluster, pch=19)

km1<-kmeans(dat, 3)
km1$tot.withinss
km2<-kmeans(dat, 3, nstart = 10)
km2$tot.withinss

## 5/05/25

set.seed(4)
# genero un cc ampiezza 100 da N(0,1)
# calcolo la media
hist(rnorm(1000, 0, 1))  

#Esercizio

sanita<-read.table("sanita.txt",header =TRUE,row.names=1)
str(sanita) # struttura dei dati
head(sanita)  # prime righe
View(sanita)  # dati aperti nella nuova finestra 
summary(sanita)  # indici sintetici

#standardizzazione e calcolo delle distanze euclidee
sanita_st<-scale(sanita)
d<-dist(sanita_st)
#dist(sanita_st, method = "euclidean")
# in alternativa
#library(cluster)
#daisy(sanita_st)

# eseguire kmeans con centers=3, nstart=1
set.seed(4)
K<-3
res<-kmeans(sanita_st, K)
res
?kmeans
# eseguire kmeans con centers=3, nstart=20
set.seed(4)
K<-3
res<-kmeans(sanita_st, K, nstart = 20)

rownames(sanita_st)
Ind1<-which(res$cluster==1)  # cluster1 
Ind2<-which(res$cluster==2)  # cluster 2
Ind3<-which(res$cluster==3)  # cluster 3

rownames(sanita_st)[Ind1]  # vettore dei nomi delle regioni in Cl1
rownames(sanita_st)[Ind2]  # vettore dei nomi delle regioni in Cl2
rownames(sanita_st)[Ind3]  # vettore dei nomi delle regioni in Cl3

### 12/05/2025

library(cluster)
dat<-read.csv("food_spending.csv", header = TRUE)
?pam
pam.out<-pam(dat, 2, metric="euclidean", stand=TRUE)
pam.out
pam.out$medoids 
pam.out$id.med
pam.out$clustering
plot(dat, col=(pam.out$clustering+1))
points(pam.out$medoids, pch=8, cex=2)
# selezione del numero di gruppi in base alla ASW
min_nc<-2
max_nc<-10
ASW<-array(0, c(max_nc-min_nc+1, 2))
ASW[,1]<-c(min_nc : max_nc)
for(nc in min_nc:max_nc){
  pam.out<-pam(dat, nc, metric="euclidean", 
               stand=TRUE)
  ASW[nc-min_nc+1, 2]<-pam.out$silinfo$avg.width
}
plot(ASW[,1], ASW[, 2], type="l",
     xlab="K", 
     ylab="ASW")

# fpc per diversi indici di validazione

## PAM & Gower distance
d<-read.csv(file="claims.csv", header = TRUE)
str(d)
d$litig<-as.factor(d$litig)
d$soft_injury<-as.factor(d$soft_injury)
d$emergency_tr<-as.factor(d$emergency_tr)
str(d)
# Gower's distance
mixd<-daisy(d, metric = "gower")
dclaims<-as.matrix(mixd)  # matrice simmetrica 1000x1000

#PAM  (k=3)
pam.claims<-pam(dclaims, k=3, diss=TRUE)
# dclaims non è la matrice dei dati ma la matrice 
# delle distanze di Gower
pam.claims$clustering
# aggiungere al dataset il vettore del clustering 
#dim(d)

dd<-cbind(d, clusterID=as.factor(pam.claims$clustering))
str(dd)
dim(dd)
View(dd)

# esploriamo il comportamento di "litig" in base al gruppo
table(dd$clusterID)
table(dd$litig, dd$clusterID)
results<-data.frame(table(dd$litig, dd$clusterID))
results
library(ggplot2)
ggplot(data=results,aes(x=Var2, y=Freq, fill=Var1))+
  geom_bar(stat="identity")+
  labs(x="\n Cluster", y="frequency\n", 
       title="litig", 
       fill="litig")

# split 
groups<-split(dd, dd$clusterID) # elementi per gruppo 
# il risultato è una lista 
# estrazione cluster
cl1<-groups[[1]] # cluster 1
#View(cl1)
cl2<-groups[[2]] # cluster 2
cl3<-groups[[3]] # cluster 3