sparse jaccard
[iramuteq] / Rscripts / simi.R
old mode 100644 (file)
new mode 100755 (executable)
index 98bacf7..2a4c8b7
@@ -1,6 +1,6 @@
 #from proxy package
 #############################################################
-#a, b, c, and d are the counts of all (TRUE, TRUE), (TRUE, FALSE), (FALSE, TRUE), and (FALSE, FALSE) 
+#a, b, c, and d are the counts of all (TRUE, TRUE), (TRUE, FALSE), (FALSE, TRUE), and (FALSE, FALSE)
 # n <- a + b + c + d = nrow(x)
 
 make.a <- function(x) {
@@ -31,20 +31,39 @@ my.jaccard <- function(x) {
     a <- make.a(x)
     b <- make.b(x)
     c <- make.c(x)
-    d <- make.d(x, a, b, c)
+    #d <- make.d(x, a, b, c)
     jac <- a / (a + b + c)
     jac
 }
 
+#Col-wise Jaccard similarity
+#http://stats.stackexchange.com/a/89947/2817
+sparse.jaccard <- function(x) {
+  A = crossprod(x)
+  ix = which(A > 0, arr.ind=TRUE)
+  b = colSums(x)
+  Aix = A[ix]
+  J = sparseMatrix(
+    i = ix[,1],
+    j = ix[,2],
+    x = Aix / (b[ix[,1]] + b[ix[,2]] - Aix),
+    dims = dim(A)
+  )
+  colnames(J) <- colnames(x)
+  rownames(J) <- row.names(x)
+  return(J)
+}
+
+
 
 prcooc <- function(x, a) {
-    prc <- (a / nrow(x)) 
+    prc <- (a / nrow(x))
     prc
 }
 
 make.bin <- function(cs, a, i, j, nb) {
     if (a[i, j] >= 1) {
-        ab <- a[i, j] - 1 
+        ab <- a[i, j] - 1
         res <- binom.test(ab, nb, (cs[i]/nb) * (cs[j]/nb), "less")
     } else {
         res <- NULL
@@ -114,7 +133,7 @@ BuildProf01<-function(x,classes) {
        mat
 }
 
-do.simi <- function(x, method = 'cooc',seuil = NULL, p.type = 'tkplot',layout.type = 'frutch', max.tree = TRUE, coeff.vertex=NULL, coeff.edge = NULL, minmaxeff=c(NULL,NULL), vcexminmax= c(NULL,NULL), cex = 1, coords = NULL) {
+do.simi <- function(x, method = 'cooc',seuil = NULL, p.type = 'tkplot',layout.type = 'frutch', max.tree = TRUE, coeff.vertex=NULL, coeff.edge = NULL, minmaxeff=c(NULL,NULL), vcexminmax= c(NULL,NULL), cex = 1, coords = NULL, communities = NULL, halo = FALSE, fromcoords=NULL, forvertex=NULL, index.word=NULL) {
        mat.simi <- x$mat
     mat.eff <- x$eff
     v.label <- colnames(mat.simi)
@@ -122,28 +141,38 @@ do.simi <- function(x, method = 'cooc',seuil = NULL, p.type = 'tkplot',layout.ty
        g.toplot<-g1
        weori<-get.edge.attribute(g1,'weight')
        if (max.tree) {
-               invw<-1/weori
+        if (method == 'cooc') {
+                   invw <- 1 / weori
+        } else {
+            invw <- 1 - weori
+        }
                E(g1)$weight<-invw
                g.max<-minimum.spanning.tree(g1)
-               E(g.max)$weight<-1/E(g.max)$weight
+        if (method == 'cooc') {
+                   E(g.max)$weight<-1 / E(g.max)$weight
+        } else {
+            E(g.max)$weight<-1 - E(g.max)$weight
+        }
                g.toplot<-g.max
        }
 
     if (!is.null(seuil)) {
-        if (seuil >= max(mat.simi)) seuil <- max(mat.simi)-1
+        if (seuil >= max(mat.simi)) seuil <- -Inf
         vec<-vector()
         w<-E(g.toplot)$weight
         tovire <- which(w<=seuil)
-        g.toplot <- delete.edges(g.toplot,(tovire-1))
-        for (i in 0:(length(V(g.toplot))-1)) {
+        g.toplot <- delete.edges(g.toplot,(tovire))
+        for (i in 1:(length(V(g.toplot)))) {
             if (length(neighbors(g.toplot,i))==0) {
                 vec<-append(vec,i)
             }
         }
         g.toplot <- delete.vertices(g.toplot,vec)
         v.label <- V(g.toplot)$name
-        if (!is.logical(vec)) mat.eff <- mat.eff[-(vec+1)]
-    }
+        if (!is.logical(vec)) mat.eff <- mat.eff[-(vec)]
+    } else {
+               vec <- NULL
+       }
 
        if (!is.null(minmaxeff[1])) {
         eff<-norm.vec(mat.eff,minmaxeff[1],minmaxeff[2])
@@ -156,41 +185,93 @@ do.simi <- function(x, method = 'cooc',seuil = NULL, p.type = 'tkplot',layout.ty
         label.cex = cex
     }
     if (!is.null(coeff.edge)) {
-        we.width <- norm.vec(abs(E(g.toplot)$weight), coeff.edge[1], coeff.edge[2]) 
+        #FIXME
+        we.width <- norm.vec(abs(E(g.toplot)$weight), coeff.edge[1], coeff.edge[2])
            #we.width <- abs((E(g.toplot)$weight/max(abs(E(g.toplot)$weight)))*coeff.edge)
     } else {
         we.width <- NULL
     }
     if (method != 'binom') {
-        we.label <- round(E(g.toplot)$weight,1)
-    } else {
         we.label <- round(E(g.toplot)$weight,3)
+    } else {
+        we.label <- round(E(g.toplot)$weight,4)
     }
-       if (p.type=='rgl') {
+       if (p.type=='rgl' || p.type=='rglweb') {
         nd<-3
     } else {
         nd<-2
     }
+    if (! is.null(fromcoords)) {
+        newfrom <- matrix(runif(nd*length(V(g.toplot)$name),min(fromcoords)),max(fromcoords),ncol=nd, nrow=length(V(g.toplot)$name))
+        for (i in 1:length(V(g.toplot)$name)) {
+            if(V(g.toplot)$name[i] %in% forvertex) {
+                newfrom[i,] <- fromcoords[which(forvertex==V(g.toplot)$name[i]),]
+            }
+        }
+       fromcoords <- newfrom
+    }
+    #print(layout.type)
     if (is.null(coords)) {
-       if (layout.type == 'frutch')
-               lo <- layout.fruchterman.reingold(g.toplot,dim=nd)#, weightsA=E(g.toplot)$weight)
-       if (layout.type == 'kawa')
-               lo <- layout.kamada.kawai(g.toplot,dim=nd)
+       if (layout.type == 'frutch') {
+            #lo <- layout_with_drl(g.toplot,dim=nd)
+            #lo <- layout_with_fr(g.toplot,dim=nd, grid="grid", niter=10000, weights=1/E(g.toplot)$weight)#, start.temp = 1)#, )
+                       if (nd==2) {
+                               library(sna)
+                               library(intergraph)
+                               lo <- gplot.layout.fruchtermanreingold(asNetwork(g.toplot), list())
+                               detach("package:intergraph", unload=TRUE)
+                               detach("package:sna", unload=TRUE)
+                               detach("package:network", unload=TRUE)
+                               library(igraph)
+                       } else {
+                               lo <- layout_with_fr(g.toplot,dim=nd)
+                       }
+        }
+       if (layout.type == 'kawa') {
+               lo <- layout_with_kk(g.toplot,dim=nd, weights=1/E(g.toplot)$weight, start=fromcoords, epsilon=0, maxiter = 10000)
+            #print(lo)
+        }
        if (layout.type == 'random')
-               lo <- layout.random(g.toplot,dim=nd)
+               lo <- layout_on_grid(g.toplot,dim=nd)
        if (layout.type == 'circle' & p.type != 'rgl')
-               lo <- layout.circle(g.toplot)
+               lo <- layout_in_circle(g.toplot)
        if (layout.type == 'circle' & p.type == 'rgl')
-               lo <- layout.sphere(g.toplot)
+               lo <- layout_on_sphere(g.toplot)
         if (layout.type == 'graphopt')
-            lo <- layout.graphopt(g.toplot)
+            lo <- layout_as_tree(g.toplot, circular = TRUE)
+               if (layout.type == 'spirale')
+                       lo <- spirale(g.toplot, E(g.toplot)$weight, index.word)
+               if (layout.type == 'spirale3D')
+                       lo <- spirale3D(g.toplot, E(g.toplot)$weight, index.word)
     } else {
         lo <- coords
     }
-       out <- list(graph = g.toplot, mat.eff = mat.eff, eff = eff, mat = mat.simi, v.label = v.label, we.width = we.width, we.label=we.label, label.cex = label.cex, layout = lo)
+    if (!is.null(communities)) {
+        if (communities == 0 ){
+            com <- edge.betweenness.community(g.toplot)
+        } else if (communities == 1) {
+            com <- fastgreedy.community(g.toplot)
+        } else if (communities == 2) {
+            com <- label.propagation.community(g.toplot)
+        } else if (communities == 3) {
+            com <- leading.eigenvector.community(g.toplot)
+        } else if (communities == 4) {
+            com <- multilevel.community(g.toplot)
+        } else if (communities == 5) {
+            com <- optimal.community(g.toplot)
+        } else if (communities == 6) {
+            com <- spinglass.community(g.toplot)
+        } else if (communities == 7) {
+            com <- walktrap.community(g.toplot)
+        }
+    } else {
+        com <- NULL
+    }
+
+       out <- list(graph = g.toplot, mat.eff = mat.eff, eff = eff, mat = mat.simi, v.label = v.label, we.width = we.width, we.label=we.label, label.cex = label.cex, layout = lo, communities = com, halo = halo, elim=vec)
 }
-       
-plot.simi <- function(graph.simi, p.type = 'tkplot',filename=NULL, vertex.col = 'red', edge.col = 'black', edge.label = TRUE, vertex.label=TRUE, vertex.label.color = 'black', vertex.label.cex= NULL, vertex.size=NULL, leg=NULL, width = 800, height = 800, alpha = 0.1, cexalpha = FALSE, movie = NULL) {
+
+plot.simi <- function(graph.simi, p.type = 'tkplot',filename=NULL, communities = NULL, vertex.col = 'red', edge.col = 'black', edge.label = TRUE, vertex.label=TRUE, vertex.label.color = 'black', vertex.label.cex= NULL, vertex.size=NULL, leg=NULL, width = 800, height = 800, alpha = 0.1, cexalpha = FALSE, movie = NULL, edge.curved = TRUE, svg = FALSE, bg='white') {
        mat.simi <- graph.simi$mat
        g.toplot <- graph.simi$graph
     if (is.null(vertex.size)) {
@@ -211,11 +292,13 @@ plot.simi <- function(graph.simi, p.type = 'tkplot',filename=NULL, vertex.col =
         we.label <- NA
     }
        lo <- graph.simi$layout
+    #rownames(lo) <- v.label
     if (!is.null(vertex.label.cex)) {
         label.cex<-vertex.label.cex
     } else {
         label.cex = graph.simi$label.cex
     }
     if (cexalpha) {
         alphas <- norm.vec(label.cex, 0.5,1)
         nvlc <- NULL
@@ -232,14 +315,24 @@ plot.simi <- function(graph.simi, p.type = 'tkplot',filename=NULL, vertex.col =
     }
     if (p.type=='nplot') {
         #print('ATTENTION - PAS OPEN FILE')
-        open_file_graph(filename, width = width, height = height)
+        open_file_graph(filename, width = width, height = height, svg = svg)
         par(mar=c(2,2,2,2))
+        par(bg=bg)
         if (!is.null(leg)) {
             layout(matrix(c(1,2),1,2, byrow=TRUE),widths=c(3,lcm(7)))
             par(mar=c(2,2,1,0))
         }
         par(pch=' ')
-        plot(g.toplot,vertex.label='', edge.width=we.width, vertex.size=vertex.size, vertex.color=vertex.col, vertex.label.color='white', edge.label=we.label, edge.label.cex=cex, edge.color=edge.col, vertex.label.cex = 0, layout=lo)#, rescale = FALSE)
+        if (is.null(graph.simi$com)) {
+            plot(g.toplot,vertex.label='', edge.width=we.width, vertex.size=vertex.size, vertex.color=vertex.col, vertex.label.color='white', edge.label=we.label, edge.label.cex=cex, edge.color=edge.col, vertex.label.cex = 0, layout=lo, edge.curved=edge.curved)#, rescale = FALSE)
+        } else {
+            if (graph.simi$halo) {
+                mark.groups <- communities(graph.simi$com)
+            } else {
+                mark.groups <- NULL
+            }
+            plot(com, g.toplot,vertex.label='', edge.width=we.width, vertex.size=vertex.size, vertex.color=vertex.col, vertex.label.color='white', edge.label=we.label, edge.label.cex=cex, edge.color=edge.col, vertex.label.cex = 0, layout=lo, mark.groups = mark.groups, edge.curved=edge.curved)
+        }
         #txt.layout <- lo
         txt.layout <- layout.norm(lo, -1, 1, -1, 1, -1, 1)
         #txt.layout <- txt.layout[order(label.cex),]
@@ -267,14 +360,18 @@ plot.simi <- function(graph.simi, p.type = 'tkplot',filename=NULL, vertex.col =
     return(coords)
        }
        
-       if (p.type == 'rgl') {
+       if (p.type == 'rgl' || p.type == 'rglweb') {
                library('rgl')
         #rgl.open()
         #par3d(cex=0.8)
-               rglplot(g.toplot,vertex.label= vire.nonascii(v.label), edge.width=we.width/10, vertex.size=0.01, vertex.color=vertex.col, vertex.label.color="black", edge.color = edge.col, layout=lo)
-        los <- layout.norm(lo, -1, 1, -1, 1, -1, 1)
-        rgl.spheres(los, col = vertex.col, radius = vertex.size/100, alpha = alpha)
-               rgl.bg(color = c('white','black'))
+        lo <- layout.norm(lo, -10, 10, -10, 10, -10, 10)
+               bg3d('white')
+               rglplot(g.toplot,vertex.label='', edge.width=we.width/10, vertex.size=0.01, vertex.color=vertex.col, vertex.label.color="black", edge.color = edge.col, layout=lo, rescale = FALSE)
+        #los <- layout.norm(lo, -1, 1, -1, 1, -1, 1)
+               text3d(lo[,1], lo[,2], lo[,3], vire.nonascii(v.label), col = vertex.label.color, alpha = 1, cex = vertex.label.cex)
+        rgl.spheres(lo, col = vertex.col, radius = vertex.size/100, alpha = alpha)
+        #rgl.bg(color = c('white','black'))
+        #bg3d('white')
         if (!is.null(movie)) {
             require(tcltk)
             ReturnVal <- tkmessageBox(title="RGL 3 D",message="Cliquez pour commencer le film",icon="info",type="ok")
@@ -283,11 +380,26 @@ plot.simi <- function(graph.simi, p.type = 'tkplot',filename=NULL, vertex.col =
             ReturnVal <- tkmessageBox(title="RGL 3 D",message="Film fini !",icon="info",type="ok")
         }
         #play3d(spin3d(axis=c(0,1,0),rpm=6))
-        require(tcltk)
-        ReturnVal <- tkmessageBox(title="RGL 3 D",message="Cliquez pour fermer",icon="info",type="ok")
+        if (p.type == 'rglweb') {
+            writeWebGL(dir = filename, width = width, height= height)
+                       #rglwidget()
+        }# else {
+            require(tcltk)
+            ReturnVal <- tkmessageBox(title="RGL 3 D",message="Cliquez pour fermer",icon="info",type="ok")
+        #}
         rgl.close()
        #       while (rgl.cur() != 0)
        #               Sys.sleep(0.5)
+       } else if (p.type == 'web') {
+               library(rgexf)
+        graph.simi$label.cex <- label.cex
+               if (length(vertex.col)==1) {
+                       vertex.col <- rep(vertex.col, length(v.label))
+               }
+        graph.simi$color <- vertex.col
+        label <- v.label
+        nodes.attr <- data.frame(label)
+               simi.to.gexf(filename, graph.simi, nodes.attr = nodes.attr)
        }
 }
 
@@ -298,3 +410,254 @@ graph.word <- function(mat.simi, index) {
     nm[index,] <- mat.simi[index,]
     nm
 }
+
+#from : 
+#http://gopalakrishna.palem.in/iGraphExport.html#GexfExport
+# Converts the given igraph object to GEXF format and saves it at the given filepath location
+#     g: input igraph object to be converted to gexf format
+#     filepath: file location where the output gexf file should be saved
+#
+saveAsGEXF = function(g, filepath="converted_graph.gexf")
+{
+  require(igraph)
+  require(rgexf)
+  
+  # gexf nodes require two column data frame (id, label)
+  # check if the input vertices has label already present
+  # if not, just have the ids themselves as the label
+  if(is.null(V(g)$label))
+    V(g)$label <- as.character(V(g))
+  
+  # similarily if edges does not have weight, add default 1 weight
+  if(is.null(E(g)$weight))
+    E(g)$weight <- rep.int(1, ecount(g))
+  
+  nodes <- data.frame(cbind(1:vcount(g), V(g)$label))
+  nodes[,1] <- as.character(nodes[,1])
+  nodes[,2] <- as.character(nodes[,2])
+  edges <- t(Vectorize(get.edge, vectorize.args='id')(g, 1:ecount(g)))
+  
+  # combine all node attributes into a matrix (and take care of & for xml)
+  vAttrNames <- setdiff(list.vertex.attributes(g), "label")
+  for (val in c("x","y","color")) {
+        vAttrNames <- setdiff(vAttrNames, val)
+  }
+  nodesAtt <- data.frame(sapply(vAttrNames, function(attr) sub("&", "&",get.vertex.attribute(g, attr))))
+  for (i in 1:ncol(nodesAtt)) {
+      nodesAtt[,i] <- as.character(nodesAtt[,i])
+  }
+  
+  # combine all edge attributes into a matrix (and take care of & for xml)
+  eAttrNames <- setdiff(list.edge.attributes(g), "weight") 
+  edgesAtt <- data.frame(sapply(eAttrNames, function(attr) sub("&", "&",get.edge.attribute(g, attr))))
+  
+  # combine all graph attributes into a meta-data
+  graphAtt <- sapply(list.graph.attributes(g), function(attr) sub("&", "&",get.graph.attribute(g, attr)))
+  ll <- length(V(g)$x)
+  cc <- t(sapply(V(g)$color, col2rgb, alpha=TRUE))
+  cc[,4] <- cc[,4]/255
+  # generate the gexf object
+  output <- write.gexf(nodes, edges, 
+                       edgesWeight=E(g)$weight,
+                       edgesAtt = edgesAtt,
+                       #edgesVizAtt = list(size=as.matrix(E(g)$weight)),
+                       nodesAtt = nodesAtt,
+                       nodesVizAtt=list(color=cc, position=cbind(V(g)$x,V(g)$y, rep(0,ll)), size=V(g)$weight),
+                       meta=c(list(creator="iramuteq", description="igraph -> gexf converted file", keywords="igraph, gexf, R, rgexf"), graphAtt))
+  
+  print(output, filepath, replace=T)
+}
+
+
+merge.graph <- function(graphs) {
+    library(colorspace)
+    ng <- graph.union(graphs, byname=T)
+    V.weight <- V(ng)$weight_1 
+    E.weight <- E(ng)$weight_1
+    cols <- rainbow(length(graphs))
+    V.color <- rep(cols[1], length(V.weight))
+    for (i in 2:length(graphs)) {
+        tw <- paste('weight_', i, sep='')
+        tocomp <- get.vertex.attribute(ng,tw)
+        totest <- intersect(which(!is.na(V.weight)), which(!is.na(tocomp)))
+        maxmat <- cbind(V.weight[totest], tocomp[totest])
+        resmax <- apply(maxmat, 1, which.max)
+        ncolor <- c(cols[(i-1)], cols[i])
+        #rbgcol1 <- col2rgb(cols[(i-1)])
+        #rbgcol1 <- rbgcol1/255
+        #rgbcol1 <- RGB(rbgcol1[1],rbgcol1[2],rbgcol1[3])
+        rbgcol2 <- col2rgb(cols[i])
+        rbgcol2 <- rbgcol2/255
+        #rgbcol2 <- RGB(rbgcol2[1],rbgcol2[2],rbgcol2[3])       
+        for (j in totest) {
+            alpha <- tocomp[j] /(V.weight[j] + tocomp[j])
+            rbgcol1 <- col2rgb(V.color[j])
+            rbgcol1 <- rbgcol1/255
+            #mix.col <- mixcolor(alpha,rbgcol1, rbgcol2)
+            mix.col <- mixcolor(alpha, RGB(rbgcol1[1],rbgcol1[2],rbgcol1[3]), RGB(rbgcol2[1],rbgcol2[2],rbgcol2[3]))
+            V.color[j] <- hex(mix.col)
+                       #V.color[j] <- adjustcolor(hex(mix.col), 0.6)
+        }
+        #to.change <- totest[which(resmax == 2)]
+        #V.color[to.change] <- cols[i]
+        V.weight[totest] <- apply(maxmat, 1, max)
+        nas <- which(is.na(V.weight))
+        nas2 <- which(is.na(tocomp))
+        fr2 <- setdiff(nas,nas2)
+        V.weight[fr2] <- tocomp[fr2]
+        V.color[fr2] <- cols[i]
+        tocomp <- get.edge.attribute(ng, tw)
+        totest <- intersect(which(!is.na(E.weight)), which(!is.na(tocomp)))
+        maxmat <- cbind(E.weight[totest], tocomp[totest])
+        resmax <- apply(maxmat, 1, which.max)
+        E.weight[totest] <- apply(maxmat, 1, max)
+        nas <- which(is.na(E.weight))
+        nas2 <- which(is.na(tocomp))
+        fr2 <- setdiff(nas,nas2)
+        E.weight[fr2] <- tocomp[fr2]        
+    }
+    V(ng)$weight <- V.weight
+    V(ng)$color <- V.color
+    E(ng)$weight <- E.weight
+       colors <- col2rgb(V(ng)$color)
+       V(ng)$r <- colors["red", ]
+       V(ng)$g <- colors["green", ]
+       V(ng)$b <- colors["blue", ]
+    ng
+}
+
+merge.graph.proto <- function(graphs) {
+    library(colorspace)
+    ng <- graph.union(graphs, byname=T)
+    V.weight <- V(ng)$weight_1 
+    E.weight <- E(ng)$weight_1
+       V.proto.color <- V(ng)$proto.color_1
+       cols <- rainbow(length(graphs))
+       V.color <- rep(cols[1], length(V.weight))
+    for (i in 2:length(graphs)) {
+        tw <- paste('weight_', i, sep='')
+        tocomp <- get.vertex.attribute(ng,tw)
+        totest <- intersect(which(!is.na(V.weight)), which(!is.na(tocomp)))
+        maxmat <- cbind(V.weight[totest], tocomp[totest])
+        resmax <- apply(maxmat, 1, which.max)
+        V.weight[totest] <- apply(maxmat, 1, max)
+        nas <- which(is.na(V.weight))
+        nas2 <- which(is.na(tocomp))
+        fr2 <- setdiff(nas,nas2)
+        V.weight[fr2] <- tocomp[fr2]
+
+               cw <- paste('proto.color_', i, sep='')
+               tocomp.col <- get.vertex.attribute(ng,cw)
+               which.sup <- which(resmax==2)
+               V.proto.color[totest[which.sup]] <- tocomp.col[totest[which.sup]]
+               V.proto.color[fr2] <- tocomp.col[fr2]
+
+               V.color[totest[which.sup]] <- cols[i]
+               V.color[fr2] <- cols[i]
+
+        tocomp <- get.edge.attribute(ng, tw)
+        totest <- intersect(which(!is.na(E.weight)), which(!is.na(tocomp)))
+        maxmat <- cbind(E.weight[totest], tocomp[totest])
+        resmax <- apply(maxmat, 1, which.max)
+        E.weight[totest] <- apply(maxmat, 1, max)
+        nas <- which(is.na(E.weight))
+        nas2 <- which(is.na(tocomp))
+        fr2 <- setdiff(nas,nas2)
+        E.weight[fr2] <- tocomp[fr2]           
+       }
+    V(ng)$weight <- V.weight
+    V(ng)$proto.color <- V.proto.color
+       V(ng)$color <- V.proto.color
+    E(ng)$weight <- E.weight
+       V(ng)$ocolor <- V.color
+       colors <- col2rgb(V(ng)$color)
+       V(ng)$r <- colors["red", ]
+       V(ng)$g <- colors["green", ]
+       V(ng)$b <- colors["blue", ]
+    ng
+}
+
+
+spirale <- function(g, weigth, center, miny=0.1) {
+       ncoord <- matrix(0, nrow=length(weigth)+1, ncol=2)
+       v.names <- V(g)$name 
+       center.name <- v.names[center]
+       first <- which.max(weigth)[1]
+       if (head_of(g, first)$name == center.name) {
+               n.name <- tail_of(g, first)
+       } else {
+               n.name <- head_of(g, first)
+       }
+       n.name <- n.name$name
+       nb <- length(weigth)
+       ncoord[which(v.names==n.name),] <- c(0,1)
+       weigth[first] <- 0
+       rs <- norm.vec(weigth,1, miny)
+       nbt <- nb %/% 50
+       if (nbt == 0) nbt <- 1
+       angler <- ((360 * nbt) / (nb- 1)) * (pi/180)
+       ang <- 90 * (pi/180)
+       rr <- (1-miny) / (nb-1)
+       r <- 1
+       while (max(weigth != 0)) {
+               first <- which.max(weigth)[1]
+               if (head_of(g, first)$name == center.name) {
+                       n.name <- tail_of(g, first)
+               } else {
+                       n.name <- head_of(g, first)
+               }
+               n.name <- n.name$name
+               #r <- rs[first]
+               r <- r - rr
+               ang <- ang + angler 
+               x <- r * cos(ang)
+               y <- r * sin(ang)
+               weigth[first] <- 0
+               ncoord[which(v.names==n.name),] <- c(x,y)
+       }
+       ncoord
+}
+
+spirale3D <- function(g, weigth, center, miny=0.1) {
+       ncoord <- matrix(0, nrow=length(weigth)+1, ncol=3)
+       v.names <- V(g)$name 
+       center.name <- v.names[center]
+       first <- which.max(weigth)[1]
+       if (head_of(g, first)$name == center.name) {
+               n.name <- tail_of(g, first)
+       } else {
+               n.name <- head_of(g, first)
+       }
+       n.name <- n.name$name
+       nb <- length(weigth)
+       ncoord[which(v.names==n.name),] <- c(0,0,1)
+       weigth[first] <- 0
+       rs <- norm.vec(weigth,1, miny)
+       nbt <- nb %/% 50
+       if (nbt == 0) nbt <- 1
+       angler <- ((360 * nbt) / (nb- 1)) * (pi/180)
+       theta <- 0
+       phi <- 90 * (pi/180)
+       rr <- (1-miny) / (nb-1)
+       r <- 1
+       while (max(weigth != 0)) {
+               first <- which.max(weigth)[1]
+               if (head_of(g, first)$name == center.name) {
+                       n.name <- tail_of(g, first)
+               } else {
+                       n.name <- head_of(g, first)
+               }
+               n.name <- n.name$name
+               #r <- rs[first]
+               r <- r - rr
+               theta <- theta + angler
+               phi <- phi + angler/2
+               x <- r * sin(theta) * cos(phi)
+               y <- r * sin(theta) * sin(phi)
+               z <- r * cos(theta)
+               weigth[first] <- 0
+               ncoord[which(v.names==n.name),] <- c(x,y,z)
+       }
+       ncoord
+}
+