## Metodi gerarchici agglomerativi (AHC)

library(cluster)  
library(fpc)
library(pdfCluster)

data(eurodist)
?eurodist
euro.mat<-as.matrix(eurodist)
#View(euro.mat)
dim(euro.mat)

#?hclust
#?agnes

# single linkage
agnes.single<-agnes(eurodist, method = "single")
# coefficiente agglomerativo
agnes.single$ac
# otteniamo il dendrogramma
pltree(agnes.single, main="Dendrogram (SL)", cex=1.2, 
       hang=-1)
abline(h=700, col=2)

# in alternativa: hclust()
hc<-hclust(eurodist, method = "single")
plot(hc, hang=-1)

# complete linkage
agnes.complete<-agnes(eurodist, method = "complete")
agnes.complete$ac
pltree(agnes.complete, main="Dendrogram (CL)",
cex=1.2, hang=-1)

# average linkage
agnes.average<-agnes(eurodist, method = "average")
agnes.average$ac
pltree(agnes.average, main="Dendrogram (AL)", 
       cex=1.2, hang=-1)

# taglio a k=4 gruppi
rect.hclust(agnes.average, k=4, border=3:6)
# k=4
res<-cutree(agnes.average, 4)
res
cnames<-colnames(euro.mat)
# estraiamo le città (i nomi) che appartengono ai clusters 
# 1, 2, 3 e 4
cnames[res==1]
cnames[res==2]
cnames[res==3]
cnames[res==4]

## per metodo di Ward si utilizza agnes()
# con method="ward"

agnes.ward<-agnes(eurodist, method = "ward")
pltree(agnes.ward, hang=-1, cex=1, main="AHC-Ward")

# confronto tra "res" e la soluzione in k=4 con Ward
res2<-cutree(agnes.ward, 4)
adj.rand.index(res, res2)
#?adj.rand.index


# ESERCIZIO
# Si consideri il dataset ‘world.cities’ nel pacchetto maps
# Si selezionano le capitali nel dataset e si estraggono 
# le colonne relative a latitudine (lat) e longitudine (long).

library(maps)
data(world.cities)
capital<-subset(world.cities,capital==1)
d0<-capital[ , c(4,5)]

# 1. effettuare AHC con average linkage e distanza euclidea 
# 2. tagliare il dendrogramma a 4 e 9 gruppi, ottenendo i due
# raggruppamenti
# 3. ottenere una mappa dei gruppi
# 4. ottenere la silhouette

ris.agnes<-agnes(d0, method="average")
# dendrogram
par(mar=c(2,2,2,2))
pltree(ris.agnes, main="", cex=0.3)

# taglio con k=4 gruppi
gruppi_4<-cutree(ris.agnes,k=4)
#table(gruppi_4)
# taglio con k=9 gruppi
gruppi_9<-cutree(ris.agnes,k=9)

# Soluzione con K=4
par(mar=c(0,0,0,0))
map("world", interior=FALSE)
points(d0[,2],d0[,1],col=gruppi_4,pch=gruppi_4)
# Soluzione con K=9
par(mar=c(0,0,0,0))
map("world",interior=FALSE)
points(d0[,2],d0[,1],col=gruppi_9,pch=gruppi_9)

# otteniamo l'indice di silhouette media al variare di 
# del numero di cluster K

ris.agnes<-agnes(d0, metric = "euclidean", 
                 method = "average")
asw<-c()
kmin<-2
kmax<-10
ran.k<-length(c(kmin:kmax))

ED<-dist(d0)  # dist. euclidee per la funzione
# cluster.stats che calcola l'ASW

for(nc in 1:ran.k){
   gruppi<-cutree(ris.agnes, nc+1)
   stat<-cluster.stats(d=ED, gruppi)
   asw[nc]<-stat$avg.silwidth  
 }

plot(2:10, asw, ylab="asw", xlab="k", type="o", pch=19)
abline(v=4, lty=2)