...
[iramuteq] / Rscripts / Rgraph.R
1 ############FIXME##################
2 #PlotDendroComp <- function(chd,filename,reso) {
3 #   jpeg(filename,res=reso)
4 #   par(cex=PARCEX)
5 #   plot(chd,which.plots=2, hang=-1)
6 #   dev.off()
7 #}
8 #
9 #PlotDendroHori <- function(dendrocutupper,filename,reso) {
10 #   jpeg(filename,res=reso)
11 #   par(cex=PARCEX)
12 #   nP <- list(col=3:2, cex=c(0.5, 0.75), pch= 21:22, bg= c('light blue', 'pink'),lab.cex = 0.75, lab.col = 'tomato')
13 #   plot(dendrocutupper,nodePar= nP, edgePar = list(col='gray', lwd=2),horiz=TRUE, center=FALSE)
14 #   dev.off()
15 #}
16
17 PlotDendroCut <- function(chd,filename,reso,clusternb) {   
18    h.chd <- as.hclust(chd)
19    memb <- cutree(h.chd, k = clusternb)
20    cent <- NULL
21    for(k in 1:clusternb){
22        cent <- rbind(cent, k)
23    }
24    h.chd1 <- hclust(dist(cent)^2, method = 'cen', members = table(memb))
25    h.chd1$labels <- sprintf('CL %02d',1:clusternb)
26    nP <- list(col=3:2, cex=c(2.0, 0.75), pch= 21:22, bg= c('light blue', 'pink'),lab.cex = 0.75, lab.col = 'tomato')
27    jpeg(filename,res=reso)
28    par(cex=PARCEX)
29    plot(h.chd1, nodePar= nP, edgePar = list(col='gray', lwd=2), horiz=TRUE, center=TRUE, hang= -1)
30    dev.off()
31 }
32
33 #PlotAfc<- function(afc, filename, width=800, height=800, quality=100, reso=200, toplot=c('all','all'), PARCEX=PARCEX) {
34 #       if (Sys.info()["sysname"]=='Darwin') {
35 #               width<-width/74.97
36 #               height<-height/74.97
37 #               quartz(file=filename,type='jpeg',width=width,height=height)
38 #       } else {
39 #               jpeg(filename,width=width,height=height,quality=quality,res=reso)
40 #       }
41 #       par(cex=PARCEX)
42 #       plot(afc,what=toplot,labels=c(1,1),contrib=c('absolute','relative'))
43 #       dev.off()
44 #}
45
46 PlotAfc2dCoul<- function(afc,chisqrtable,filename, what='coord',col=FALSE, axetoplot=c(1,2), deb=0,fin=0, width=900, height=900, quality=100, reso=200, parcex=PARCEX, xlab = NULL, ylab = NULL, xmin=NULL, xmax=NULL, ymin=NULL, ymax=NULL, active = TRUE) {
47         if (col) {
48                 if (what == 'coord') {
49                         rowcoord <- as.matrix(afc$colcoord)
50                 } else {
51                         rowcoord <- as.matrix(afc$colcrl)
52                 }
53         } else {
54                 if (what == 'coord') {
55                         rowcoord <- as.matrix(afc$rowcoord)
56                 } else {
57                         rowcoord <- as.matrix(afc$rowcrl)
58                 }
59         }
60         x <- axetoplot[1]
61         y <- axetoplot[2]
62         if (col)
63                 rownames(rowcoord) <- afc$colnames
64         if (!col){
65                 rownames(rowcoord) <- afc$rownames
66                 rowcoord <- as.matrix(rowcoord[deb:fin,])
67                 chitable<- as.matrix(chisqrtable[deb:fin,])
68                 #row_keep <- select_point_nb(chitable,15)
69         }
70         if (ncol(rowcoord) == 1) {
71                 rowcoord <- t(rowcoord)
72         }
73         clnb <- ncol(chisqrtable)
74         
75         if (!col) {
76         classes <- as.matrix(apply(chitable,1,which.max))
77         cex.par <- norm.vec(apply(chitable,1,max), 0.8,3)
78         row.keep <- select.chi.classe(chitable, 80, active=active)
79         rowcoord <- rowcoord[row.keep,]
80         classes <- classes[row.keep]
81         cex.par <- cex.par[row.keep]
82         } else {
83         classes <- 1:clnb
84         cex.par <- rep(1,clnb)
85     }
86     if (is.null(xmin)) {
87         table.in <- rowcoord
88         xminmax <- c(min(table.in[,1], na.rm = TRUE) + ((max(cex.par)/10) * min(table.in[,1], na.rm = TRUE)), max(table.in[,1], na.rm = TRUE) + ((max(cex.par)/10) * max(table.in[,1], na.rm = TRUE)))
89         xmin <- xminmax[1]
90         xmax <- xminmax[2]
91         yminmax <- c(min(table.in[,2], na.rm = TRUE) + ((max(cex.par)/10) * min(table.in[,2], na.rm = TRUE)), max(table.in[,2], na.rm = TRUE) + ((max(cex.par)/10) * max(table.in[,2], na.rm = TRUE)))
92         ymin <- yminmax[1]
93         ymax <- yminmax[2]
94      }
95         #ntabtot <- cbind(rowcoord, classes)
96         #if (!col) ntabtot <- ntabtot[row_keep,]
97     xlab <- paste('facteur ', x, ' -')
98     ylab <- paste('facteur ', y, ' -')
99     xlab <- paste(xlab,round(afc_table$facteur[x,2],2),sep = ' ')
100     xlab <- paste(xlab,' %%',sep = '')
101     ylab <- paste(ylab,round(afc_table$facteur[y,2],2),sep = ' ')
102     ylab <- paste(ylab,' %%',sep = '')
103
104         open_file_graph(filename, width = width, height = height)
105         par(cex=PARCEX)
106     table.in <- rowcoord[order(cex.par, decreasing = TRUE),]
107     classes <- classes[order(cex.par, decreasing = TRUE)]
108     cex.par <- cex.par[order(cex.par, decreasing = TRUE)]
109     table.out <- stopoverlap(table.in, cex.par=cex.par, xlim = c(xmin,xmax), ylim = c(ymin,ymax))
110     table.in <- table.out$toplot
111     notplot <- table.out$notplot
112     if (! is.null(notplot)) {
113         write.csv2(notplot, file = paste(filename, '_notplotted.csv', sep = ''))
114     }
115     classes <- classes[table.in[,4]]
116     cex.par <- cex.par[table.in[,4]]
117     make_afc_graph(table.in, classes, clnb, xlab, ylab, cex.txt = cex.par, xminmax=c(xmin,xmax), yminmax=c(ymin,ymax))
118     xyminmax <- list(yminmax = c(ymin,ymax), xminmax = c(xmin,xmax))
119     xyminmax 
120         #plot(rowcoord[,x],rowcoord[,y], pch='', xlab = xlab, ylab = ylab)
121         #abline(h=0,v=0)
122         #for (i in 1:clnb) {
123         #       ntab <- subset(ntabtot,ntabtot[,ncol(ntabtot)] == i)
124         #       if (nrow(ntab) != 0)
125         #               text(ntab[,x],ntab[,y],rownames(ntab),col=rainbow(clnb)[i])
126         #}
127         #dev.off()
128 }
129
130 filename.to.svg <- function(filename) {
131     filename <- gsub('.png', '.svg', filename)
132     return(filename)
133 }
134
135 open_file_graph <- function (filename, width=800, height = 800, quality = 100, svg = FALSE) {
136         if (Sys.info()["sysname"] == 'Darwin') {
137         width <- width/74.97
138         height <- height/74.97
139         if (!svg) {
140                     quartz(file = filename, type = 'png', width = width, height = height)
141         } else {
142             svg(filename.to.svg(filename), width=width, height=height)
143         }
144         } else {
145         if (svg) {
146             svg(filename.to.svg(filename), width=width/74.97, height=height/74.97)
147         } else {
148                     png(filename, width=width, height=height)#, quality = quality)
149         }
150         }
151 }
152
153 #################################################@@
154 #from wordcloud
155 overlap <- function(x1, y1, sw1, sh1, boxes) {
156     use.r.layout <- FALSE
157         if(!use.r.layout)
158                 return(.overlap(x1,y1,sw1,sh1,boxes))
159         s <- 0
160         if (length(boxes) == 0) 
161                 return(FALSE)
162         for (i in c(last,1:length(boxes))) {
163                 bnds <- boxes[[i]]
164                 x2 <- bnds[1]
165                 y2 <- bnds[2]
166                 sw2 <- bnds[3]
167                 sh2 <- bnds[4]
168                 if (x1 < x2) 
169                         overlap <- x1 + sw1 > x2-s
170                 else 
171                         overlap <- x2 + sw2 > x1-s
172                 
173                 if (y1 < y2) 
174                         overlap <- overlap && (y1 + sh1 > y2-s)
175                 else 
176                         overlap <- overlap && (y2 + sh2 > y1-s)
177                 if(overlap){
178                         last <<- i
179                         return(TRUE)
180                 }
181         }
182         FALSE
183 }
184
185 .overlap <- function(x11,y11,sw11,sh11,boxes1){
186     .Call("is_overlap",x11,y11,sw11,sh11,boxes1)
187 }
188 ########################################################
189 stopoverlap <- function(x, cex.par = NULL, xlim = NULL, ylim = NULL) {
190 #from wordcloud
191     library(wordcloud)
192     tails <- "g|j|p|q|y"
193     rot.per <- 0 
194     last <- 1
195     thetaStep <- .1
196     rStep <- .5
197     toplot <- NULL
198     notplot <- NULL
199
200 #    plot.new()
201     plot(x[,1],x[,2], pch='', xlim = xlim, ylim = ylim)
202
203     words <- rownames(x)
204     if  (is.null(cex.par))  {
205         size <- rep(0.9, nrow(x))
206     } else {
207         size <- cex.par
208     }
209     #cols <- rainbow(clnb)
210     boxes <- list()
211     for (i in 1:nrow(x)) {
212         rotWord <- runif(1)<rot.per
213         r <-0
214                 theta <- runif(1,0,2*pi)
215                 x1<- x[i,1] 
216                 y1<- x[i,2]
217                 wid <- strwidth(words[i],cex=size[i])
218                 ht <- strheight(words[i],cex=size[i])
219                 isOverlaped <- TRUE
220                 while(isOverlaped){
221                         if(!overlap(x1-.5*wid,y1-.5*ht,wid,ht, boxes)) { #&&
222                 toplot <- rbind(toplot, c(x1, y1, size[i], i)) 
223                                 #text(x1,y1,words[i],cex=size[i],offset=0,srt=rotWord*90,
224                                 #               col=cols[classes[i]])
225                                 boxes[[length(boxes)+1]] <- c(x1-.5*wid,y1-.5*ht,wid,ht)
226                                 isOverlaped <- FALSE
227                         } else {
228                                 if(r>sqrt(.5)){
229                                         print(paste(words[i], "could not be fit on page. It will not be plotted."))
230                     notplot <- rbind(notplot,c(words[i], x[i,1], x[i,2]))
231                                         isOverlaped <- FALSE
232                                 }
233                                 theta <- theta+thetaStep
234                                 r <- r + rStep*thetaStep/(2*pi)
235                 x1 <- x[i,1]+r*cos(theta)
236                                 y1 <- x[i,2]+r*sin(theta)
237                         }
238                 }
239     }
240     row.names(toplot) <- words[toplot[,4]]
241     return(list(toplot = toplot, notplot = notplot))
242 }
243 ###############################################################################
244
245 getwordcloudcoord <- function(words,freq,scale=c(4,.5),min.freq=3,max.words=Inf,random.order=TRUE,random.color=FALSE,
246                 rot.per=.1,colors="black",ordered.colors=FALSE,use.r.layout=FALSE,fixed.asp=TRUE,...) { 
247         tails <- "g|j|p|q|y"
248         last <- 1
249         
250         overlap <- function(x1, y1, sw1, sh1) {
251                 if(!use.r.layout)
252                         return(.overlap(x1,y1,sw1,sh1,boxes))
253                 s <- 0
254                 if (length(boxes) == 0) 
255                         return(FALSE)
256                 for (i in c(last,1:length(boxes))) {
257                         bnds <- boxes[[i]]
258                         x2 <- bnds[1]
259                         y2 <- bnds[2]
260                         sw2 <- bnds[3]
261                         sh2 <- bnds[4]
262                         if (x1 < x2) 
263                                 overlap <- x1 + sw1 > x2-s
264                         else 
265                                 overlap <- x2 + sw2 > x1-s
266                         
267                         if (y1 < y2) 
268                                 overlap <- overlap && (y1 + sh1 > y2-s)
269                         else 
270                                 overlap <- overlap && (y2 + sh2 > y1-s)
271                         if(overlap){
272                                 last <<- i
273                                 return(TRUE)
274                         }
275                 }
276                 FALSE
277         }
278         
279         ord <- rank(-freq, ties.method = "random")
280         words <- words[ord<=max.words]
281         freq <- freq[ord<=max.words]
282
283
284         ord <- order(freq,decreasing=TRUE)
285         words <- words[ord]
286         freq <- freq[ord]
287         words <- words[freq>=min.freq]
288         freq <- freq[freq>=min.freq]
289         if (ordered.colors) {
290                 colors <- colors[ord][freq>=min.freq]
291         }
292         
293         thetaStep <- .1
294         rStep <- .05
295         plot.new()
296
297         normedFreq <- freq/max(freq)
298         size <- (scale[1]-scale[2])*normedFreq + scale[2]
299         boxes <- list()
300         toplot <- NULL  
301         
302         
303         for(i in 1:length(words)){
304                 rotWord <- runif(1)<rot.per
305                 r <-0
306                 theta <- runif(1,0,2*pi)
307                 x1<-.5
308                 y1<-.5
309                 wid <- strwidth(words[i],cex=size[i],...)
310                 ht <- strheight(words[i],cex=size[i],...)
311                 #mind your ps and qs
312                 if(grepl(tails,words[i]))
313                         ht <- ht + ht*.2
314                 if(rotWord){
315                         tmp <- ht
316                         ht <- wid
317                         wid <- tmp      
318                 }
319                 isOverlaped <- TRUE
320                 while(isOverlaped){
321                         if(!overlap(x1-.5*wid,y1-.5*ht,wid,ht) &&
322                                         x1-.5*wid>0 && y1-.5*ht>0 &&
323                                         x1+.5*wid<1 && y1+.5*ht<1){
324                                 toplot <- rbind(toplot, c(x1,y1,size[i], i))
325                                 boxes[[length(boxes)+1]] <- c(x1-.5*wid,y1-.5*ht,wid,ht)
326                                 isOverlaped <- FALSE
327                         }else{
328                                 if(r>sqrt(.5)){
329                                         warning(paste(words[i],
330                                                                         "could not be fit on page. It will not be plotted."))
331                                         isOverlaped <- FALSE
332                                 }
333                                 theta <- theta+thetaStep
334                                 r <- r + rStep*thetaStep/(2*pi)
335                                 x1 <- .5+r*cos(theta)
336                                 y1 <- .5+r*sin(theta)
337                         }
338                 }
339         }
340         toplot <- cbind(toplot,norm.vec(freq[toplot[,4]], 1, 50))
341         row.names(toplot) <- words[toplot[,4]]
342         toplot <- toplot[,-4]
343         return(toplot)
344 }
345
346 make_tree_tot <- function (chd) {
347         library(ape)
348         lf<-chd$list_fille
349         clus<-'a1a;'
350         for (i in 1:length(lf)) {
351                 if (!is.null(lf[[i]])) {
352                         clus<-gsub(paste('a',i,'a',sep=''),paste('(','a',lf[[i]][1],'a',',','a',lf[[i]][2],'a',')',sep=''),clus)
353         }
354         }
355         dendro_tuple <- clus
356         clus <- gsub('a','',clus)
357         tree.cl <- read.tree(text = clus)
358         res<-list(tree.cl = tree.cl, dendro_tuple = dendro_tuple)
359         res
360 }
361
362 make_dendro_cut_tuple <- function(dendro_in, coordok, classeuce, x, nbt = 9) {
363         library(ape)
364         dendro<-dendro_in
365         i <- 0
366         for (cl in coordok[,x]) {
367                 i <- i + 1
368                 fcl<-fille(cl,classeuce)
369                 for (fi in fcl) {
370                         dendro <- gsub(paste('a',fi,'a',sep=''),paste('b',i,'b',sep=''),dendro)
371                 }
372         }
373         clnb <- nrow(coordok)
374     tcl=((nbt+1) *2) - 2
375         for (i in 1:(tcl + 1)) {
376                 dendro <- gsub(paste('a',i,'a',sep=''),paste('b',0,'b',sep=''),dendro)
377         }
378         dendro <- gsub('b','',dendro)
379         dendro <- gsub('a','',dendro)
380         dendro_tot_cl <- read.tree(text = dendro)
381         #FIXME
382         for (i in 1:100) {
383                 for (cl in 1:clnb) {
384                         dendro <- gsub(paste('\\(',cl,',',cl,'\\)',sep=''),cl,dendro)
385                 }
386         }
387         for (i in 1:100) {
388                 dendro <- gsub(paste('\\(',0,',',0,'\\)',sep=''),0,dendro)
389                 for (cl in 1:clnb) {
390                         dendro <- gsub(paste('\\(',0,',',cl,'\\)',sep=''),cl,dendro)
391                         dendro <- gsub(paste('\\(',cl,',',0,'\\)',sep=''),cl,dendro)
392                 }
393         }
394         print(dendro)
395         tree.cl <- read.tree(text = dendro)
396     lab <- tree.cl$tip.label
397     if ("0" %in% lab) {
398         tovire <- which(lab == "0")
399         tree.cl <- drop.tip(tree.cl, tip = tovire)
400     }
401         res <- list(tree.cl = tree.cl, dendro_tuple_cut = dendro, dendro_tot_cl = dendro_tot_cl)
402         res
403 }
404
405 select_point_nb <- function(tablechi, nb) {
406         chimax<-as.matrix(apply(tablechi,1,max))
407         chimax<-cbind(chimax,1:nrow(tablechi))
408         order_chi<-as.matrix(chimax[order(chimax[,1],decreasing = TRUE),])
409         row_keep <- order_chi[,2][1:nb]
410         row_keep
411 }
412
413 select_point_chi <- function(tablechi, chi_limit) {
414         chimax<-as.matrix(apply(tablechi,1,max))
415         row_keep <- which(chimax >= chi_limit)
416         row_keep
417 }
418
419 select.chi.classe <- function(tablechi, nb, active = TRUE) {
420     rowkeep <- NULL
421     if (active & !is.null(debsup)) {
422         print(debsup)
423         print('###############################################################@')
424         tablechi <- tablechi[1:(debsup-1),]
425     }
426     if (nb > nrow(tablechi)) {
427         nb <- nrow(tablechi)
428     }
429     for (i in 1:ncol(tablechi)) {
430         rowkeep <- append(rowkeep,order(tablechi[,i], decreasing = TRUE)[1:nb])
431     }
432     rowkeep <- unique(rowkeep)
433     rowkeep
434 }
435
436 #from summary.ca
437 summary.ca.dm <- function(object, scree = TRUE, ...){
438   obj <- object
439   nd  <- obj$nd
440   if (is.na(nd)){
441     nd <- 2
442     } else {
443     if (nd > length(obj$sv)) nd <- length(obj$sv)
444     }  
445  # principal coordinates:
446   K   <- nd
447   I   <- dim(obj$rowcoord)[1] ; J <- dim(obj$colcoord)[1]
448   svF <- matrix(rep(obj$sv[1:K], I), I, K, byrow = TRUE)
449   svG <- matrix(rep(obj$sv[1:K], J), J, K, byrow = TRUE)
450   rpc <- obj$rowcoord[,1:K] * svF
451   cpc <- obj$colcoord[,1:K] * svG
452
453  # rows:
454   r.names <- obj$rownames
455   sr      <- obj$rowsup
456   if (!is.na(sr[1])) r.names[sr] <- paste("(*)", r.names[sr], sep = "")
457   r.mass <- obj$rowmass
458   r.inr  <- obj$rowinertia / sum(obj$rowinertia, na.rm = TRUE)
459   r.COR  <- matrix(NA, nrow = length(r.names), ncol = nd)
460   colnames(r.COR) <- paste('COR -facteur', 1:nd, sep=' ')
461   r.CTR  <- matrix(NA, nrow = length(r.names), ncol = nd)
462   colnames(r.CTR) <- paste('CTR -facteur', 1:nd, sep=' ')
463   for (i in 1:nd){
464     r.COR[,i] <- obj$rowmass * rpc[,i]^2 / obj$rowinertia
465     r.CTR[,i] <- obj$rowmass * rpc[,i]^2 / obj$sv[i]^2
466     }
467  # cor and quality for supplementary rows
468   if (length(obj$rowsup) > 0){
469     i0 <- obj$rowsup
470     for (i in 1:nd){
471       r.COR[i0,i] <- obj$rowmass[i0] * rpc[i0,i]^2
472       r.CTR[i0,i] <- NA
473     }
474     }
475
476  # columns:
477   c.names <- obj$colnames
478   sc      <- obj$colsup
479   if (!is.na(sc[1])) c.names[sc] <- paste("(*)", c.names[sc], sep = "")
480   c.mass  <- obj$colmass
481   c.inr   <- obj$colinertia / sum(obj$colinertia, na.rm = TRUE)
482   c.COR   <- matrix(NA, nrow = length(c.names), ncol = nd)
483   colnames(c.COR) <- paste('COR -facteur', 1:nd, sep=' ')
484   c.CTR   <- matrix(NA, nrow = length(c.names), ncol = nd)
485   colnames(c.CTR) <- paste('CTR -facteur', 1:nd, sep=' ')
486   for (i in 1:nd)
487     {
488     c.COR[,i] <- obj$colmass * cpc[,i]^2 / obj$colinertia
489     c.CTR[,i] <- obj$colmass * cpc[,i]^2 / obj$sv[i]^2
490     }
491   if (length(obj$colsup) > 0){
492     i0 <- obj$colsup
493     for (i in 1:nd){
494       c.COR[i0,i] <- obj$colmass[i0] * cpc[i0,i]^2
495       c.CTR[i0,i] <- NA
496       }
497     }
498
499  # scree plot:
500   if (scree) {
501     values     <- obj$sv^2
502     values2    <- 100*(obj$sv^2)/sum(obj$sv^2)
503     values3    <- cumsum(100*(obj$sv^2)/sum(obj$sv^2))
504     scree.out  <- cbind(values, values2, values3)
505     } else {
506     scree.out <- NA
507     }
508
509   obj$r.COR <- r.COR
510   obj$r.CTR <- r.CTR
511   obj$c.COR <- c.COR
512   obj$c.CTR <- c.CTR
513   obj$facteur <- scree.out
514   return(obj)
515   }
516
517 create_afc_table <- function(x) {
518    #x = afc
519         facteur.table <- as.matrix(x$facteur)
520     nd <- ncol(x$colcoord)
521         rownames(facteur.table) <- paste('facteur',1:nrow(facteur.table),sep = ' ')
522     colnames(facteur.table) <- c('Valeurs propres', 'Pourcentages', 'Pourcentage cumules')
523         ligne.table <- as.matrix(x$rowcoord)
524         rownames(ligne.table) <- x$rownames
525         colnames(ligne.table) <- paste('Coord. facteur', 1:nd, sep=' ')
526     tmp <- as.matrix(x$rowcrl)
527         colnames(tmp) <- paste('Corr. facteur', 1:nd, sep=' ')
528         ligne.table <- cbind(ligne.table,tmp)
529         ligne.table <- cbind(ligne.table, x$r.COR)
530         ligne.table <- cbind(ligne.table, x$r.CTR)
531         ligne.table <- cbind(ligne.table, mass = x$rowmass)
532         ligne.table <- cbind(ligne.table, chi.distance = x$rowdist)
533         ligne.table <- cbind(ligne.table, inertie = x$rowinertia)
534     colonne.table <- x$colcoord
535         rownames(colonne.table) <- paste('classe', 1:(nrow(colonne.table)),sep=' ')
536         colnames(colonne.table) <- paste('Coord. facteur', 1:nd, sep=' ')
537     tmp <- as.matrix(x$colcrl)
538         colnames(tmp) <- paste('Corr. facteur', 1:nd, sep=' ')
539         colonne.table <- cbind(colonne.table, tmp)
540         colonne.table <- cbind(colonne.table, x$c.COR)
541         colonne.table <- cbind(colonne.table, x$c.CTR)
542         colonne.table <- cbind(colonne.table, mass = x$colmass)
543         colonne.table <- cbind(colonne.table, chi.distance = x$coldist)
544         colonne.table <- cbind(colonne.table, inertie = x$colinertia)
545     res <- list(facteur = facteur.table, ligne = ligne.table, colonne = colonne.table)
546         res
547 }
548
549 is.yellow <- function(my.color) {
550     if ((my.color[1] > 200) & (my.color[2] > 200) & (my.color[3] < 20)) {
551         return(TRUE)
552     } else {
553         return(FALSE)
554     }
555 }
556
557 del.yellow <- function(colors) {
558     rgbs <- col2rgb(colors)
559     tochange <- apply(rgbs, 2, is.yellow)
560     tochange <- which(tochange)
561     if (length(tochange)) {
562         gr.col <- grey.colors(length(tochange), start = 0.5, end = 0.8)
563     }
564     compt <- 1
565     for (val in tochange) {
566         colors[val] <- gr.col[compt]
567         compt <- compt + 1
568     }
569     colors
570 }
571
572 make_afc_graph <- function(toplot, classes, clnb, xlab, ylab, cex.txt = NULL, leg = FALSE, cmd = FALSE, black = FALSE, xminmax=NULL, yminmax=NULL) {
573     
574     rain <- rainbow(clnb)
575     compt <- 1
576     tochange <- NULL
577     #for (my.color in rain) {
578     #    my.color <- col2rgb(my.color)
579     #    if ((my.color[1] > 200) & (my.color[2] > 200) & (my.color[3] < 20)) {
580     #       tochange <- append(tochange, compt)   
581     #    }
582     #    compt <- compt + 1
583     #}
584     #if (!is.null(tochange)) {
585     #    gr.col <- grey.colors(length(tochange))
586     #    compt <- 1
587     #    for (val in tochange) {
588     #        rain[val] <- gr.col[compt]
589     #        compt <- compt + 1
590     #    }
591     #}
592         rain <- del.yellow(rain)
593     cl.color <- rain[classes]
594     if (black) {
595         cl.color <- 'black'
596     }
597         plot(toplot[,1],toplot[,2], pch='', xlab = xlab, ylab = ylab, xlim=xminmax, ylim = yminmax)
598         abline(h=0, v=0, lty = 'dashed')
599         if (is.null(cex.txt))
600                 text(toplot[,1],toplot[,2],rownames(toplot),col=cl.color, offset=0)
601         else 
602                 #require(wordcloud)
603                 #textplot(toplot[,1],toplot[,2],rownames(toplot),col=cl.color, cex = cex.txt, xlim=xminmax, ylim = yminmax)
604         text(toplot[,1],toplot[,2],rownames(toplot),col=cl.color, cex = cex.txt, offset=0)
605
606     if (!cmd) {    
607             dev.off()
608     }
609 }
610
611 plot.dendro.prof <- function(tree, classes, chisqtable, nbbycl = 60, type.dendro = "phylogram", from.cmd = FALSE, bw = FALSE, lab = NULL) {
612     library(ape)
613     library(wordcloud)
614     classes<-classes[classes!=0]
615         classes<-as.factor(classes)
616         sum.cl<-as.matrix(summary(classes, maxsum=1000000))
617         sum.cl<-(sum.cl/colSums(sum.cl)*100)
618         sum.cl<-round(sum.cl,2)
619         sum.cl<-cbind(sum.cl,as.matrix(100-sum.cl[,1]))
620     sum.cl <- sum.cl[,1]
621     tree.order<- as.numeric(tree$tip.label)
622         vec.mat<-NULL
623     row.keep <- select.chi.classe(chisqtable, nbbycl)
624     toplot <- chisqtable[row.keep,]
625     lclasses <- list()
626     for (classe in 1:length(sum.cl)) {
627        ntoplot <- toplot[,classe]
628        names(ntoplot) <- rownames(toplot)
629        ntoplot <- ntoplot[order(ntoplot, decreasing = TRUE)]
630        ntoplot <- round(ntoplot, 0)
631        if (length(toplot) > nbbycl) {
632            ntoplot <- ntoplot[1:nbbycl]
633        }       
634        ntoplot <- ntoplot[which(ntoplot > 0)]
635        #ntoplot <- ntoplot[order(ntoplot)]
636        #ntoplot <- ifelse(length(ntoplot) > nbbycl, ntoplot[1:nbbycl], ntoplot)
637        lclasses[[classe]] <- ntoplot
638     }
639     vec.mat <- matrix(1, nrow = 3, ncol = length(sum.cl))
640         vec.mat[2,] <- 2
641     vec.mat[3,] <- 3:(length(sum.cl)+2)
642     layout(matrix(vec.mat, nrow=3, ncol=length(sum.cl)),heights=c(2,1,6))
643     if (! bw) {
644         col <- rainbow(length(sum.cl))
645         col <- del.yellow(col)
646         col <- col[as.numeric(tree$tip.label)]
647         colcloud <- rainbow(length(sum.cl))
648         colcloud <- del.yellow(colcloud)
649     }
650     label.ori<-tree[[2]]
651     if (!is.null(lab)) {
652         tree$tip.label <- lab
653     } else {
654             tree[[2]]<-paste('classe ',tree[[2]])
655     }
656         par(mar=c(2,1,0,1))
657         plot.phylo(tree,label.offset=0, tip.col=col, type=type.dendro, direction = 'downwards', srt=90, adj = 0.5, cex = 1.4, y.lim=c(-0.3,tree$Nnode))
658         par(mar=c(0,0,0,0))
659         d <- barplot(-sum.cl[tree.order], col=col, names.arg='', axes=FALSE, axisname=FALSE)
660         text(x=d, y=(-sum.cl[tree.order]+3), label=paste(round(sum.cl[tree.order],1),'%'), cex=1.4)
661     for (i in tree.order) {
662         par(mar=c(0,0,1,0),cex=0.7)
663         #wordcloud(names(lclasses[[i]]), lclasses[[i]], scale = c(1.5, 0.2), random.order=FALSE, colors = colcloud[i])
664         yval <- 1.1
665         plot(0,0,pch='', axes = FALSE)
666         vcex <- norm.vec(lclasses[[i]], 1, 2)
667         for (j in 1:length(lclasses[[i]])) {
668             yval <- yval-(strheight( names(lclasses[[i]])[j],cex=vcex[j])+0.02)
669             text(-0.9, yval, names(lclasses[[i]])[j], cex = vcex[j], col = colcloud[i], adj=0)
670         }
671     }
672     if (!from.cmd) {
673         dev.off()
674     }
675     
676 }
677
678 plot.dendro.cloud <- function(tree, classes, chisqtable, nbbycl = 60, type.dendro = "phylogram", from.cmd = FALSE, bw = FALSE, lab = NULL) {
679     library(wordcloud)
680     library(ape)
681     classes<-classes[classes!=0]
682         classes<-as.factor(classes)
683         sum.cl<-as.matrix(summary(classes, maxsum=1000000))
684         sum.cl<-(sum.cl/colSums(sum.cl)*100)
685         sum.cl<-round(sum.cl,2)
686         sum.cl<-cbind(sum.cl,as.matrix(100-sum.cl[,1]))
687     sum.cl <- sum.cl[,1]
688     tree.order<- as.numeric(tree$tip.label)
689         vec.mat<-NULL
690     row.keep <- select.chi.classe(chisqtable, nbbycl)
691     toplot <- chisqtable[row.keep,]
692     lclasses <- list()
693     for (classe in 1:length(sum.cl)) {
694        ntoplot <- toplot[,classe]
695        names(ntoplot) <- rownames(toplot)
696        ntoplot <- ntoplot[order(ntoplot, decreasing = TRUE)]
697        ntoplot <- round(ntoplot, 0)
698        if (length(toplot) > nbbycl) {
699             ntoplot <- ntoplot[1:nbbycl]
700        }
701        ntoplot <- ntoplot[order(ntoplot)]
702        ntoplot <- ntoplot[which(ntoplot > 0)]
703        #ntoplot <- ifelse(length(ntoplot) > nbbycl, ntoplot[1:nbbycl], ntoplot)
704        lclasses[[classe]] <- ntoplot
705     }
706         for (i in 1:length(sum.cl)) vec.mat<-append(vec.mat,1)
707         v<-2
708         for (i in 1:length(sum.cl)) {
709                 vec.mat<-append(vec.mat,v)
710                 v<-v+1
711         }    
712     layout(matrix(vec.mat,length(sum.cl),2),widths=c(1,2))
713     if (! bw) {
714         col <- rainbow(length(sum.cl))[as.numeric(tree$tip.label)]
715         colcloud <- rainbow(length(sum.cl))
716     }
717     par(mar=c(0,0,0,0))
718     label.ori<-tree[[2]]
719     if (!is.null(lab)) {
720         tree$tip.label <- lab
721     } else {
722             tree[[2]]<-paste('classe ',tree[[2]])
723     }
724         plot.phylo(tree,label.offset=0.1,tip.col=col, type=type.dendro)
725     for (i in rev(tree.order)) {
726         par(mar=c(0,0,1,0),cex=0.9)
727         wordcloud(names(lclasses[[i]]), lclasses[[i]], scale = c(2.5, 0.5), random.order=FALSE, colors = colcloud[i])
728     }
729 }
730
731 plot.dendropr <- function(tree, classes, type.dendro="phylogram", histo=FALSE, from.cmd=FALSE, bw=FALSE, lab = NULL, tclasse=TRUE) {
732         classes<-classes[classes!=0]
733         classes<-as.factor(classes)
734         sum.cl<-as.matrix(summary(classes, maxsum=1000000))
735         sum.cl<-(sum.cl/colSums(sum.cl)*100)
736         sum.cl<-round(sum.cl,2)
737         sum.cl<-cbind(sum.cl,as.matrix(100-sum.cl[,1]))
738     tree.order<- as.numeric(tree$tip.label)
739
740
741     if (! bw) {
742         col <- rainbow(nrow(sum.cl))[as.numeric(tree$tip.label)]
743         col <- del.yellow(col)
744         col.bars <- col
745         col.pie <- rainbow(nrow(sum.cl))
746         col.pie <- del.yellow(col.pie)
747             #col.vec<-rainbow(nrow(sum.cl))[as.numeric(tree[[2]])]
748     } else {
749         col = 'black'
750         col.bars = 'grey'
751         col.pie <- rep('grey',nrow(sum.cl))
752     }
753         vec.mat<-NULL
754         for (i in 1:nrow(sum.cl)) vec.mat<-append(vec.mat,1)
755         v<-2
756         for (i in 1:nrow(sum.cl)) {
757                 vec.mat<-append(vec.mat,v)
758                 v<-v+1
759         }
760         par(mar=c(0,0,0,0))
761     if (tclasse) {
762         if (! histo) {
763                 layout(matrix(vec.mat,nrow(sum.cl),2),widths=c(3,1))
764         } else {
765             layout(matrix(c(1,2),1,byrow=TRUE), widths=c(3,2),TRUE)
766         }
767     }
768         par(mar=c(0,0,0,0),cex=1)
769         label.ori<-tree[[2]]
770     if (!is.null(lab)) {
771         tree$tip.label <- lab
772     } else {
773             tree[[2]]<-paste('classe ',tree[[2]])
774     }
775         plot.phylo(tree,label.offset=0.1,tip.col=col, type=type.dendro)
776     #cl.order <- as.numeric(label.ori)
777     #sum.cl[cl.order,1]
778         #for (i in 1:nrow(sum.cl)) {
779     if (tclasse) {
780         if (! histo) {
781             for (i in rev(tree.order)) {
782                 par(mar=c(0,0,1,0),cex=0.7)
783                     pie(sum.cl[i,],col=c(col.pie[i],'white'),radius = 1, labels='', clockwise=TRUE, main = paste('classe ',i,' - ',sum.cl[i,1],'%' ))
784             }
785         } else {
786             par(cex=0.7)
787             par(mar=c(0,0,0,1))
788             to.plot <- sum.cl[tree.order,1]
789             d <- barplot(to.plot,horiz=TRUE, col=col.bars, names.arg='', axes=FALSE, axisname=FALSE)
790             text(x=to.plot, y=d[,1], label=paste(round(to.plot,1),'%'), adj=1.2)
791         }
792     }
793     if (!from.cmd) dev.off()
794         tree[[2]]<-label.ori
795 }
796 #tree <- tree.cut1$tree.cl
797 #to.plot <- di
798 plot.dendro.lex <- function(tree, to.plot, bw=FALSE, lab=NULL, lay.width=c(3,3,2), colbar=NULL, classes=NULL, cmd=FALSE) {
799     tree.order<- as.numeric(tree$tip.label)
800         if (!is.null(classes)) {
801                 classes<-classes[classes!=0]
802                 classes<-as.factor(classes)
803                 sum.cl<-as.matrix(summary(classes, maxsum=1000000))
804                 sum.cl<-(sum.cl/colSums(sum.cl)*100)
805                 sum.cl<-round(sum.cl,2)
806                 sum.cl<-cbind(sum.cl,as.matrix(100-sum.cl[,1]))
807         }
808     par(mar=c(0,0,0,0))
809         if (!is.null(classes)) {
810                 matlay <- matrix(c(1,2,3,4),1,byrow=TRUE)
811                 lay.width <- c(3,2,3,2)
812         } else {
813                 matlay <- matrix(c(1,2,3),1,byrow=TRUE)
814         }
815     layout(matlay, widths=lay.width,TRUE)
816         par(mar=c(3,0,2,4),cex=1)
817         label.ori<-tree[[2]]
818     if (!is.null(lab)) {
819         tree$tip.label <- lab[tree.order]
820     } else {
821             tree[[2]]<-paste('classe ',tree[[2]])
822     }
823     to.plot <- matrix(to.plot[,tree.order], nrow=nrow(to.plot), dimnames=list(rownames(to.plot), colnames(to.plot)))
824     if (!bw) {
825                 col <- rainbow(ncol(to.plot))
826                 col <- del.yellow(col)
827                 if (is.null(colbar)) {
828                 col.bars <- rainbow(nrow(to.plot))
829                 col.bars <- del.yellow(col.bars)
830                 } else {
831                         col.bars <- colbar
832                 }
833     } else {
834         col <- 'black'
835         col.bars <- grey.colors(nrow(to.plot),0,0.8)
836     }
837     col <- col[tree.order]
838         plot.phylo(tree,label.offset=0.2,tip.col=col)
839         if (!is.null(classes)) {
840                 par(cex=0.7)
841                 par(mar=c(3,0,2,1))
842                 to.plota <- sum.cl[tree.order,1]
843                 d <- barplot(to.plota,horiz=TRUE, col=col, names.arg='', axes=FALSE, axisname=FALSE)
844                 text(x=to.plota, y=d[,1], label=paste(round(to.plota,1),'%'), adj=1.2)
845         }
846     par(mar=c(3,0,2,1))
847     d <- barplot(to.plot,horiz=TRUE, col=col.bars, beside=TRUE, names.arg='', space = c(0.1,0.6), axisname=FALSE)
848     c <- colMeans(d)
849     c1 <- c[-1]
850     c2 <- c[-length(c)]
851     cc <- cbind(c1,c2)
852     lcoord <- apply(cc, 1, mean)
853     abline(h=lcoord)
854     if (min(to.plot) < 0) {
855         amp <- abs(max(to.plot) - min(to.plot))
856     } else {
857         amp <- max(to.plot)
858     }
859     if (amp < 10) {
860         d <- 2
861     } else {
862         d <- signif(amp%/%10,1)
863     }
864     mn <- round(min(to.plot))
865     mx <- round(max(to.plot))
866     for (i in mn:mx) {
867         if ((i/d) == (i%/%d)) { 
868             abline(v=i,lty=3)
869         }
870     }    
871     par(mar=c(0,0,0,0))
872     plot(0, axes = FALSE, pch = '')
873     legend(x = 'center' , rownames(to.plot), fill = col.bars)
874     if (!cmd) {
875         dev.off()
876     }
877         tree[[2]]<-label.ori
878 }
879
880 plot.alceste.graph <- function(rdata,nd=3,layout='fruke', chilim = 2) {
881     load(rdata)
882     if (is.null(debsup)) {
883         tab.toplot<-afctable[1:(debet+1),]
884         chitab<-chistabletot[1:(debet+1),]
885     } else {
886         tab.toplot<-afctable[1:(debsup+1),]
887         chitab<-chistabletot[1:(debsup+1),]
888     }
889     rkeep<-select_point_chi(chitab,chilim)
890     tab.toplot<-tab.toplot[rkeep,]
891     chitab<-chitab[rkeep,]
892     dm<-dist(tab.toplot,diag=TRUE,upper=TRUE)
893     cn<-rownames(tab.toplot)
894     cl.toplot<-apply(chitab,1,which.max)
895     col<-rainbow(ncol(tab.toplot))[cl.toplot]
896     library(igraph)
897     g1 <- graph.adjacency(as.matrix(dm), mode = 'lower', weighted = TRUE)
898     g.max<-minimum.spanning.tree(g1)
899     we<-(rowSums(tab.toplot)/max(rowSums(tab.toplot)))*2
900     #lo <- layout.fruchterman.reingold(g.max,dim=nd)
901     lo<- layout.kamada.kawai(g.max,dim=nd)
902     print(nrow(tab.toplot))
903     print(nrow(chitab))
904     print(length(we))
905     print(length(col))
906     print(length(cn))
907     if (nd == 3) {
908         rglplot(g.max, vertex.label = cn, vertex.size = we*3, edge.width = 0.5, edge.color='black', vertex.label.color = col,vertex.color = col, layout = lo, vertex.label.cex = 1)
909     } else if (nd == 2) {
910         plot(g.max, vertex.label = cn, vertex.size = we, edge.width = 0.5, edge.color='black', vertex.label.color = col,vertex.color = col, layout = lo, vertex.label.cex = 0.8)
911     }
912
913 }
914
915 make.simi.afc <- function(x,chitable,lim=0, alpha = 0.1, movie = NULL) {
916     library(igraph)
917     chimax<-as.matrix(apply(chitable,1,max))
918     chimax<-as.matrix(chimax[,1][1:nrow(x)])
919     chimax<-cbind(chimax,1:nrow(x))
920     order_chi<-as.matrix(chimax[order(chimax[,1],decreasing = TRUE),])
921     if ((lim == 0) || (lim>nrow(x))) lim <- nrow(x)
922     x<-x[order_chi[,2][1:lim],]
923     maxchi <- chimax[order_chi[,2][1:lim],1]
924     #-------------------------------------------------------
925     limit<-nrow(x)
926     distm<-dist(x,diag=TRUE)
927     distm<-as.matrix(distm)
928     g1<-graph.adjacency(distm,mode='lower',weighted=TRUE)
929     g1<-minimum.spanning.tree(g1)
930     lo<-layout.kamada.kawai(g1,dim=3)
931     lo <- layout.norm(lo, -3, 3, -3, 3, -3, 3)
932     mc<-rainbow(ncol(chistabletot))
933     chitable<-chitable[order_chi[,2][1:lim],]
934     cc <- apply(chitable, 1, which.max)
935     cc<-mc[cc]
936     #mass<-(rowSums(x)/max(rowSums(x))) * 5
937     maxchi<-norm.vec(maxchi, 0.03, 0.3)
938     rglplot(g1,vertex.label = vire.nonascii(rownames(x)),vertex.label.color= cc,vertex.label.cex = maxchi, vertex.size = 0.1, layout=lo, rescale=FALSE)
939     text3d(lo[,1], lo[,2],lo[,3], rownames(x), cex=maxchi, col=cc)
940     #rgl.spheres(lo, col = cc, radius = maxchi, alpha = alpha)
941     rgl.bg(color = c('white','black'))
942     if (!is.null(movie)) {
943         require(tcltk)
944         ReturnVal <- tkmessageBox(title="RGL 3 D",message="Cliquez pour commencer le film",icon="info",type="ok")
945
946         movie3d(spin3d(axis=c(0,1,0),rpm=6), movie = 'film_graph', frames = "tmpfilm", duration=10, clean=TRUE, top = TRUE, dir = movie)
947         ReturnVal <- tkmessageBox(title="RGL 3 D",message="Film fini !",icon="info",type="ok")
948     }
949         while (rgl.cur() != 0)
950                 Sys.sleep(1)
951
952 }
953
954 # from igraph
955 norm.vec <- function(v, min, max) {
956
957   vr <- range(v)
958   if (vr[1]==vr[2]) {
959     fac <- 1
960   } else {
961     fac <- (max-min)/(vr[2]-vr[1])
962   }
963   (v-vr[1]) * fac + min
964 }
965
966
967 vire.nonascii <- function(rnames) {
968     print('vire non ascii')
969     couple <- list(c('é','e'),
970                 c('è','e'),
971                 c('ê','e'),
972                 c('ë','e'),
973                 c('î','i'),
974                 c('ï','i'),
975                 c('ì','i'),
976                 c('à','a'),
977                 c('â','a'),
978                 c('ä','a'),
979                 c('á','a'),
980                 c('ù','u'),
981                 c('û','u'),
982                 c('ü','u'),
983                 c('ç','c'),
984                 c('ò','o'),
985                 c('ô','o'),
986                 c('ö','o'),
987                 c('ñ','n')
988                 )
989
990     for (c in couple) {
991         rnames<-gsub(c[1],c[2], rnames)
992     }
993     rnames
994 }
995
996
997
998 #par(mar=c(0,0,0,0))
999 #layout(matrix(c(1,2),1,byrow=TRUE), widths=c(3,2),TRUE)
1000 #par(mar=c(1,0,1,0), cex=1)
1001 #plot.phylo(tree,label.offset=0.1)
1002 #par(mar=c(0,0,0,1))
1003 #to.plot <- sum.cl[cl.order,1]
1004 #d <- barplot(to.plot,horiz=TRUE, names.arg='', axes=FALSE, axisname=FALSE)
1005 #text(x=to.plot, y=d[,1], label=round(to.plot,1), adj=1.2)
1006
1007 make.afc.attributes <- function(rn, afc.table, contafc, clnb, column = FALSE, x=1, y=2) {
1008     if (!column){
1009         nd <- clnb - 1
1010         afc.res <- afc.table$ligne
1011         #tokeep <- which(row.names(afc.res) %in% rn)
1012         afc.res <- afc.res[rn,]
1013         debcor <- (nd*2) + 1
1014         cor <- afc.res[,debcor:(debcor+nd-1)][,c(x,y)]
1015         debctr <- (nd*3) + 1
1016         ctr <- afc.res[,debctr:(debctr+nd-1)][,c(x,y)]
1017         massdeb <- (nd*4) + 1
1018         mass <- afc.res[,massdeb]
1019         chideb <- massdeb + 1
1020         chi <- afc.res[,chideb]
1021         inertiadeb <- chideb + 1
1022         inertia <- afc.res[,inertiadeb]
1023         frequence <- rowSums(contafc[rn,])
1024     }
1025     res <- list(frequence=frequence, cor, ctr, mass = mass, chi=chi, inertia=inertia)
1026     return(res)
1027 }
1028
1029
1030 afctogexf <- function(fileout, toplot, classes, clnb, sizes, nodes.attr=NULL) {
1031     toplot <- toplot[,1:3]
1032     toplot[,3] <- 0
1033     #toplot <- afc$rowcoord[1:100,1:3]
1034     #toplot[,3] <- 0
1035     #rownames(toplot)<-afc$rownames[1:100]
1036     cc <- rainbow(clnb)[classes]
1037     cc <- t(sapply(cc, col2rgb, alpha=TRUE))
1038     #sizes <- apply(chistabletot[1:100,], 1, max)
1039     
1040     nodes <- data.frame(cbind(1:nrow(toplot), rownames(toplot)))
1041     colnames(nodes) <- c('id', 'label')
1042     nodes[,1] <- as.character(nodes[,1])
1043     nodes[,2] <- as.character(nodes[,2])
1044     #nodes attributs
1045     if (! is.null(nodes.attr)) {
1046         nodesatt <- as.data.frame(nodes.attr)
1047     } else {
1048         nodesatt <- data.frame(cbind(toplot[,1],toplot[,2]))
1049     }
1050     #make axes
1051     edges<-matrix(c(1,1),ncol=2)
1052     xmin <- min(toplot[,1])
1053     xmax <- max(toplot[,1])
1054     ymin <- min(toplot[,2])
1055     ymax <- max(toplot[,2])
1056     nodes<-rbind(nodes, c(nrow(nodes)+1, 'F1'))
1057     nodes<-rbind(nodes, c(nrow(nodes)+1, 'F1'))
1058     nodes<-rbind(nodes, c(nrow(nodes)+1, 'F2'))
1059     nodes<-rbind(nodes, c(nrow(nodes)+1, 'F2'))
1060     nodesatt<-rbind(nodesatt, c(0,0))
1061     nodesatt<-rbind(nodesatt, c(0,0))
1062     nodesatt<-rbind(nodesatt, c(0,0))
1063     nodesatt<-rbind(nodesatt, c(0,0))
1064     toplot <- rbind(toplot, c(xmin, 0,0))
1065     toplot <- rbind(toplot, c(xmax,0,0))
1066     toplot <- rbind(toplot, c(0,ymin,0))
1067     toplot <- rbind(toplot, c(0,ymax,0))
1068     cc <- rbind(cc, c(255,255,255,1))
1069     cc <- rbind(cc, c(255,255,255,1))
1070     cc <- rbind(cc, c(255,255,255,1))
1071     cc <- rbind(cc, c(255,255,255,1))
1072     sizes <- c(sizes, c(0.5, 0.5, 0.5, 0.5))
1073     edges <- rbind(edges, c(nrow(nodes)-3, nrow(nodes)-2))
1074     edges <- rbind(edges, c(nrow(nodes)-1, nrow(nodes)))
1075     write.gexf(nodes, edges, output=fileout, nodesAtt=nodesatt, nodesVizAtt=list(color=cc, position=toplot, size=sizes))
1076 }
1077
1078 simi.to.gexf <- function(fileout, graph.simi, nodes.attr = NULL) {
1079         lo <- graph.simi$layout
1080         if (ncol(lo) == 3) {
1081                 lo[,3] <- 0
1082         } else {
1083                 lo <- cbind(lo, rep(0,nrow(lo)))
1084         }
1085         g <- graph.simi$graph
1086         nodes <- data.frame(cbind(1:nrow(lo), V(g)$name))
1087         colnames(nodes) <- c('id', 'label')
1088         print(nodes)
1089         if (! is.null(nodes.attr)) {
1090                 nodesatt <- as.data.frame(nodes.attr)
1091         } else {
1092                 nodesatt <- data.frame(cbind(lo[,1],lo[,2]))
1093         }
1094         edges <- as.data.frame(get.edges(g, c(1:ecount(g))))
1095         col <- rep('red', nrow(lo))
1096         col <- t(sapply(col, col2rgb, alpha=TRUE))
1097         write.gexf(nodes, edges, output=fileout, nodesAtt=nodesatt, nodesVizAtt=list(color=col,position=lo))
1098 }
1099
1100
1101 graph.to.file <- function(grah.simi, nodesfile = NULL, edgesfile = NULL, community = FALSE, color = NULL, sweight = NULL) {
1102         require(igraph)
1103         g <- graph.simi$graph
1104         V(g)$weight <- graph.simi$eff
1105         V(g)$x <- graph.simi$layout[,1]
1106         V(g)$y <- graph.simi$layout[,2]
1107         if (ncol(graph.simi$layout) == 3) {
1108                 V(g)$z <- graph.simi$layout[,3]
1109         }
1110         if (community) {
1111                 member <- graph.simi$communities$membership
1112                 col <- rainbow(max(member))
1113                 v.colors <- col[member]
1114                 v.colors <- col2rgb(v.colors)
1115                 V(g)$r <- v.colors[1,]
1116                 V(g)$g <- v.colors[2,]
1117                 V(g)$b <- v.colors[3,]
1118         }
1119         if (!is.null(color)) {
1120                 v.colors <- col2rgb(color)
1121                 V(g)$r <- v.colors[1,]
1122                 V(g)$g <- v.colors[2,]
1123                 V(g)$b <- v.colors[3,]          
1124         }
1125         if (!is.null(sweight)) {
1126                 V(g)$sweight <- sweight
1127         }
1128         df <- get.data.frame(g, what='both')
1129         if (!is.null(nodesfile)) {
1130                 write.table(df$vertices, nodesfile, sep='\t')
1131         }
1132         if (!is.null(edgesfile)) {
1133                 write.table(df$edges, edgesfile, sep='\t')
1134         }
1135         if (is.null(edgesfile) & is.null(nodesfile)) {
1136                 df
1137         }
1138 }