-
Notifications
You must be signed in to change notification settings - Fork 0
/
Copy patholeg_rnaseq_volcanoPlot.R
70 lines (63 loc) · 2.86 KB
/
oleg_rnaseq_volcanoPlot.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
data=x[,c('PPEE','PostFC','X')]
data = data[data$X!="N/A",]
data = data[complete.cases(data),]
data[is.na(data[,1]),1] = 1
data[data[,1]==0,1] = 1e-16
data[,2] = log2(data[,2])
dds_res = data.frame(log2FoldChange=data[,2],padj=data[,1])
library(graphics)
library(gplots)
library(ggplot2)
pdf("oleg_rnaseq_volcanoPlot.pdf")
plot(dds_res$log2FoldChange,-log10(dds_res$padj),xlab=expression('Log'[2]*' Fold Change ( shNT / shH2AFV )'),
ylab=expression('-Log'[10]*' Q-values'),col=alpha("grey",.5))
abline(v=-.5,lty = 2,col="grey")
abline(v=.5,lty = 2,col="grey")
abline(h=-log10(0.05),lty = 2,col="grey")
points(dds_res$log2FoldChange[abs(dds_res$log2FoldChange)>.6 & dds_res$padj<0.05],
-log10(dds_res$padj)[abs(dds_res$log2FoldChange)>.6 & dds_res$padj<0.05],
col="red")
legend("topright", paste("shNT:",length(which(dds_res$log2FoldChange>.6 & dds_res$padj<0.05))), bty="n")
legend("topleft", paste("shH2AFV:",length(which(dds_res$log2FoldChange<(-.6) & dds_res$padj<0.05))), bty="n")
dev.off()
########
#REDO .37
x = read.csv("Master_tab_o_202_143_400_R.csv")
data=x[,c('PPEE','PostFC','X')]
data = data[data$X!="N/A",]
data = data[complete.cases(data),]
data[is.na(data[,1]),1] = 1
data[data[,1]==0,1] = 1e-16
data[,2] = log2(data[,2])
dds_res = data.frame(log2FoldChange=data[,2],padj=data[,1])
library(graphics)
library(gplots)
library(ggplot2)
pdf("oleg_rnaseq_volcanoPlot_30perc_log2FC37_withNUMBERS.pdf")
plot(dds_res$log2FoldChange,-log10(dds_res$padj),xlab=expression('Log'[2]*' Fold Change ( shNT / shH2AFV )'),
ylab=expression('-Log'[10]*' Q-values'),col=alpha("grey",.5))
abline(v=-.37,lty = 2,col="grey")
abline(v=.37,lty = 2,col="grey")
abline(h=-log10(0.05),lty = 2,col="grey")
points(dds_res$log2FoldChange[abs(dds_res$log2FoldChange)>.37 & dds_res$padj<0.05],
-log10(dds_res$padj)[abs(dds_res$log2FoldChange)>.37 & dds_res$padj<0.05],
col="red")
legend("topright", paste("shNT:",length(which(dds_res$log2FoldChange>.37 & dds_res$padj<0.05))), bty="n")
legend("topleft", paste("shH2AFV:",length(which(dds_res$log2FoldChange<(-.37) & dds_res$padj<0.05))), bty="n")
abline(v=-.37,lty = 2,col="grey")
abline(v=.37,lty = 2,col="grey")
abline(h=-log10(0.05),lty = 2,col="grey")
dev.off()
pdf("oleg_rnaseq_volcanoPlot_30perc_log2FC37_withoutNUMBERS.pdf")
plot(dds_res$log2FoldChange,-log10(dds_res$padj),xlab=expression('Log'[2]*' Fold Change ( shNT / shH2AFV )'),
ylab=expression('-Log'[10]*' Q-values'),col=alpha("grey",.5))
abline(v=-.37,lty = 2,col="grey")
abline(v=.37,lty = 2,col="grey")
abline(h=-log10(0.05),lty = 2,col="grey")
points(dds_res$log2FoldChange[abs(dds_res$log2FoldChange)>.37 & dds_res$padj<0.05],
-log10(dds_res$padj)[abs(dds_res$log2FoldChange)>.37 & dds_res$padj<0.05],
col="red")
abline(v=-.37,lty = 2,col="grey")
abline(v=.37,lty = 2,col="grey")
abline(h=-log10(0.05),lty = 2,col="grey")
dev.off()