-
Notifications
You must be signed in to change notification settings - Fork 13
/
Copy pathVisulz_bulkRNA_MAplot.R
75 lines (60 loc) · 1.91 KB
/
Visulz_bulkRNA_MAplot.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
# MESSAGE -----------------------------------------------------------------
#
# author: Yulin Lyu
# email: lvyulin@pku.edu.cn
#
# require: R whatever
#
# ---
# 1. Load packages --------------------------------------------------------
setwd("exampleData/RNA")
# grammar
library(tidyverse)
library(magrittr)
library(glue)
library(data.table)
# analysis
library(DESeq2)
# graphics
library(ggplot2)
library(ggrepel)
library(ggsci)
library(scales)
library(latex2exp)
# 2. Load data ------------------------------------------------------------
diffData <- fread("mid/ES_vs_Fib.DEG.csv")
colnames(diffData)[1] <- "gene"
diffData[is.na(padj), padj := 1][]
diffData[, p := -log10(padj)][]
diffData[, baseMean := log2(baseMean)][]
# 3. Plot -----------------------------------------------------------------
diffData[, type := "ns"][]
diffData[log2FoldChange > 3 & padj < 0.05, type := "up"][log2FoldChange < -3 & padj < 0.05, type := "down"][]
labelGene <- diffData[order(p, decreasing = T)][type == "up"][1:10]
pal_nejm()(8) %>% show_col()
typeColor <- structure(
c(pal_nejm()(2), "gray80"),
names = c("up", "down", "ns")
)
ggplot(diffData, aes(x = baseMean, y = log2FoldChange)) +
geom_point(aes(color = type, size = p), show.legend = F) +
geom_hline(yintercept = 0, color = "black") +
geom_hline(yintercept = 3, color = "gray60", linetype = "dashed") +
geom_hline(yintercept = -3, color = "gray60", linetype = "dashed") +
geom_text_repel(
data = labelGene, aes(label = gene),
size = 3, fontface = 3,
nudge_x = .5, nudge_y = .5,
max.overlaps = 100) +
scale_radius(range = c(.1, 2)) +
scale_color_manual(values = typeColor) +
labs(
x = TeX("$log_{2}(base\\,Mean)$"),
y = TeX("$log_{2}(Fold\\,Change)$")) +
theme(
aspect.ratio = 1,
panel.background = element_blank(),
panel.grid = element_blank(),
axis.line = element_line()
)
ggsave("plot/ES_vs_Fib.MAplot.pdf", width = 4, height = 4)