Skip to content
Snippets Groups Projects
Commit 45256da9 authored by UMEC Mathieu's avatar UMEC Mathieu
Browse files

Merge branch 'tesselation_voronoi' into mapping_ma_r

parents 5b3f2f09 07e7e5cf
No related branches found
No related tags found
No related merge requests found
Hide whitespace changes
Inline Side-by-side
Tesselation_de_voronoi_modif.ipynb 0 → 100644
+ 402
− 0
View file @ 45256da9
%% Cell type:code id:d5fbe906-9395-4c84-ac8d-5a2d774bf39b tags:
``` R
#author: "Clément Frainay"
#output: html_document
#knitr::opts_chunk$set(echo = TRUE)
```
%% Cell type:code id:8c08fd21-ec30-4f74-83e8-768ba4540915 tags:
``` R
library(igraph)
library(scales)
library(openxlsx)
```
%% Output
Warning message:
"package 'igraph' was built under R version 3.6.3"
Attaching package: 'igraph'
The following objects are masked from 'package:stats':
decompose, spectrum
The following object is masked from 'package:base':
union
Warning message:
"package 'openxlsx' was built under R version 3.6.3"
%% Cell type:code id:449d0f84-38f5-4b22-adff-fad36d012f10 tags:
``` R
#import networks in gml ( graph format)
g <- read_graph("C:/Users/mumec/Desktop/Mini_codes/CompoundNet_KEGG2SBML",format = "gml")
```
%% Output
Error in match.arg(arg = arg, choices = choices, several.ok = several.ok): 'arg' should be one of "edgelist", "pajek", "ncol", "lgl", "graphml", "dimacs", "graphdb", "gml", "dl"
Traceback:
1. read_graph("C:/Users/mumec/Desktop/Mini_codes/CompoundNet_KEGG2SBML",
. format = "txt")
2. igraph.match.arg(format)
3. match.arg(arg = arg, choices = choices, several.ok = several.ok)
4. stop(gettextf("'arg' should be one of %s", paste(dQuote(choices),
. collapse = ", ")), domain = NA)
%% Cell type:code id:a975c1b3-2cc6-47d8-a744-79d5165779ba tags:
``` R
corresp <- read.csv("C:/Users/mumec/Desktop/Mini_codes/100_metabolites_KEGG_ID.csv", sep=";")
```
%% Cell type:code id:704eb832-7948-43d0-8a4e-5e9ca49de549 tags:
``` R
Name_cor <- as.character(corresp$ NAME)
ID_cor <- as.character(corresp$ KEGG)
```
%% Cell type:code id:ccaf5941-4624-428b-8a53-beb73af67f93 tags:
``` R
#extract main connected component (otherwise, inf in the matrix)
g <- decompose(g, mode="weak", max.comps = 1)[[1]]
m <- distances(g) # Distance matrix without ID
rownames(m) <- V(g)$label # recup of name
colnames(m) <- V(g)$label
```
%% Cell type:code id:f48cf29d-e94f-4ba0-8b71-d54c21bf7e65 tags:
``` R
## Center-Nodes computation
grid_size = c(10,10)# put the dimension of the grid
nclusters = grid_size[1]*grid_size[2]
centers <- data.frame(ID_cor)
```
%% Cell type:code id:2e5d1af7-b4a9-4897-8ec1-26d8663cdc3b tags:
``` R
## Neighborhood assignments
submed <- m[rownames(m) %in% centers[,1],]
minDist <- apply(submed,2, min)
minDistMatrix <- matrix(rep(minDist, each = nrow(submed)), nrow = nrow(submed), ncol = ncol(submed))
assignMatrix <- (submed==minDistMatrix)
assign <- apply(assignMatrix, 1, function(row) colnames(submed)[row])
vornoiBridge<-colnames(assignMatrix)[colSums(assignMatrix) > 1]
```
%% Cell type:code id:cab448f4-8181-4fcd-a030-787b277e40f5 tags:
``` R
## Subgraph building
Name_of_out_file="tesselation_test_modif"
df <- data.frame(matrix(ncol = 5, nrow = length(assign)))
colnames(df) <- c("center", "nodes", "edges", "outEdges","nbCC")
pdf(file=paste0(Name_of_out_file,nclusters,".pdf"),width=50, height=50)
#par(mfrow=c(grid_size[1],grid_size[2]), cex.main = 4)
par(mfrow=c(7,7), cex.main = 4)
for(i in 1:length(assign)){
medoidId = names(assign[i]) # on récupére l'id e question ici on peux modifier l'ID par son nom usuel
#get neighborhood vertices
clusterVertex = V(g)[V(g)$label%in%assign[[i]]]
indice <- which(ID_cor == medoidId)
Name_ID = Name_cor[indice][1]
#build induced subgraph
clusterGraph = induced_subgraph(g,clusterVertex,impl = "auto")
#get edges that goes outside the cluster
incidentEdges <- incident_edges(g, clusterVertex, mode = "all")
incidentEdges <- setdiff(incidentEdges,E(clusterGraph))
#cluser data
df$center[i] = Name_ID # dzdzdzdazdazdazdazdazdazdazd
df$nodes[i] = length(V(clusterGraph)) # number of graph vertices
df$edges[i] = length(E(clusterGraph)) # Number of graph edges
df$outEdges[i] = length(incidentEdges) # Number of out edges
df$nbCC[i] = count_components(clusterGraph, mode = "weak") # Counting related components
df$diameter[i]= diameter(clusterGraph)
#clusters visualization
V(clusterGraph)$size <- 6
V(clusterGraph)$shape <- "circle"
V(clusterGraph)$color <- "aquamarine"
V(clusterGraph)[V(clusterGraph)$label == medoidId]$size <- 15
V(clusterGraph)[V(clusterGraph)$label == medoidId]$color <- "khaki1"
V(clusterGraph)[V(clusterGraph)$label %in% vornoiBridge]$shape <- "square"
V(clusterGraph)[V(clusterGraph)$label %in% vornoiBridge]$size <- 4
V(clusterGraph)[V(clusterGraph)$label %in% vornoiBridge]$color <- "mediumorchid1"
plot(clusterGraph,vertex.label=NA, main =paste(Name_ID,collapse = '\n'))
box(col = "black",lwd = 7)
}
dev.off()
```
%% Output
**png:** 2
\textbf{png:} 2
%% Cell type:code id:b6c838d4-3bdf-4a2c-a987-f6e313e8d7a6 tags:
``` R
# Sub-graphs analysis
print(df)
```
%% Output
center nodes edges outEdges nbCC diameter
1 bilirubin 244 945 243 1 13
2 sphinganine 176 695 173 1 13
3 pipecolinic acid 31 121 31 1 6
4 pyruvic acid 364 1054 363 1 16
5 L-glutamic acid 475 1783 475 1 13
6 glycine 305 782 305 1 14
7 succinic acid 267 1085 265 1 11
8 L-aspartic acid 247 878 247 1 13
9 L-lysine 99 210 99 1 11
10 L-serine 251 823 248 1 14
11 L-arginine 101 214 101 1 9
12 L-tryptophan 114 247 114 1 14
13 L-ornithine 5 10 5 1 2
14 L-methionine 232 617 232 1 15
15 L-tyrosine 126 316 126 1 12
16 beta-alanine 18 45 18 1 6
17 L-cysteine 335 1198 332 1 13
18 serotonin 50 112 50 1 10
19 trans-urocanic acid 2 2 2 1 1
20 creatine 51 120 51 1 6
21 spermidine 8 24 8 1 4
22 L-citrulline 96 203 96 1 9
23 glycine betaine 51 124 51 1 6
24 oleic acid 3 2 3 1 1
25 1-methylnicotinamide 51 120 51 1 6
26 4-acetamidobutanoic acid 2 2 2 1 1
27 9H-xanthine 10 36 10 1 4
28 5-hydroxyindoleacetic acid 40 95 40 1 5
29 octadecanoic acid 4 3 4 1 3
30 cholic acid 96 199 96 1 9
31 sarcosine 116 393 116 1 9
32 4-methyl-2-oxovaleric acid 39 85 39 1 8
33 hexadecanoic acid 172 475 172 1 13
34 taurine 25 49 25 1 6
35 riboflavin 164 551 164 1 14
36 D-glyceric acid 20 56 20 1 10
37 all-cis-5_8_11_14_17-icosapentaenoic acid 22 33 22 1 6
38 octanoic acid 37 144 37 1 7
39 hypoxanthine 15 27 15 1 6
40 linoleic acid 172 661 169 1 12
41 N6_N6_N6-trimethyl-L-lysine 2 2 2 1 1
42 uridine 148 650 148 1 13
43 inosine 26 60 26 1 10
44 Nn-methyl-L-histidine 2 4 2 1 1
45 trans-4-hydroxy-L-proline 37 75 37 1 6
46 5-oxo-L-proline 97 223 97 1 8
47 uracil 12 26 12 1 5
48 phosphocholine 65 243 65 1 10
49 choline 255 853 252 1 16
50 L-leucine 94 192 94 1 9
51 fumaric acid 27 50 27 1 7
52 L-histidine 102 212 102 1 10
53 indole-3-acetic acid 2 4 2 1 1
54 L-thyroxine 3 4 3 1 2
55 N-acetylputrescine 274 1055 271 1 16
56 3-methyl-2-oxobutanoic acid 43 112 43 1 8
57 (S)-malic acid 27 61 27 1 7
58 L-proline 111 232 111 1 13
59 adenine 107 231 107 1 10
60 beta-D-fructofuranose 37 121 37 1 9
61 nicotinamide 49 125 49 1 7
62 L-asparagine 159 434 159 1 12
63 citric acid 108 238 107 1 13
64 (2-hydroxyphenyl)acetic acid 3 4 3 1 2
65 N_N-dimethylglycine 70 177 70 1 9
66 L-threonine 97 204 97 1 9
67 cholesterol 352 1261 348 1 15
68 L-valine 102 211 102 1 11
69 4-guanidinobutanoic acid 2 4 2 1 1
%% Cell type:code id:d0a0f33e-9f25-45f6-a4e2-c0863d1ff8bf tags:
``` R
outfile <- "C:/Users/mumec/Desktop/Mini_codes/data_out_tesselation.xlsx"
write.xlsx(df, outfile, sheetName = "Tesselation Voronoï", row.names = FALSE)
```
%% Cell type:code id:262f1c2c-cd8a-4e5b-9503-0c1fa36b4bbf tags:
``` R
df$density <- (df$diameter/(df$nodes*(df$nodes-1)))
plot(df$nodes ~ df$diameter,
ylab="Number of cell nodes",
xlab="Cell diameter",
main="Voronoï cells topology",
cex = rescale(df$density, to = c(1, 5)))
```
%% Output
%% Cell type:code id:70da4ec7-437d-4ee4-9d5a-98f633cbcf9c tags:
``` R
df$density <- (df$edges/(df$nodes*(df$nodes-1)))
plot(df$nodes ~ df$edges,
ylab="Number of cell nodes",
xlab="Total edges number",
main="Voronoï cells topology",
cex = rescale(df$density, to = c(1, 5)))
```
%% Output
%% Cell type:code id:2f3c9af6-56c6-4190-996a-29089c35e5cb tags:
``` R
boxplot(df$nodes,
ylab="Number of nodes",
main="BoxPlot of number of nodes"
)
```
%% Output
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment