-
Notifications
You must be signed in to change notification settings - Fork 1
/
Copy pathlong_read_TGS_AS_plot.R
119 lines (100 loc) · 3.37 KB
/
long_read_TGS_AS_plot.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
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
---
title: "Untitled"
author: "xumengying"
date: "2020/11/27"
output: html_document
---
```{r setup, include=FALSE}
knitr::opts_chunk$set(echo = TRUE)
```
# 1. load data
```{r}
files <- list.files("/mnt/data5/xumengying/proj/NSC_yan_2020_10_30/analysis/04.AS_events/02.SUPP2/minimap/v2", full.names = T, pattern = ".ioe")
supp2 <- lapply(1:length(files), function(x){
a <- read.table(
files[x],
header = T,
sep = "\t",
stringsAsFactors = F)
a$type <- unlist(strsplit(basename(files[x]), "_"))[2]
a$stage <- unlist(strsplit(basename(files[x]), "_"))[1]
return(a)
} )
supp <- do.call(rbind ,supp2)
supp$stage <- factor(supp$stage,
levels = c(paste0("barcode0", rep(1:5)), "TGS"),
labels = c("E15.5", "È17.5", "P1.5", "P8", "Adult", "Complex"))
```
# 2. find overlap
```{r fig.width=8.5,fig.height=5}
sample <- unique(supp$stage)
veen_list <- lapply(1:5, function(x){
do.call(rbind ,strsplit(supp[supp$stage == sample[x],]$event_id, ";"))[,2]
})
names(veen_list) <- sample[1:5]
library(UpSetR)
pdf("./upset.pdf",8,5)
upset(fromList(veen_list), nsets = 9, order.by = "freq")
dev.off()
```
## 2.1 veen plot
```{r }
futile.logger::flog.threshold(futile.logger::ERROR, name = "VennDiagramLogger")
library(VennDiagram)
g = venn.diagram(
x = list( E15_5 = veen_list[[1]],
E17_5 = veen_list[[2]],
P1_5 = veen_list[[3]],
P8 = veen_list[[4]],
P60 = veen_list[[5]]),
category.names = c("E15.5", "E17.5", "P1.5","P8", "Adult"),
filename = NULL,
output=FALSE,
# Circles
lwd = 2,
lty = 'blank',
fill = c("#F8766D", "#A3A500", "#00BF7D", "#00B0F6", "#E76BF3"),#c("dodgerblue", "goldenrod1", "darkorange1", "seagreen3", "orchid3"),
#
# # Numbers
cex = .5,
fontface = "bold",
# fontfamily = "sans",
#
# # Set names
cat.cex = 1,
cat.fontface = "bold",
cat.default.pos = "outer",
# cat.pos = c(-27, 27, 135),
# cat.dist = c(0.055, 0.055, 0.085),
# cat.fontfamily = "sans",
# rotation = 1,
#
margin = 0.07,
print.mode=c("percent")
)
pdf("./veen.pdf",width = 4, height = 4)
grid::grid.draw(g)
dev.off()
```
# 3. sum as events
```{r}
supp_anno <- supp
supp_anno %>% dplyr::group_by(stage) %>% mutate(total_event = length(unique(event_id)))-> supp_anno
supp_anno %>% dplyr::group_by(stage,type) %>% mutate(Events.1 = length(unique(event_id)), Events = round(length(unique(event_id))/total_event *100 , 2)) -> supp_anno
plot.events <- unique(supp_anno[,c("type", "stage","Events.1","Events","total_event")])
plot.events <- plot.events[plot.events$stage != "Complex", ]
```
```{r fig.width=10, fig.height=4}
color <- "Dark2"
pdf("./supp2_as_events.pdf", 10,5)
ggplot(plot.events, aes(stage, Events.1, fill = type)) +
geom_col(position=position_dodge())+
labs(x = "Stage", y = "Number")+
theme_classic() +
scale_fill_brewer(palette=color)+
# scale_fill_brewer()+
geom_text(aes(label =round(Events,0), y = Events.1+ 25), position=position_dodge(0.9), vjust = 0.1, size = 4)+
theme(axis.title = element_text(size = 16),
axis.text = element_text(size = 12))
dev.off()
```