Isomap algorithm demo

isomap.R

require(RColorBrewer); require(igraph); set.seed(100)

make_spiral <- function(step,v,pts) {
x <- y <- c()
a <- b <- 1
for(i in 1:pts){
  theta = step * i
  x <- append(x, (a + b * theta) * cos(theta) + runif(1,-v,v))
  y <- append(y, (a + b * theta) * sin(theta) + runif(1,-v,v))
}
color <- rainbow(pts * 1.2)[1:pts]
plot(y~x, pch=20, col=color, main= 'Spiral')
return(data.frame(x,y,color))
}

nearest <- function(data,K){
  data <- as.matrix(dist(spiral[,1:2]))
  out <- matrix(0, ncol = ncol(data), nrow = nrow(data), dimnames = list(rownames(data),rownames(data)))
  for(i in 1:nrow(data)){
    neighbours <- as.integer(names(sort(data[i,])[2:K+1]))
    out[i,neighbours] = 1
     }
  return(out)
  }

get_graph <- function(s){
  graph <- graph_from_adjacency_matrix(s, weighted =T, mode = 'undirected')
  set_vertex_attr(graph, "label", value = 1:nrow(spiral))
  #plot(graph, vertex.size = 1, vertex.label = NA)
  geo <- distances(graph, algorithm = 'dijkstra')
  return(geo)
}

X = 1000
spiral <- make_spiral(10/X,1,X)
text(spiral[1,1:2], label='A', pos =4)
text(spiral[1000,1:2], label='B', pos=3)

pc_spiral <- prcomp(spiral[,1:2])
plot(data.frame(pc_spiral$x[,1],1), col = as.character(spiral$color),pch=20,
     xlab='Principal Component 1', ylab='', yaxt= 'n', main='PCA of Spiral')
text(data.frame(pc_spiral$x[1,1],1), label='A', pos =3)
text(data.frame(pc_spiral$x[1000,1],1), label='B', pos =1)
nearest_neighbours <- nearest(spiral, 5)
g <- get_graph(nearest_neighbours)
md <- data.frame('scaled' = cmdscale(g,1), 'color' = spiral$color)
plot(data.frame(md$scaled,1), col = as.character(md$color), pch = 20,
     xlab='Dimension 1', ylab='', yaxt= 'n', main='Isomap')
text(md$scaled[1],1, label='A', pos =3)
text(md$scaled[1000],1, label='B', pos =1)