data(wine,package="rattle")

wine.stand <- scale(wine[-1])  # To standarize the variables

# K-Means
k.means.fit <- kmeans(wine.stand, 3) # k = 3

attributes(k.means.fit)

# Centroids:
k.means.fit$centers

# Clusters:
k.means.fit$cluster

# Cluster size:
k.means.fit$size

wssplot <- function(data, nc=15, seed=1234){
  wss <- (nrow(data)-1)*sum(apply(data,2,var))
  for (i in 2:nc){
    set.seed(seed)
    wss[i] <- sum(kmeans(data, centers=i)$withinss)}
  plot(1:nc, wss, type="b", xlab="Number of Clusters",
       ylab="Within groups sum of squares")}

library(cluster)
clusplot(wine.stand, k.means.fit$cluster, main='2D representation of the Cluster solution',
         color=TRUE, shade=TRUE,
         labels=2, lines=0)


table(wine[,1],k.means.fit$cluster)



##Hierarchical clustering
d <- dist(wine.stand, method = "euclidean") # Euclidean distance matrix.

H.fit <- hclust(d, method="ward")

plot(H.fit) # display dendogram
groups <- cutree(H.fit, k=3) # cut tree into 5 clusters
# draw dendogram with red borders around the 5 clusters
rect.hclust(H.fit, k=3, border="red") 

table(wine[,1],groups)

##Gaussian mixtures

library(mclust)
fit = Mclust(wine.stand, G=4,modelNames="EEE")
show(fit$bic)
p=13
n=nrow(wine.stand)
G=4
show(-2*fit$loglik+log(n)*(G*p+G-1+p*(p+1)/2))
#plot(fit,what="density")

BIC = mclustBIC(wine.stand,modelNames="EEE")
plot(BIC)

fit = Mclust(wine.stand, G=7,modelNames="EEE")
plot(fit,what="classification")

##Principal components
wine.pr = prcomp(wine.stand,retx=TRUE,scale=TRUE)
BIC = mclustBIC(wine.pr$x[,1:2],modelNames="EEE")
plot(BIC)

fit = Mclust(wine.pr$x[,1:2], G=4,modelNames="EEE")
plot(wine.pr$x[,1:2])
plot(fit,what="classification")

BIC = mclustBIC(wine.pr$x[,1:5],modelNames="EEE")
plot(BIC)

fit = Mclust(wine.pr$x[,c(1,5)], G=5,modelNames="EEE")
plot(fit,what="classification")

fit = Mclust(wine.pr$x[,1:2], G=5,modelNames="EEE")
plot(fit,what="classification")

##Density estimation
fit = densityMclust(wine.stand[,1:2],G=1)
plot(fit,what="density",data=wine.stand[,1:2])

fit = densityMclust(wine.stand[,1:2],G=3)
plot(fit,what="density",data=wine.stand[,1:2])

fit = densityMclust(wine.stand[,1:2])
plot(fit,what="density",data=wine.stand[,1:2])

fit = densityMclust(wine.stand[,1:3])
plot(fit,what="density",data=wine.stand[,1:3])