一个看似简单的需求——修改富集分析的dotplot图
刘小泽写于2020.2.6
最近再一次做起了转录组，但这一次需求有点小改变，需要自己定制一下，具体原因看本文吧。其中要特别表扬花花💏同学，帮了个大忙！
问题由来
我们一般进行富集分析，一般的做法都是：
# 先得到上调的差异基因gs_up，可以先不进行pvalue、qvalue的过滤
up <- enrichGO(gene = gs_up,
keyType = 'ENSEMBL',
OrgDb = org.Hs.eg.db,
ont = "BP",
pAdjustMethod = "BH",
pvalueCutoff = 1,
qvalueCutoff = 1,
readable = TRUE)
dotplot(up, showCategory=10,title="EnrichmentGO_BP")
得到的结果默认是：
现在，想把横坐标进行改变，改成：Fold Enrichment，像这种图
什么是Fold Enrichment？和原来的GeneRatio啥关系？
看这里：https://yulab-smu.github.io/clusterProfiler-book/chapter13.html
实际上，Fold Enrichment也是一个比值，是在原来GeneRatio的基础上，再除以BgRatio。
如下图，以第一行为例，默认情况下，dotplot给出的横坐标上，它会显示：7/797=0.0087
如果要改成Fold Enrichment，就是：(7/797)/(20/20610)=90.51
那么如何改变横坐标？
感觉很简单的样子，只要自己计算出来Fold Enrichment，然后加到原来数据中就好了。但是dotplot不太容易定制，它只支持GeneRatio or Count
看说明书：https://bioconductor.org/packages/devel/bioc/manuals/enrichplot/man/enrichplot.pdf
不过，好在Y叔的函数都是基于ggplot2的，于是我们可以自己拿到数据去作图
尝试第一步 准备作图数据
# 之前计算的结果
up <- enrichGO(gene = gs_up,
keyType = 'ENSEMBL',
OrgDb = org.Hs.eg.db,
ont = "BP",
pAdjustMethod = "BH",
pvalueCutoff = 1,
qvalueCutoff = 1,
readable = TRUE)
# 根据显著性取前10个
test=as.data.frame(up)
test=test[1:10,]
# 然后自己计算Fold Enrichment，并按照Fold Enrichment升序排序
library(stringr)
gr1 <- as.numeric(str_split(test$GeneRatio,"/",simplify = T)[,1])
gr2 <- as.numeric(str_split(test$GeneRatio,"/",simplify = T)[,2])
bg1 <- as.numeric(str_split(test$BgRatio,"/",simplify = T)[,1])
bg2 <- as.numeric(str_split(test$BgRatio,"/",simplify = T)[,2])
test$fold <- (gr1/gr2)/(bg1/bg2)
test <- arrange(test,fold)
尝试第二步 作图
使用最简单的作图代码
ggplot(test,aes(x = fold,y = Description))+
geom_point(aes(color = p.adjust,
size = Count))+
scale_color_gradient(low = "red", high = "blue")+
xlab("Fold Enrichment")+
theme_bw()
结果出来是这样：
尝试修改
虽然结果出来了，但总感觉怪怪的，没有按照一个顺序排列。
这里花花给出了最大的帮助，她一眼（可能是两眼）就看到了，虽然X轴是从小到大排列，但Y轴不是，它是按照首字母的顺序进行的排列，所以结果是乱的
于是她帮我加了一行代码，将Y轴的Description设成有序型因子变量，
test$Description = factor(test$Description,levels = test$Description,ordered = T)
实际上将我已经按照Fold Enrichment排序好的数据框保持原样，并且其中的Description同样按照从上到下的顺序进行了排列
比如原来的Description是字符型变量，作图默认按照字母顺序，也就是carboxylic这一条会是第一个，但现在按照Fold Enrichment排序后，它到了第3位。我们就顺势把这个carboxylic设成有序型因子，它现在就不再是一个字符串了，而是按照一个数3进行后面的排列【这里看不懂也没啥关系，反之因子型的概念就是有点绕】：
这个修改完之后，再次运行作图代码：
完事了吗？好像还差一点！
好像图例中这个p.adjust的顺序有点...强迫症犯了，最好是从小到大排列
于是，搜索了一下，找到了答案(https://stackoverflow.com/questions/22458970/how-to-reverse-legend-labels-and-color-so-high-value-starts-downstairs)：
ggplot(test,aes(x = fold,y = Description))+
geom_point(aes(color = p.adjust,
size = Count))+
scale_color_gradient(low = "red", high = "blue")+
xlab("Fold Enrichment")+
theme_bw()+
#edit legends
guides(
#reverse color order (higher value on top)
color = guide_colorbar(reverse = TRUE))
#reverse size order (higher diameter on top)
#size = guide_legend(reverse = TRUE))
嗯，这回感觉清爽多了
欢迎关注我们的公众号～_～
我们是两个农转生信的小硕，打造生信星球，想让它成为一个不拽术语、通俗易懂的生信知识平台。需要帮助或提出意见请后台留言或发送邮件到jieandze1314@gmail.com
Welcome to our bioinfoplanet!
