####
#
# Animate the k-means algorithm (Amesian code)

library(cluster)

data(ruspini)

x<-ruspini
par(mfrow=c(1,1))
stepKMean(x, 3)
stepKMean(x, 4)
stepKMean(x, 5)


toCluster <- function(x) {
   return(which.min(as.matrix(dist(rbind(x,centers)))[1,-1]))
}


stepKMean <- function(x,K) {
   centers <<- vector()
   centers <<- x[sample(1:dim(x)[[1]], K), , drop = FALSE]
   plot(x)
   points(centers,col=c(2:(K+1)),pch="x")

   ch <- scan()

   while(length(ch) == 0) {
      cluster <- apply(x, MARGIN=1,FUN=toCluster)
      if (dim(x)[2]>2)
         plot(x,col=cluster+1)
      else 
         points(x,col=cluster+1)

      ch <- scan()

      centersNew <- t(sapply(1:K,FUN=function(i){ return(apply(x[cluster==i,],MARGIN=2,FUN=mean))}))

      points(centers,col=1,pch="x")
      points(centersNew,col=c(2:(K+1)),pch="x")
      for(i in 1:K) {
         lines(c(centers[i,1],centersNew[i,1]),c(centers[i,2],centersNew[i,2]),col=i+1)
      }

      centers <<- centersNew
      ch <- scan()
   }

   result <- list()
   result[[1]] <- cluster
   result[[2]] <- centers
   names(result) <- c("cluster","centers")

   return(result)
}