cell.jpg
跟着CELL学作图|1.火山图
“实践是检验真理的唯一标准。”
“复现是学习R语言的最好办法。”
2021.4.12_1
这篇2020年发表在cell上关于新冠的组学文章里面有大量的生信内容。今天带大家复现其中的一个Supplemental Figure：火山图。
22
本文代码及示例数据领取：后台回复“210412”
2021.4.12_2
这图确实比一般的火山图美观且简洁。
火山图的意义
火山图可用于展示两组样本间基因表达水平差异的分布状况。
横轴log2 fold change差异表达倍数（Fold Change值，简称FC），差异越大的基因分布X轴在两端。
纵坐标用-log10 p-value表示，对P值进行-log10的转化。转化后，值越大就表示差异越显著。
数据格式
2021.4.12_3
绘制
setwd(".../data")#设置目标路径,自己修改library(RColorBrewer)#配色用
df <- read.csv("df.csv",row.names = 1) #导入数据，第一列作为行名
fd <- 0.25 #设置foldchange阈值
cut.fd <- 0.25
pvalue <- 0.05 #设置p阈值
pdf( "df_volcano.pdf") #打开画板
plot(df$fd, -log10(df$P_value_adjust), col="#00000033", pch=19,
xlab=paste("log2 (fold change)"),
ylab="-log10 (P_value_adjust)")
#筛选上下调
up <- subset(df, df$P_value_adjust < pvalue & df$fd > cut.fd)
down <- subset(df, df$P_value_adjust< pvalue & df$fd< as.numeric(cut.fd*(-1)))
#绘制上下调
points(up$fd, -log10(up$P_value_adjust), col=1, bg = brewer.pal(9, "YlOrRd")[6], pch=21, cex=1.5)
points(down$fd, -log10(down$P_value_adjust), col = 1, bg = brewer.pal(11,"RdBu")[9], pch = 21,cex=1.5)
#加上线p、fd阈值线
abline(h=-log10(pvalue),v=c(-1*fd,fd),lty=2,lwd=1)
dev.off()#关闭
注：也可以用ggplot来绘制。
2021.4.12_4
大功告成！
往期内容：
