跟着CELL学作图之火山图
“实践是检验真理的唯一标准。” “复现是学习R语言的最好办法。”
DOI: 10.1016/j.cell.2020.05.032
这篇2020年发表在cell上关于新冠的组学文章里面有大量的生信内容。今天带大家复现其中的一个Supplemental Figure：火山图。
本文代码及示例数据领取：后台回复“210412” （这个是我转载的，需要到 木舟笔记这个公众号留言获取示例数据和代码）
这图确实比一般的火山图美观且简洁。
火山图的意义
火山图可用于展示两组样本间基因表达水平差异的分布状况。
横轴log2 fold change差异表达倍数（Fold Change值，简称FC），差异越大的基因分布X轴在两端。
纵坐标用-log10 p-value表示，对P值进行-log10的转化。转化后，值越大就表示差异越显著。
数据格式
绘制
代码语言：javascript复制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来绘制。
大功告成！
