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)
6 make.a <- function(x) {
11 make.b <- function(x) {
16 make.c <- function(x) {
21 make.d <- function(x, a, b, c) {
23 d <- ncol(x) - a - b - c
27 ###########################################
29 ###########################################
30 my.jaccard <- function(x) {
34 #d <- make.d(x, a, b, c)
35 jac <- a / (a + b + c)
40 prcooc <- function(x, a) {
45 make.bin <- function(cs, a, i, j, nb) {
48 res <- binom.test(ab, nb, (cs[i]/nb) * (cs[j]/nb), "less")
53 #res <- binom.test(ab, nb, (cs[i]/nb) * (cs[j]/nb), "less")
57 binom.sim <- function(x) {
61 mat <- matrix(0,ncol(x),ncol(x))
62 colnames(mat)<-colnames(a)
63 rownames(mat)<-rownames(a)
64 for (i in 1:(ncol(x)-1)) {
65 for (j in (i+1):ncol(x)) {
66 mat[j,i] <- make.bin(cs, a, i, j , n)
74 ############################################
76 ############################################
77 # jaccard a, b, c a / (a + b + c)
78 # Kulczynski1 a, b, c a / (b + c)
79 # Kulczynski2 a, b, c [a / (a + b) + a / (a + c)] / 2
80 # Mountford a, b, c 2a / (ab + ac + 2bc)
81 # Fager, McGowan a, b, c a / sqrt((a + b)(a + c)) - 1 / 2 sqrt(a + c)
83 # Dice, Czekanowski, Sorensen a, b, c 2a / (2a + b + c)
84 # Mozley, Margalef a, b, c an / (a + b)(a + c)
85 # Ochiai a, b, c a / sqrt[(a + b)(a + c)]
86 # Simpson a, b, c a / min{(a + b), (a + c)}
87 # Braun-Blanquet a, b, c a / max{(a + b), (a + c)}
89 #simple matching, Sokal/Michener a, b, c, d, ((a + d) /n)
90 # Hamman, a, b, c, d, ([a + d] - [b + c]) / n
91 # Faith , a, b, c, d, (a + d/2) / n
92 # Tanimoto, Rogers a, b, c, d, (a + d) / (a + 2b + 2c + d)
93 # Phi a, b, c, d (ad - bc) / sqrt[(a + b)(c + d)(a + c)(b + d)]
94 # Stiles a, b, c, d log(n(|ad-bc| - 0.5n)^2 / [(a + b)(c + d)(a + c)(b + d)])
95 # Michael a, b, c, d 4(ad - bc) / [(a + d)^2 + (b + c)^2]
96 # Yule a, b, c, d (ad - bc) / (ad + bc)
97 # Yule2 a, b, c, d (sqrt(ad) - sqrt(bc)) / (sqrt(ad) + sqrt(bc))
99 BuildProf01<-function(x,classes) {
101 #classes : classes de chaque lignes de x
102 dm<-cbind(x,cl=classes)
103 clnb=length(summary(as.data.frame(as.character(classes)),max=100))
105 print(summary(as.data.frame(as.character(classes)),max=100))
106 mat<-matrix(0,ncol(x),clnb)
107 rownames(mat)<-colnames(x)
109 dtmp<-dm[which(dm$cl==i),]
110 for (j in 1:(ncol(dtmp)-1)) {
111 mat[j,i]<-sum(dtmp[,j])
117 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) {
120 v.label <- colnames(mat.simi)
121 g1<-graph.adjacency(mat.simi,mode="lower",weighted=TRUE)
123 weori<-get.edge.attribute(g1,'weight')
125 if (method == 'cooc') {
131 g.max<-minimum.spanning.tree(g1)
132 if (method == 'cooc') {
133 E(g.max)$weight<-1 / E(g.max)$weight
135 E(g.max)$weight<-1 - E(g.max)$weight
140 if (!is.null(seuil)) {
141 if (seuil >= max(mat.simi)) seuil <- -Inf
143 w<-E(g.toplot)$weight
144 tovire <- which(w<=seuil)
145 g.toplot <- delete.edges(g.toplot,(tovire))
146 for (i in 1:(length(V(g.toplot)))) {
147 if (length(neighbors(g.toplot,i))==0) {
151 g.toplot <- delete.vertices(g.toplot,vec)
152 v.label <- V(g.toplot)$name
153 if (!is.logical(vec)) mat.eff <- mat.eff[-(vec)]
158 if (!is.null(minmaxeff[1])) {
159 eff<-norm.vec(mat.eff,minmaxeff[1],minmaxeff[2])
163 if (!is.null(vcexminmax[1])) {
164 label.cex = norm.vec(mat.eff, vcexminmax[1], vcexminmax[2])
168 if (!is.null(coeff.edge)) {
170 we.width <- norm.vec(abs(E(g.toplot)$weight), coeff.edge[1], coeff.edge[2])
171 #we.width <- abs((E(g.toplot)$weight/max(abs(E(g.toplot)$weight)))*coeff.edge)
175 if (method != 'binom') {
176 we.label <- round(E(g.toplot)$weight,3)
178 we.label <- round(E(g.toplot)$weight,4)
180 if (p.type=='rgl' || p.type=='rglweb') {
185 if (! is.null(fromcoords)) {
186 newfrom <- matrix(runif(nd*length(V(g.toplot)$name),min(fromcoords)),max(fromcoords),ncol=nd, nrow=length(V(g.toplot)$name))
187 for (i in 1:length(V(g.toplot)$name)) {
188 if(V(g.toplot)$name[i] %in% forvertex) {
189 newfrom[i,] <- fromcoords[which(forvertex==V(g.toplot)$name[i]),]
192 fromcoords <- newfrom
195 if (is.null(coords)) {
196 if (layout.type == 'frutch') {
197 #lo <- layout_with_drl(g.toplot,dim=nd)
198 #lo <- layout_with_fr(g.toplot,dim=nd, grid="grid", niter=10000, weights=1/E(g.toplot)$weight)#, start.temp = 1)#, )
202 lo <- gplot.layout.fruchtermanreingold(asNetwork(g.toplot), list())
203 detach("package:intergraph", unload=TRUE)
204 detach("package:sna", unload=TRUE)
205 detach("package:network", unload=TRUE)
208 lo <- layout_with_fr(g.toplot,dim=nd)
211 if (layout.type == 'kawa') {
212 lo <- layout_with_kk(g.toplot,dim=nd, weights=1/E(g.toplot)$weight, start=fromcoords, epsilon=0, maxiter = 10000)
215 if (layout.type == 'random')
216 lo <- layout_on_grid(g.toplot,dim=nd)
217 if (layout.type == 'circle' & p.type != 'rgl')
218 lo <- layout_in_circle(g.toplot)
219 if (layout.type == 'circle' & p.type == 'rgl')
220 lo <- layout_on_sphere(g.toplot)
221 if (layout.type == 'graphopt')
222 lo <- layout_as_tree(g.toplot, circular = TRUE)
223 if (layout.type == 'spirale')
224 lo <- spirale(g.toplot, E(g.toplot)$weight, index.word)
225 if (layout.type == 'spirale3D')
226 lo <- spirale3D(g.toplot, E(g.toplot)$weight, index.word)
230 if (!is.null(communities)) {
231 if (communities == 0 ){
232 com <- edge.betweenness.community(g.toplot)
233 } else if (communities == 1) {
234 com <- fastgreedy.community(g.toplot)
235 } else if (communities == 2) {
236 com <- label.propagation.community(g.toplot)
237 } else if (communities == 3) {
238 com <- leading.eigenvector.community(g.toplot)
239 } else if (communities == 4) {
240 com <- multilevel.community(g.toplot)
241 } else if (communities == 5) {
242 com <- optimal.community(g.toplot)
243 } else if (communities == 6) {
244 com <- spinglass.community(g.toplot)
245 } else if (communities == 7) {
246 com <- walktrap.community(g.toplot)
252 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)
255 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') {
256 mat.simi <- graph.simi$mat
257 g.toplot <- graph.simi$graph
258 if (is.null(vertex.size)) {
259 vertex.size <- graph.simi$eff
261 vertex.size <- vertex.size
263 we.width <- graph.simi$we.width
265 #v.label <- vire.nonascii(graph.simi$v.label)
266 v.label <- graph.simi$v.label
271 we.label <- graph.simi$we.label
275 lo <- graph.simi$layout
276 #rownames(lo) <- v.label
277 if (!is.null(vertex.label.cex)) {
278 label.cex<-vertex.label.cex
280 label.cex = graph.simi$label.cex
284 alphas <- norm.vec(label.cex, 0.5,1)
286 if (length(vertex.label.color) == 1) {
287 for (i in 1:length(alphas)) {
288 nvlc <- append(nvlc, adjustcolor(vertex.label.color, alpha=alphas[i]))
291 for (i in 1:length(alphas)) {
292 nvlc <- append(nvlc, adjustcolor(vertex.label.color[i], alpha=alphas[i]))
295 vertex.label.color <- nvlc
297 if (p.type=='nplot') {
298 #print('ATTENTION - PAS OPEN FILE')
299 open_file_graph(filename, width = width, height = height, svg = svg)
303 layout(matrix(c(1,2),1,2, byrow=TRUE),widths=c(3,lcm(7)))
307 if (is.null(graph.simi$com)) {
308 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)
310 if (graph.simi$halo) {
311 mark.groups <- communities(graph.simi$com)
315 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)
318 txt.layout <- layout.norm(lo, -1, 1, -1, 1, -1, 1)
319 #txt.layout <- txt.layout[order(label.cex),]
320 #vertex.label.color <- vertex.label.color[order(label.cex)]
321 #v.label <- v.label[order(label.cex)]
322 #label.cex <- label.cex[order(label.cex)]
323 text(txt.layout[,1], txt.layout[,2], v.label, cex=label.cex, col=vertex.label.color)
326 plot(0, axes = FALSE, pch = '')
327 legend(x = 'center' , leg$unetoile, fill = leg$gcol)
332 if (p.type=='tkplot') {
333 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)
334 coords = tkplot.getcoords(id)
335 ok <- try(coords <- tkplot.getcoords(id), TRUE)
336 while (is.matrix(ok)) {
337 ok <- try(coords <- tkplot.getcoords(id), TRUE)
344 if (p.type == 'rgl' || p.type == 'rglweb') {
348 lo <- layout.norm(lo, -10, 10, -10, 10, -10, 10)
350 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)
351 #los <- layout.norm(lo, -1, 1, -1, 1, -1, 1)
352 text3d(lo[,1], lo[,2], lo[,3], vire.nonascii(v.label), col = vertex.label.color, alpha = 1, cex = vertex.label.cex)
353 rgl.spheres(lo, col = vertex.col, radius = vertex.size/100, alpha = alpha)
354 #rgl.bg(color = c('white','black'))
356 if (!is.null(movie)) {
358 ReturnVal <- tkmessageBox(title="RGL 3 D",message="Cliquez pour commencer le film",icon="info",type="ok")
360 movie3d(spin3d(axis=c(0,1,0),rpm=6), movie = 'film_graph', frames = "tmpfilm", duration=10, clean=TRUE, top = TRUE, dir = movie)
361 ReturnVal <- tkmessageBox(title="RGL 3 D",message="Film fini !",icon="info",type="ok")
363 #play3d(spin3d(axis=c(0,1,0),rpm=6))
364 if (p.type == 'rglweb') {
365 writeWebGL(dir = filename, width = width, height= height)
368 ReturnVal <- tkmessageBox(title="RGL 3 D",message="Cliquez pour fermer",icon="info",type="ok")
371 # while (rgl.cur() != 0)
373 } else if (p.type == 'web') {
375 graph.simi$label.cex <- label.cex
376 if (length(vertex.col)==1) {
377 vertex.col <- rep(vertex.col, length(v.label))
379 graph.simi$color <- vertex.col
381 nodes.attr <- data.frame(label)
382 simi.to.gexf(filename, graph.simi, nodes.attr = nodes.attr)
387 graph.word <- function(mat.simi, index) {
388 nm <- matrix(0, ncol = ncol(mat.simi), nrow=nrow(mat.simi), dimnames=list(row.names(mat.simi), colnames(mat.simi)))
389 nm[,index] <- mat.simi[,index]
390 nm[index,] <- mat.simi[index,]
395 #http://gopalakrishna.palem.in/iGraphExport.html#GexfExport
396 # Converts the given igraph object to GEXF format and saves it at the given filepath location
397 # g: input igraph object to be converted to gexf format
398 # filepath: file location where the output gexf file should be saved
400 saveAsGEXF = function(g, filepath="converted_graph.gexf")
405 # gexf nodes require two column data frame (id, label)
406 # check if the input vertices has label already present
407 # if not, just have the ids themselves as the label
408 if(is.null(V(g)$label))
409 V(g)$label <- as.character(V(g))
411 # similarily if edges does not have weight, add default 1 weight
412 if(is.null(E(g)$weight))
413 E(g)$weight <- rep.int(1, ecount(g))
415 nodes <- data.frame(cbind(1:vcount(g), V(g)$label))
416 nodes[,1] <- as.character(nodes[,1])
417 nodes[,2] <- as.character(nodes[,2])
418 edges <- t(Vectorize(get.edge, vectorize.args='id')(g, 1:ecount(g)))
420 # combine all node attributes into a matrix (and take care of & for xml)
421 vAttrNames <- setdiff(list.vertex.attributes(g), "label")
422 for (val in c("x","y","color")) {
423 vAttrNames <- setdiff(vAttrNames, val)
425 nodesAtt <- data.frame(sapply(vAttrNames, function(attr) sub("&", "&",get.vertex.attribute(g, attr))))
426 for (i in 1:ncol(nodesAtt)) {
427 nodesAtt[,i] <- as.character(nodesAtt[,i])
430 # combine all edge attributes into a matrix (and take care of & for xml)
431 eAttrNames <- setdiff(list.edge.attributes(g), "weight")
432 edgesAtt <- data.frame(sapply(eAttrNames, function(attr) sub("&", "&",get.edge.attribute(g, attr))))
434 # combine all graph attributes into a meta-data
435 graphAtt <- sapply(list.graph.attributes(g), function(attr) sub("&", "&",get.graph.attribute(g, attr)))
437 cc <- t(sapply(V(g)$color, col2rgb, alpha=TRUE))
439 # generate the gexf object
440 output <- write.gexf(nodes, edges,
441 edgesWeight=E(g)$weight,
443 #edgesVizAtt = list(size=as.matrix(E(g)$weight)),
445 nodesVizAtt=list(color=cc, position=cbind(V(g)$x,V(g)$y, rep(0,ll)), size=V(g)$weight),
446 meta=c(list(creator="iramuteq", description="igraph -> gexf converted file", keywords="igraph, gexf, R, rgexf"), graphAtt))
448 print(output, filepath, replace=T)
452 merge.graph <- function(graphs) {
454 ng <- graph.union(graphs, byname=T)
455 V.weight <- V(ng)$weight_1
456 E.weight <- E(ng)$weight_1
457 cols <- rainbow(length(graphs))
459 V.color <- rep(cols[1], length(V.weight))
460 for (i in 2:length(graphs)) {
461 tw <- paste('weight_', i, sep='')
462 tocomp <- get.vertex.attribute(ng,tw)
463 totest <- intersect(which(!is.na(V.weight)), which(!is.na(tocomp)))
464 maxmat <- cbind(V.weight[totest], tocomp[totest])
465 resmax <- apply(maxmat, 1, which.max)
466 ncolor <- c(cols[(i-1)], cols[i])
467 #rbgcol1 <- col2rgb(cols[(i-1)])
468 #rbgcol1 <- rbgcol1/255
469 #rgbcol1 <- RGB(rbgcol1[1],rbgcol1[2],rbgcol1[3])
470 rbgcol2 <- col2rgb(cols[i])
471 rbgcol2 <- rbgcol2/255
472 #rgbcol2 <- RGB(rbgcol2[1],rbgcol2[2],rbgcol2[3])
474 alpha <- tocomp[j] /(V.weight[j] + tocomp[j])
475 rbgcol1 <- col2rgb(V.color[j])
476 rbgcol1 <- rbgcol1/255
477 #mix.col <- mixcolor(alpha,rbgcol1, rbgcol2)
478 mix.col <- mixcolor(alpha, RGB(rbgcol1[1],rbgcol1[2],rbgcol1[3]), RGB(rbgcol2[1],rbgcol2[2],rbgcol2[3]))
479 V.color[j] <- hex(mix.col)
480 #V.color[j] <- adjustcolor(hex(mix.col), 0.6)
482 #to.change <- totest[which(resmax == 2)]
483 #V.color[to.change] <- cols[i]
484 V.weight[totest] <- apply(maxmat, 1, max)
485 nas <- which(is.na(V.weight))
486 nas2 <- which(is.na(tocomp))
487 fr2 <- setdiff(nas,nas2)
488 V.weight[fr2] <- tocomp[fr2]
489 V.color[fr2] <- cols[i]
490 tocomp <- get.edge.attribute(ng, tw)
491 totest <- intersect(which(!is.na(E.weight)), which(!is.na(tocomp)))
492 maxmat <- cbind(E.weight[totest], tocomp[totest])
493 resmax <- apply(maxmat, 1, which.max)
494 E.weight[totest] <- apply(maxmat, 1, max)
495 nas <- which(is.na(E.weight))
496 nas2 <- which(is.na(tocomp))
497 fr2 <- setdiff(nas,nas2)
498 E.weight[fr2] <- tocomp[fr2]
500 V(ng)$weight <- V.weight
501 V(ng)$color <- V.color
502 E(ng)$weight <- E.weight
503 colors <- col2rgb(V(ng)$color)
504 V(ng)$r <- colors["red", ]
505 V(ng)$g <- colors["green", ]
506 V(ng)$b <- colors["blue", ]
510 spirale <- function(g, weigth, center, miny=0.1) {
511 ncoord <- matrix(0, nrow=length(weigth)+1, ncol=2)
513 center.name <- v.names[center]
514 first <- which.max(weigth)[1]
515 if (head_of(g, first)$name == center.name) {
516 n.name <- tail_of(g, first)
518 n.name <- head_of(g, first)
520 n.name <- n.name$name
522 ncoord[which(v.names==n.name),] <- c(0,1)
524 rs <- norm.vec(weigth,1, miny)
526 if (nbt == 0) nbt <- 1
527 angler <- ((360 * nbt) / (nb- 1)) * (pi/180)
529 rr <- (1-miny) / (nb-1)
531 while (max(weigth != 0)) {
532 first <- which.max(weigth)[1]
533 if (head_of(g, first)$name == center.name) {
534 n.name <- tail_of(g, first)
536 n.name <- head_of(g, first)
538 n.name <- n.name$name
545 ncoord[which(v.names==n.name),] <- c(x,y)
550 spirale3D <- function(g, weigth, center, miny=0.1) {
551 ncoord <- matrix(0, nrow=length(weigth)+1, ncol=3)
553 center.name <- v.names[center]
554 first <- which.max(weigth)[1]
555 if (head_of(g, first)$name == center.name) {
556 n.name <- tail_of(g, first)
558 n.name <- head_of(g, first)
560 n.name <- n.name$name
562 ncoord[which(v.names==n.name),] <- c(0,0,1)
564 rs <- norm.vec(weigth,1, miny)
566 if (nbt == 0) nbt <- 1
567 angler <- ((360 * nbt) / (nb- 1)) * (pi/180)
570 rr <- (1-miny) / (nb-1)
572 while (max(weigth != 0)) {
573 first <- which.max(weigth)[1]
574 if (head_of(g, first)$name == center.name) {
575 n.name <- tail_of(g, first)
577 n.name <- head_of(g, first)
579 n.name <- n.name$name
582 theta <- theta + angler
583 phi <- phi + angler/2
584 x <- r * sin(theta) * cos(phi)
585 y <- r * sin(theta) * sin(phi)
588 ncoord[which(v.names==n.name),] <- c(x,y,z)