Laplacian eigenmap demo

eigenmap.R

require('geigen'); 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
  }
  for(j in 1:nrow(out)){
    for(k in 1:ncol(out)){
      if(out[j,k] == 1){
        out[k,j] = 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)
}

spiral <- make_spiral(0.02,1,1000)
text(spiral[1,1:2], label='A', pos =4)
text(spiral[1000,1:2], label='B', pos=2)
nearest_neighbours <- nearest(spiral, 8)
g <- get_graph(nearest_neighbours)
W <- nearest_neighbours 
D <- diag(nrow(W))

for(i in 1:nrow(D)){
  D[i,i]=as.vector(colSums(W))[i]
}

L = D - W
eL <- geigen(L,D)
eL$values[1:10]
embedding <- eL$vectors[,1:2]
plot(-data.frame(embedding[,1],1) ,col =as.character(spiral$color),pch=20,
     xlab='Dimension 1', ylab='', yaxt= 'n', main='Laplacian Eigenmapping', xaxt='n')
text(-data.frame(embedding[1,1],1), label='A', pos =3)
text(-data.frame(embedding[1000,1],1), label='B', pos =1)