RNA-seq

差异分析后画火山图

2021-12-02  本文已影响0人  LELELEANNA

> library(ggplot2)

> library(ggrepel)

> res = na.omit(res)#去NA

> res$type = ifelse(res$log2FoldChange >= 1.5 & res$padj <=0.001,"up_regulated", ifelse(res$log2FoldChange<= -1.5 & res$padj <=0.001,"down_regulated","none siginificant"))

> table(res$type)#列出上调、下调基因的数目

> title = paste0('The number of up gene is ', nrow(res[res$group == 'up_regulated',]),'\nThe number of down gene is ', nrow(res[res$group == 'down_regulated',]))

> res = as.data.frame(res)#不转换会出现报错:`data` must be a data frame

> ggplot(res, aes(x = log2FoldChange, y = -log10(padj)))+geom_point(aes(color = group)) + scale_color_manual(values = c("purple","grey","navy blue"),limits = c('up_regulated','none siginificant'," down_regulated")) + theme_bw(base_size = 20) + theme(plot.title = element_text(size=15,hjust = 0.5) + theme_classic())+ggtitle(title)

上述代码可画出最简单的火山图,如果想加点花里胡哨的东西↓

参考:TCGA数据分析系列之火山图 (qq.com)

> p = ggplot(res, aes(log2FoldChange, -log10(padj),color  = type))

> p +geom_point()

> x_lim = max(res$log2FoldChange, - res$log2FoldChange)#设置x轴坐标

> gg = p + geom_point(size =1) + xlim(-x_lim, x_lim) +labs(x= "log2FoldChange", y = "-log10(padj)") +scale_color_manual(values = c("#A52A2A","grey","#f8766d")) + geom_hline(aes(yintercept=-1*log10(0.05)),colour="black", linetype="dashed") +     geom_vline(xintercept=c(-2,2),colour="black", linetype="dashed")

> print(gg)

上一篇下一篇

猜你喜欢

热点阅读