-
Notifications
You must be signed in to change notification settings - Fork 13
/
Copy pathVisulz_bulkATAC_network.R
85 lines (66 loc) · 2.28 KB
/
Visulz_bulkATAC_network.R
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
# MESSAGE -----------------------------------------------------------------
#
# author: Yulin Lyu
# email: lvyulin@pku.edu.cn
#
# require: R whatever
#
# ---
# * 1. Load packages ------------------------------------------------------
setwd("exampleData/ATAC")
# grammar
library(tidyverse)
library(magrittr)
library(glue)
library(data.table)
# for graph
library(XML)
# graphics
library(ggplot2)
library(ggsci)
library(RColorBrewer)
library(viridis)
library(ggrepel)
library(scales)
library(ggnewscale)
# * 2. Load data ----------------------------------------------------------
eData <- fread("network/4_TFgene.csv")
vData <- fread("network/4_node.csv")
setkey(vData, Id)
graphData <- xmlParse("network/g4.gexf") %>% xmlToList()
vData[map(graphData$graph$nodes, ".attrs") %>% map("id") %>% unlist() %>% unname(), `:=` (
x = map(graphData$graph$nodes, "position") %>% map("x") %>% unlist() %>% unname() %>% as.numeric(),
y = map(graphData$graph$nodes, "position") %>% map("y") %>% unlist() %>% unname() %>% as.numeric()
)][]
eData[, `:=` (x = vData[Source, x], y = vData[Source, y], xend = vData[Target, x], yend = vData[Target, y])][]
# * 3. Plot ---------------------------------------------------------------
vData[, color := "other"][type == "tf", color := cut(Weight, 5)][]
colorName <- vData[type == "tf"][order(-Weight), color] %>% unique() %>% c("other")
geneColor <- structure(
c(viridis(5, begin = .2, direction = -1), "gray80"),
names = colorName
)
ggplot() +
geom_segment(
data = eData, aes(x, y, xend = xend, yend = yend, size = Weight),
color = "gray90", alpha = .5, show.legend = F) +
scale_radius(range = c(0.1, 0.5)) +
new_scale("size") +
geom_point(
data = vData, aes(x, y, color = color, size = Weight),
show.legend = F) +
scale_color_manual(values = geneColor) +
scale_radius(range = c(1, 12)) +
new_scale("size") +
geom_text(
data = vData[color %in% colorName[1:3]], aes(x, y, label = Id, size = Weight),
fontface = 3, show.legend = F) +
scale_size(range = c(2, 4)) +
theme(
aspect.ratio = 1,
panel.background = element_blank(),
panel.grid = element_blank(),
axis.ticks = element_blank(),
axis.text = element_blank(),
axis.title = element_blank())
ggsave("graphics/g4_network.png", width = 6, height = 6, scale = 1.5)