...
[iramuteq] / Rscripts / simi.R
1 #from proxy package
2 #############################################################
3 #a, b, c, and d are the counts of all (TRUE, TRUE), (TRUE, FALSE), (FALSE, TRUE), and (FALSE, FALSE)
4 # n <- a + b + c + d = nrow(x)
5
6 make.a <- function(x) {
7     a  <- t(x) %*% (x)
8     a
9 }
10
11 make.b <- function(x) {
12     b <- t(x) %*% (1 - x)
13     b
14 }
15
16 make.c <- function(x) {
17     c <- (1-t(x)) %*% x
18     c
19 }
20
21 make.d <- function(x, a, b, c) {
22 #??????????? ncol ?
23     d <- ncol(x) - a - b - c
24     d
25 }
26
27 ###########################################
28 #x, a
29 ###########################################
30 my.jaccard <- function(x) {
31     a <- make.a(x)
32     b <- make.b(x)
33     c <- make.c(x)
34     #d <- make.d(x, a, b, c)
35     jac <- a / (a + b + c)
36     jac
37 }
38
39 #Col-wise Jaccard similarity
40 #http://stats.stackexchange.com/a/89947/2817
41 sparse.jaccard <- function(x) {
42   A = crossprod(x)
43   ix = which(A > 0, arr.ind=TRUE)
44   b = colSums(x)
45   Aix = A[ix]
46   J = sparseMatrix(
47     i = ix[,1],
48     j = ix[,2],
49     x = Aix / (b[ix[,1]] + b[ix[,2]] - Aix),
50     dims = dim(A)
51   )
52   colnames(J) <- colnames(x)
53   rownames(J) <- row.names(x)
54   return(J)
55 }
56
57
58
59 prcooc <- function(x, a) {
60     prc <- (a / nrow(x))
61     prc
62 }
63
64 make.bin <- function(cs, a, i, j, nb) {
65     if (a[i, j] >= 1) {
66         ab <- a[i, j] - 1
67         res <- binom.test(ab, nb, (cs[i]/nb) * (cs[j]/nb), "less")
68     } else {
69         res <- NULL
70         res$p.value <- 0
71     }
72     #res <- binom.test(ab, nb, (cs[i]/nb) * (cs[j]/nb), "less")
73     res$p.value
74     }
75
76 binom.sim <- function(x) {
77     a <- make.a(x)
78     n <- nrow(x)
79     cs <- colSums(x)
80     mat <- matrix(0,ncol(x),ncol(x))
81     colnames(mat)<-colnames(a)
82     rownames(mat)<-rownames(a)
83     for (i in 1:(ncol(x)-1)) {
84         for (j in (i+1):ncol(x)) {
85             mat[j,i] <- make.bin(cs, a, i, j , n)
86         }
87     }
88 #    print(mat)
89     mat
90 }
91
92
93 ############################################
94 # a, b, c
95 ############################################
96 # jaccard a, b, c   a / (a + b + c)
97 # Kulczynski1 a, b, c   a / (b + c)
98 # Kulczynski2 a, b, c   [a / (a + b) + a / (a + c)] / 2
99 # Mountford a, b, c    2a / (ab + ac + 2bc)
100 # Fager, McGowan a, b, c   a / sqrt((a + b)(a + c)) - 1 / 2 sqrt(a + c)
101 # Russel, Rao a (a/n)
102 # Dice, Czekanowski, Sorensen a, b, c   2a / (2a + b + c)
103 # Mozley, Margalef a, b, c  an / (a + b)(a + c)
104 # Ochiai a, b, c  a / sqrt[(a + b)(a + c)]
105 # Simpson a, b, c   a / min{(a + b), (a + c)}
106 # Braun-Blanquet a, b, c  a / max{(a + b), (a + c)}
107
108 #simple matching, Sokal/Michener a, b, c, d, ((a + d) /n)
109 # Hamman, a, b, c, d, ([a + d] - [b + c]) / n
110 # Faith , a, b, c, d, (a + d/2) / n
111 # Tanimoto, Rogers a, b, c, d, (a + d) / (a + 2b + 2c + d)
112 # Phi  a, b, c, d   (ad - bc) / sqrt[(a + b)(c + d)(a + c)(b + d)]
113 # Stiles a, b, c, d  log(n(|ad-bc| - 0.5n)^2 / [(a + b)(c + d)(a + c)(b + d)])
114 # Michael   a, b, c, d   4(ad - bc) / [(a + d)^2 + (b + c)^2]
115 # Yule a, b, c, d  (ad - bc) / (ad + bc)
116 # Yule2  a, b, c, d  (sqrt(ad) - sqrt(bc)) / (sqrt(ad) + sqrt(bc))
117
118 BuildProf01<-function(x,classes) {
119         #x : donnees en 0/1
120         #classes : classes de chaque lignes de x
121         dm<-cbind(x,cl=classes)
122         clnb=length(summary(as.data.frame(as.character(classes)),max=100))
123         print(clnb)
124         print(summary(as.data.frame(as.character(classes)),max=100))
125         mat<-matrix(0,ncol(x),clnb)
126         rownames(mat)<-colnames(x)
127         for (i in 1:clnb) {
128                 dtmp<-dm[which(dm$cl==i),]
129                 for (j in 1:(ncol(dtmp)-1)) {
130                         mat[j,i]<-sum(dtmp[,j])
131                 }
132         }
133         mat
134 }
135
136 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) {
137         mat.simi <- x$mat
138     mat.eff <- x$eff
139     v.label <- colnames(mat.simi)
140         g1<-graph.adjacency(mat.simi,mode="lower",weighted=TRUE)
141         g.toplot<-g1
142         weori<-get.edge.attribute(g1,'weight')
143         if (max.tree) {
144         if (method == 'cooc') {
145                     invw <- 1 / weori
146         } else {
147             invw <- 1 - weori
148         }
149                 E(g1)$weight<-invw
150                 g.max<-minimum.spanning.tree(g1)
151         if (method == 'cooc') {
152                     E(g.max)$weight<-1 / E(g.max)$weight
153         } else {
154             E(g.max)$weight<-1 - E(g.max)$weight
155         }
156                 g.toplot<-g.max
157         }
158
159     if (!is.null(seuil)) {
160         if (seuil >= max(mat.simi)) seuil <- -Inf
161         vec<-vector()
162         w<-E(g.toplot)$weight
163         tovire <- which(w<=seuil)
164         g.toplot <- delete.edges(g.toplot,(tovire))
165         for (i in 1:(length(V(g.toplot)))) {
166             if (length(neighbors(g.toplot,i))==0) {
167                 vec<-append(vec,i)
168             }
169         }
170         g.toplot <- delete.vertices(g.toplot,vec)
171         v.label <- V(g.toplot)$name
172         if (!is.logical(vec)) mat.eff <- mat.eff[-(vec)]
173     } else {
174                 vec <- NULL
175         }
176
177         if (!is.null(minmaxeff[1])) {
178         eff<-norm.vec(mat.eff,minmaxeff[1],minmaxeff[2])
179     } else {
180         eff<-coeff.vertex
181     }
182     if (!is.null(vcexminmax[1])) {
183         label.cex = norm.vec(mat.eff, vcexminmax[1], vcexminmax[2])
184     } else {
185         label.cex = cex
186     }
187     if (!is.null(coeff.edge)) {
188         #FIXME
189         we.width <- norm.vec(abs(E(g.toplot)$weight), coeff.edge[1], coeff.edge[2])
190             #we.width <- abs((E(g.toplot)$weight/max(abs(E(g.toplot)$weight)))*coeff.edge)
191     } else {
192         we.width <- NULL
193     }
194     if (method != 'binom') {
195         we.label <- round(E(g.toplot)$weight,3)
196     } else {
197         we.label <- round(E(g.toplot)$weight,4)
198     }
199         if (p.type=='rgl' || p.type=='rglweb') {
200         nd<-3
201     } else {
202         nd<-2
203     }
204     if (! is.null(fromcoords)) {
205         newfrom <- matrix(runif(nd*length(V(g.toplot)$name),min(fromcoords)),max(fromcoords),ncol=nd, nrow=length(V(g.toplot)$name))
206         for (i in 1:length(V(g.toplot)$name)) {
207             if(V(g.toplot)$name[i] %in% forvertex) {
208                 newfrom[i,] <- fromcoords[which(forvertex==V(g.toplot)$name[i]),]
209             }
210         }
211        fromcoords <- newfrom
212     }
213     #print(layout.type)
214     if (is.null(coords)) {
215         if (layout.type == 'frutch') {
216             #lo <- layout_with_drl(g.toplot,dim=nd)
217             #lo <- layout_with_fr(g.toplot,dim=nd, grid="grid", niter=10000, weights=1/E(g.toplot)$weight)#, start.temp = 1)#, )
218                         if (nd==2) {
219                                 library(sna)
220                                 library(intergraph)
221                                 lo <- gplot.layout.fruchtermanreingold(asNetwork(g.toplot), list())
222                                 detach("package:intergraph", unload=TRUE)
223                                 detach("package:sna", unload=TRUE)
224                                 detach("package:network", unload=TRUE)
225                                 library(igraph)
226                         } else {
227                                 lo <- layout_with_fr(g.toplot,dim=nd)
228                         }
229         }
230         if (layout.type == 'kawa') {
231                 lo <- layout_with_kk(g.toplot,dim=nd, weights=1/E(g.toplot)$weight, start=fromcoords, epsilon=0, maxiter = 10000)
232             #print(lo)
233         }
234         if (layout.type == 'random')
235                 lo <- layout_on_grid(g.toplot,dim=nd)
236         if (layout.type == 'circle' & p.type != 'rgl')
237                 lo <- layout_in_circle(g.toplot)
238         if (layout.type == 'circle' & p.type == 'rgl')
239                 lo <- layout_on_sphere(g.toplot)
240         if (layout.type == 'graphopt')
241             lo <- layout_as_tree(g.toplot, circular = TRUE)
242                 if (layout.type == 'spirale')
243                         lo <- spirale(g.toplot, E(g.toplot)$weight, index.word)
244                 if (layout.type == 'spirale3D')
245                         lo <- spirale3D(g.toplot, E(g.toplot)$weight, index.word)
246     } else {
247         lo <- coords
248     }
249     if (!is.null(communities)) {
250         if (communities == 0 ){
251             com <- edge.betweenness.community(g.toplot)
252         } else if (communities == 1) {
253             com <- fastgreedy.community(g.toplot)
254         } else if (communities == 2) {
255             com <- label.propagation.community(g.toplot)
256         } else if (communities == 3) {
257             com <- leading.eigenvector.community(g.toplot)
258         } else if (communities == 4) {
259             com <- multilevel.community(g.toplot)
260         } else if (communities == 5) {
261             com <- optimal.community(g.toplot)
262         } else if (communities == 6) {
263             com <- spinglass.community(g.toplot)
264         } else if (communities == 7) {
265             com <- walktrap.community(g.toplot)
266         }
267     } else {
268         com <- NULL
269     }
270
271         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)
272 }
273
274 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') {
275         mat.simi <- graph.simi$mat
276         g.toplot <- graph.simi$graph
277     if (is.null(vertex.size)) {
278         vertex.size <- graph.simi$eff
279     } else {
280         vertex.size <- vertex.size
281     }
282         we.width <- graph.simi$we.width
283     if (vertex.label) {
284         #v.label <- vire.nonascii(graph.simi$v.label)
285         v.label <- graph.simi$v.label
286     } else {
287         v.label <- NA
288     }
289     if (edge.label) {
290         we.label <- graph.simi$we.label
291     } else {
292         we.label <- NA
293     }
294         lo <- graph.simi$layout
295     #rownames(lo) <- v.label
296     if (!is.null(vertex.label.cex)) {
297         label.cex<-vertex.label.cex
298     } else {
299         label.cex = graph.simi$label.cex
300     }
301  
302     if (cexalpha) {
303         alphas <- norm.vec(label.cex, 0.5,1)
304         nvlc <- NULL
305         if (length(vertex.label.color) == 1) {
306             for (i in 1:length(alphas)) {
307              nvlc <- append(nvlc, adjustcolor(vertex.label.color, alpha=alphas[i]))
308             }
309         } else {
310             for (i in 1:length(alphas)) {
311                 nvlc <- append(nvlc, adjustcolor(vertex.label.color[i], alpha=alphas[i]))
312             }
313         }
314         vertex.label.color <- nvlc  
315     }
316     if (p.type=='nplot') {
317         #print('ATTENTION - PAS OPEN FILE')
318         open_file_graph(filename, width = width, height = height, svg = svg)
319         par(mar=c(2,2,2,2))
320         par(bg=bg)
321         if (!is.null(leg)) {
322             layout(matrix(c(1,2),1,2, byrow=TRUE),widths=c(3,lcm(7)))
323             par(mar=c(2,2,1,0))
324         }
325         par(pch=' ')
326         if (is.null(graph.simi$com)) {
327             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)
328         } else {
329             if (graph.simi$halo) {
330                 mark.groups <- communities(graph.simi$com)
331             } else {
332                 mark.groups <- NULL
333             }
334             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)
335         }
336         #txt.layout <- lo
337         txt.layout <- layout.norm(lo, -1, 1, -1, 1, -1, 1)
338         #txt.layout <- txt.layout[order(label.cex),]
339         #vertex.label.color <- vertex.label.color[order(label.cex)]
340         #v.label <- v.label[order(label.cex)]
341         #label.cex <- label.cex[order(label.cex)]
342         text(txt.layout[,1], txt.layout[,2], v.label, cex=label.cex, col=vertex.label.color)
343         if (!is.null(leg)) {
344             par(mar=c(0,0,0,0))
345             plot(0, axes = FALSE, pch = '')
346             legend(x = 'center' , leg$unetoile, fill = leg$gcol)
347         }
348         dev.off()
349         return(lo)
350     }
351         if (p.type=='tkplot') {
352                 id <- tkplot(g.toplot,vertex.label=v.label, edge.width=we.width, vertex.size=vertex.size, vertex.color=vertex.col, vertex.label.color=vertex.label.color, edge.label=we.label, edge.color=edge.col, layout=lo)
353         coords = tkplot.getcoords(id)
354         ok <- try(coords <- tkplot.getcoords(id), TRUE)
355                 while (is.matrix(ok)) {
356             ok <- try(coords <- tkplot.getcoords(id), TRUE)
357                         Sys.sleep(0.5)
358         }
359         tkplot.off()
360     return(coords)
361         }
362         
363         if (p.type == 'rgl' || p.type == 'rglweb') {
364                 library('rgl')
365         #rgl.open()
366         #par3d(cex=0.8)
367         lo <- layout.norm(lo, -10, 10, -10, 10, -10, 10)
368                 bg3d('white')
369                 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)
370         #los <- layout.norm(lo, -1, 1, -1, 1, -1, 1)
371                 text3d(lo[,1], lo[,2], lo[,3], vire.nonascii(v.label), col = vertex.label.color, alpha = 1, cex = vertex.label.cex)
372         rgl.spheres(lo, col = vertex.col, radius = vertex.size/100, alpha = alpha)
373         #rgl.bg(color = c('white','black'))
374         #bg3d('white')
375         if (!is.null(movie)) {
376             require(tcltk)
377             ReturnVal <- tkmessageBox(title="RGL 3 D",message="Cliquez pour commencer le film",icon="info",type="ok")
378
379             movie3d(spin3d(axis=c(0,1,0),rpm=6), movie = 'film_graph', frames = "tmpfilm", duration=10, clean=TRUE, top = TRUE, dir = movie)
380             ReturnVal <- tkmessageBox(title="RGL 3 D",message="Film fini !",icon="info",type="ok")
381         }
382         #play3d(spin3d(axis=c(0,1,0),rpm=6))
383         if (p.type == 'rglweb') {
384             writeWebGL(dir = filename, width = width, height= height)
385                         #rglwidget()
386         }# else {
387             require(tcltk)
388             ReturnVal <- tkmessageBox(title="RGL 3 D",message="Cliquez pour fermer",icon="info",type="ok")
389         #}
390         rgl.close()
391         #       while (rgl.cur() != 0)
392         #               Sys.sleep(0.5)
393         } else if (p.type == 'web') {
394                 library(rgexf)
395         graph.simi$label.cex <- label.cex
396                 if (length(vertex.col)==1) {
397                         vertex.col <- rep(vertex.col, length(v.label))
398                 }
399         graph.simi$color <- vertex.col
400         label <- v.label
401         nodes.attr <- data.frame(label)
402                 simi.to.gexf(filename, graph.simi, nodes.attr = nodes.attr)
403         }
404 }
405
406
407 graph.word <- function(mat.simi, index) {
408     nm <- matrix(0, ncol = ncol(mat.simi), nrow=nrow(mat.simi), dimnames=list(row.names(mat.simi), colnames(mat.simi)))
409     nm[,index] <- mat.simi[,index]
410     nm[index,] <- mat.simi[index,]
411     nm
412 }
413
414 #from : 
415 #http://gopalakrishna.palem.in/iGraphExport.html#GexfExport
416 # Converts the given igraph object to GEXF format and saves it at the given filepath location
417 #     g: input igraph object to be converted to gexf format
418 #     filepath: file location where the output gexf file should be saved
419 #
420 saveAsGEXF = function(g, filepath="converted_graph.gexf")
421 {
422   require(igraph)
423   require(rgexf)
424   
425   # gexf nodes require two column data frame (id, label)
426   # check if the input vertices has label already present
427   # if not, just have the ids themselves as the label
428   if(is.null(V(g)$label))
429     V(g)$label <- as.character(V(g))
430   
431   # similarily if edges does not have weight, add default 1 weight
432   if(is.null(E(g)$weight))
433     E(g)$weight <- rep.int(1, ecount(g))
434   
435   nodes <- data.frame(cbind(1:vcount(g), V(g)$label))
436   nodes[,1] <- as.character(nodes[,1])
437   nodes[,2] <- as.character(nodes[,2])
438   edges <- t(Vectorize(get.edge, vectorize.args='id')(g, 1:ecount(g)))
439   
440   # combine all node attributes into a matrix (and take care of & for xml)
441   vAttrNames <- setdiff(list.vertex.attributes(g), "label")
442   for (val in c("x","y","color")) {
443         vAttrNames <- setdiff(vAttrNames, val)
444   }
445   nodesAtt <- data.frame(sapply(vAttrNames, function(attr) sub("&", "&",get.vertex.attribute(g, attr))))
446   for (i in 1:ncol(nodesAtt)) {
447       nodesAtt[,i] <- as.character(nodesAtt[,i])
448   }
449   
450   # combine all edge attributes into a matrix (and take care of & for xml)
451   eAttrNames <- setdiff(list.edge.attributes(g), "weight") 
452   edgesAtt <- data.frame(sapply(eAttrNames, function(attr) sub("&", "&",get.edge.attribute(g, attr))))
453   
454   # combine all graph attributes into a meta-data
455   graphAtt <- sapply(list.graph.attributes(g), function(attr) sub("&", "&",get.graph.attribute(g, attr)))
456   ll <- length(V(g)$x)
457   cc <- t(sapply(V(g)$color, col2rgb, alpha=TRUE))
458   cc[,4] <- cc[,4]/255
459   # generate the gexf object
460   output <- write.gexf(nodes, edges, 
461                        edgesWeight=E(g)$weight,
462                        edgesAtt = edgesAtt,
463                        #edgesVizAtt = list(size=as.matrix(E(g)$weight)),
464                        nodesAtt = nodesAtt,
465                        nodesVizAtt=list(color=cc, position=cbind(V(g)$x,V(g)$y, rep(0,ll)), size=V(g)$weight),
466                        meta=c(list(creator="iramuteq", description="igraph -> gexf converted file", keywords="igraph, gexf, R, rgexf"), graphAtt))
467   
468   print(output, filepath, replace=T)
469 }
470
471
472 merge.graph <- function(graphs) {
473     library(colorspace)
474     ng <- graph.union(graphs, byname=T)
475     V.weight <- V(ng)$weight_1 
476     E.weight <- E(ng)$weight_1
477     cols <- rainbow(length(graphs))
478     V.color <- rep(cols[1], length(V.weight))
479     for (i in 2:length(graphs)) {
480         tw <- paste('weight_', i, sep='')
481         tocomp <- get.vertex.attribute(ng,tw)
482         totest <- intersect(which(!is.na(V.weight)), which(!is.na(tocomp)))
483         maxmat <- cbind(V.weight[totest], tocomp[totest])
484         resmax <- apply(maxmat, 1, which.max)
485         ncolor <- c(cols[(i-1)], cols[i])
486         #rbgcol1 <- col2rgb(cols[(i-1)])
487         #rbgcol1 <- rbgcol1/255
488         #rgbcol1 <- RGB(rbgcol1[1],rbgcol1[2],rbgcol1[3])
489         rbgcol2 <- col2rgb(cols[i])
490         rbgcol2 <- rbgcol2/255
491         #rgbcol2 <- RGB(rbgcol2[1],rbgcol2[2],rbgcol2[3])       
492         for (j in totest) {
493             alpha <- tocomp[j] /(V.weight[j] + tocomp[j])
494             rbgcol1 <- col2rgb(V.color[j])
495             rbgcol1 <- rbgcol1/255
496             #mix.col <- mixcolor(alpha,rbgcol1, rbgcol2)
497             mix.col <- mixcolor(alpha, RGB(rbgcol1[1],rbgcol1[2],rbgcol1[3]), RGB(rbgcol2[1],rbgcol2[2],rbgcol2[3]))
498             V.color[j] <- hex(mix.col)
499                         #V.color[j] <- adjustcolor(hex(mix.col), 0.6)
500         }
501         #to.change <- totest[which(resmax == 2)]
502         #V.color[to.change] <- cols[i]
503         V.weight[totest] <- apply(maxmat, 1, max)
504         nas <- which(is.na(V.weight))
505         nas2 <- which(is.na(tocomp))
506         fr2 <- setdiff(nas,nas2)
507         V.weight[fr2] <- tocomp[fr2]
508         V.color[fr2] <- cols[i]
509         tocomp <- get.edge.attribute(ng, tw)
510         totest <- intersect(which(!is.na(E.weight)), which(!is.na(tocomp)))
511         maxmat <- cbind(E.weight[totest], tocomp[totest])
512         resmax <- apply(maxmat, 1, which.max)
513         E.weight[totest] <- apply(maxmat, 1, max)
514         nas <- which(is.na(E.weight))
515         nas2 <- which(is.na(tocomp))
516         fr2 <- setdiff(nas,nas2)
517         E.weight[fr2] <- tocomp[fr2]        
518     }
519     V(ng)$weight <- V.weight
520     V(ng)$color <- V.color
521     E(ng)$weight <- E.weight
522         colors <- col2rgb(V(ng)$color)
523         V(ng)$r <- colors["red", ]
524         V(ng)$g <- colors["green", ]
525         V(ng)$b <- colors["blue", ]
526     ng
527 }
528
529 merge.graph.proto <- function(graphs) {
530     library(colorspace)
531     ng <- graph.union(graphs, byname=T)
532     V.weight <- V(ng)$weight_1 
533     E.weight <- E(ng)$weight_1
534         V.proto.color <- V(ng)$proto.color_1
535         cols <- rainbow(length(graphs))
536         V.color <- rep(cols[1], length(V.weight))
537     for (i in 2:length(graphs)) {
538         tw <- paste('weight_', i, sep='')
539         tocomp <- get.vertex.attribute(ng,tw)
540         totest <- intersect(which(!is.na(V.weight)), which(!is.na(tocomp)))
541         maxmat <- cbind(V.weight[totest], tocomp[totest])
542         resmax <- apply(maxmat, 1, which.max)
543         V.weight[totest] <- apply(maxmat, 1, max)
544         nas <- which(is.na(V.weight))
545         nas2 <- which(is.na(tocomp))
546         fr2 <- setdiff(nas,nas2)
547         V.weight[fr2] <- tocomp[fr2]
548
549                 cw <- paste('proto.color_', i, sep='')
550                 tocomp.col <- get.vertex.attribute(ng,cw)
551                 which.sup <- which(resmax==2)
552                 V.proto.color[totest[which.sup]] <- tocomp.col[totest[which.sup]]
553                 V.proto.color[fr2] <- tocomp.col[fr2]
554
555                 V.color[totest[which.sup]] <- cols[i]
556                 V.color[fr2] <- cols[i]
557
558         tocomp <- get.edge.attribute(ng, tw)
559         totest <- intersect(which(!is.na(E.weight)), which(!is.na(tocomp)))
560         maxmat <- cbind(E.weight[totest], tocomp[totest])
561         resmax <- apply(maxmat, 1, which.max)
562         E.weight[totest] <- apply(maxmat, 1, max)
563         nas <- which(is.na(E.weight))
564         nas2 <- which(is.na(tocomp))
565         fr2 <- setdiff(nas,nas2)
566         E.weight[fr2] <- tocomp[fr2]            
567         }
568     V(ng)$weight <- V.weight
569     V(ng)$proto.color <- V.proto.color
570         V(ng)$color <- V.proto.color
571     E(ng)$weight <- E.weight
572         V(ng)$ocolor <- V.color
573         colors <- col2rgb(V(ng)$color)
574         V(ng)$r <- colors["red", ]
575         V(ng)$g <- colors["green", ]
576         V(ng)$b <- colors["blue", ]
577     ng
578 }
579
580
581 spirale <- function(g, weigth, center, miny=0.1) {
582         ncoord <- matrix(0, nrow=length(weigth)+1, ncol=2)
583         v.names <- V(g)$name 
584         center.name <- v.names[center]
585         first <- which.max(weigth)[1]
586         if (head_of(g, first)$name == center.name) {
587                 n.name <- tail_of(g, first)
588         } else {
589                 n.name <- head_of(g, first)
590         }
591         n.name <- n.name$name
592         nb <- length(weigth)
593         ncoord[which(v.names==n.name),] <- c(0,1)
594         weigth[first] <- 0
595         rs <- norm.vec(weigth,1, miny)
596         nbt <- nb %/% 50
597         if (nbt == 0) nbt <- 1
598         angler <- ((360 * nbt) / (nb- 1)) * (pi/180)
599         ang <- 90 * (pi/180)
600         rr <- (1-miny) / (nb-1)
601         r <- 1
602         while (max(weigth != 0)) {
603                 first <- which.max(weigth)[1]
604                 if (head_of(g, first)$name == center.name) {
605                         n.name <- tail_of(g, first)
606                 } else {
607                         n.name <- head_of(g, first)
608                 }
609                 n.name <- n.name$name
610                 #r <- rs[first]
611                 r <- r - rr
612                 ang <- ang + angler 
613                 x <- r * cos(ang)
614                 y <- r * sin(ang)
615                 weigth[first] <- 0
616                 ncoord[which(v.names==n.name),] <- c(x,y)
617         }
618         ncoord
619 }
620
621 spirale3D <- function(g, weigth, center, miny=0.1) {
622         ncoord <- matrix(0, nrow=length(weigth)+1, ncol=3)
623         v.names <- V(g)$name 
624         center.name <- v.names[center]
625         first <- which.max(weigth)[1]
626         if (head_of(g, first)$name == center.name) {
627                 n.name <- tail_of(g, first)
628         } else {
629                 n.name <- head_of(g, first)
630         }
631         n.name <- n.name$name
632         nb <- length(weigth)
633         ncoord[which(v.names==n.name),] <- c(0,0,1)
634         weigth[first] <- 0
635         rs <- norm.vec(weigth,1, miny)
636         nbt <- nb %/% 50
637         if (nbt == 0) nbt <- 1
638         angler <- ((360 * nbt) / (nb- 1)) * (pi/180)
639         theta <- 0
640         phi <- 90 * (pi/180)
641         rr <- (1-miny) / (nb-1)
642         r <- 1
643         while (max(weigth != 0)) {
644                 first <- which.max(weigth)[1]
645                 if (head_of(g, first)$name == center.name) {
646                         n.name <- tail_of(g, first)
647                 } else {
648                         n.name <- head_of(g, first)
649                 }
650                 n.name <- n.name$name
651                 #r <- rs[first]
652                 r <- r - rr
653                 theta <- theta + angler
654                 phi <- phi + angler/2
655                 x <- r * sin(theta) * cos(phi)
656                 y <- r * sin(theta) * sin(phi)
657                 z <- r * cos(theta)
658                 weigth[first] <- 0
659                 ncoord[which(v.names==n.name),] <- c(x,y,z)
660         }
661         ncoord
662 }
663