统计多个基因列表交集并画韦恩图和Upset图

有时候,我们从GEO数据库挖掘了好几组数据集的差异表达基因,我们需要统计这些基因列表的交集,一般可以用韦恩图和Upset图。而韦恩图一般适合小样本的交集,而多个样本的话可以用Upset图来表示。

比如,我们通过GEO2R工具,获得了差异基因列表的结果,假设有6个基因列表,我们可以把显著性差异表达基因提取处理,用excel制作一个genelist,表头是数据集,每列是基因,保存为list.csv,导入进来以后,见表 [1]

list<-read.csv("/list.csv")
head(list)

Table 1: 6个基因列表

A B C D E F
SLC2A3 TMC4 CDC25B C15orf48 SYTL1 KIF2C
UNC5B CTHRC1 C1orf132 E2F1 CDCA5 AURKA
MT2A MMD ST3GAL5 C1orf132 BIRC5 AURKA
ST3GAL5 VSIG2 TMEM154 MTHFD2 MCM10 CDC20
TUBB6 CTHRC1 ST3GAL5 KDELR3 MTHFD2 C17orf53
CTHRC1 MMP11 TMEM154 NSG1 MELK CDCA5

但是这个表其实是有问题的,那就是这些列表里有些基因是重复的,而且还有NA值,我们可以把各自的重复值和NA值删除。记住,一定是删除各列的重复值,而不是全部的表格。

我们可以把列表单独提取处理,然后分别去除重复值和NA值,然后指定成list的形式(当然表格形式也是可以的,不过不同行数的list怎么合并成一个表格我不会,除非导出来用excel)。

## 单独提出各个list
listA<-list$A
listB<-list$B
listC<-list$C
listD<-list$D
listE<-list$E
listF<-list$F
## 去除各自列表的重复值
listA<-listA[!duplicated(listA)]
listB<-listB[!duplicated(listB)]
listC<-listC[!duplicated(listC)]
listD<-listD[!duplicated(listD)]
listE<-listE[!duplicated(listE)]
listF<-listF[!duplicated(listF)]
## 去除各自列表的NA值
listA<-na.omit(listA)
listB<-na.omit(listB)
listC<-na.omit(listC)
listD<-na.omit(listD)
listE<-na.omit(listE)
listF<-na.omit(listF)
## 合并list
list<-list(A=listA,B=listB,C=listC,D=listD,E=listE,F=listF)

不过,经过我后期验证,其实不去除重复值和NA值,最后的结果是一样的,可能已经自动去除了,具体是否这样我没去深入研究。


1 绘制韦恩图

一般而言,主流画韦恩图的包是VennDiagram,教程非常多,比如使用VennDiagram包绘制韦恩图

但是这里我想推荐一下ggVennDiagram这个包,有一天在Y数的公众号瞎溜达,看到一篇的ggVennDiagram 诞生记文章,知道这个包可以用ggplot2的语句,而且不限基因列表限制。

我以前喜欢用联川生物的在线工具画韦恩图,我承认他家可以画很漂亮的韦恩图,并且还是在线交互式绘图,我看过很多在线工具,但是他家的确实是最好看的。而这个工具利用的就是VennDiagram这个包,但是缺点是最多只能画5个,所以还是推荐ggVennDiagram这个包,可能颜值差了点,但是效果还是不错的。

先看最基本的函数,见图 [1]所示。

library(ggVennDiagram)
ggVennDiagram(list)
Figure 1: 6个基因列表的韦恩图

可以看到默认显示的是交集数字和比例,然后按照集合的数字填充颜色,对于很多个交集的话,这样的图确实很难看,可是这个包符合ggolot2语句,所以是可以高度DIY的。我们看一下help的函数

ggVennDiagram(
  x,
  category.names = names(x), #自定义数据集的名称
  show_intersect = FALSE, #是否交互式演示
  set_color = "black", # 数据集的颜色,默认即可
  set_size = NA, #数据集标签的大小,默认即可
  label = c("both", "count", "percent", "none"), # 显示数值、比例还是不显示
  label_alpha = 0.5, #各个标签集合的外框透明度,0为全透明,默认0.5
  label_geom = c("label", "text"), #显示数字和标签框
  label_color = "black",#标签颜色
  label_size = NA, #标签字体大小
  label_percent_digit = 0, #标签的百分位数
  label_txtWidth = 40, #标签文本的宽度
  edge_lty = "solid", #标签边缘类型,默认实心
  edge_size = 1, #标签边缘大小
  ...
)

可以看到可支持的函数很多,作者也写了一些引导说明,大家可以自己尝试跑一些。

比如我们只显示数字,设置为虚线显示,见图 [2]所示。

ggVennDiagram(
  list,
  label = "count",
  edge_lty = "dashed",
  edge_size = 0.5
)
Figure 2: 韦恩图

然而,ggVennDiagram是支持ggplot2语句的,也就是说我们可以继续叠加修改,比如换个颜色,加个主题,加个标题什么的,还有把那个count的标签去掉什么的,最终效果见图 [3]所示,也可以换个配色,见图[4]所示。

library(ggplot2)
ggVennDiagram(
  list,
  label = "count",
  edge_lty = "dashed",
  edge_size = 0.5)+
  scale_fill_gradient(low="steelblue",high = "brown")+
  theme_bw()+
  ggtitle("Veen plot for six sets")+
  theme(legend.position = 'none')
Figure 3: 韦恩图
ggVennDiagram(
  list,
  label = "count",
  edge_lty = "dashed",
  edge_size = 0.5)+
  scale_fill_distiller(palette = "RdBu")+
  theme_void()+
  ggtitle("Veen plot for six sets")+
  theme(legend.position = 'none')
Figure 4: 韦恩图

可是这依然不够ggplot2,更多的可以看ggVennDiagram 的新生这篇文章,我们可以把数据转换,比如给个华丽丽的标签,见图 [5]所示。

venn = Venn(list)
data = process_data(venn)
ggplot() +
  geom_sf(aes(fill = count), data = venn_region(data)) +
  geom_sf(aes(color = id), data = venn_setedge(data), show.legend = FALSE) +
  geom_sf_text(aes(label = name), data = venn_setlabel(data)) +
  geom_sf_label(aes(label = count), data = venn_region(data)) +
  theme_void()+
    theme(legend.position = 'none')+
   ggtitle("Veen plot for six sets")
Figure 5: 韦恩图

而我们喜欢按数据集填充颜色,把fill改成id即可,最后效果见图 [6]所示

ggplot() +
  geom_sf(aes(fill = id), data = venn_region(data)) +
  geom_sf(aes(color = id), data = venn_setedge(data), show.legend = FALSE) +
  geom_sf_text(aes(label = name), data = venn_setlabel(data)) +
  geom_sf_label(aes(label = count), data = venn_region(data)) +
  theme_void()+
    theme(legend.position = 'none')+
   ggtitle("Veen plot for six sets")
Figure 6: 韦恩图

2 绘制Upset图

当数据集多于5组的适合,画韦恩图就吃力了,因为圈圈太多,所以很不好看。这个时候就可以用Upset这种图来表示了。而推荐用UpsetR这个包,目前在CRAN上也收录了,很方便的可以安装。

install.package('UpSetR)

UpsetR这个包默认的画图其实不是list的形式,而是各个交集以0和1的expression形式,用来表示有还是无。这个时候就需要需要我们先把数据处理好,这个可以使用UpSetR自带的fromList()函数实现,很简单的一句话,最后的表格见表 [2]

library(UpSetR)
data2<-fromList(list)

Table 2: 6个基因列表的集合

A B C D E F
1 1 0 1 0 0
1 0 0 0 0 0
1 1 1 1 1 0
1 1 1 1 1 1
1 0 0 0 0 0
1 1 1 1 1 0

当然我们也可以不转换数据,只要加上一句fromList(list)也可以画图,现在我们用最基本的函数进行下Upset图的绘制,效果见图 [7]所示。

upset(data2) # 要先加载UpSetR
Figure 7: 数据集的Upset图
### 如果不想对list进行转换,也可以直接用下面这个代码
# upset(fromList(list))

很明显的,这里有个很大的bug,那就是明明有6个数据集,那就是默认只显示了5个数据集(nset=5),既然默认的参数统计不出来,我们就手动添加进来。

另外这个函数的功能相当的多,可以各种排序(按频率排序,按频数排序,还可以按数据集分组),还可以给点填色等等。具体的可以在示例说明找到,一切都等你发现

比如我们按freq排序,设置6个数据集,见图 [8]所示。

upset(fromList(list),nsets = 6,order.by = "freq")
Figure 8: 6个数据列表的Upset图
## upset(fromList(list),sets = c("A","B","C","D","E","F")) #也可以按需设置

我们想知道6个数据集都交集的情况,还想标个颜色,可以这样,见图 [9]所示。

upset(fromList(list),nsets = 6,order.by = "freq",,queries = list(list(query = intersects, params = list("A","B","C","D","E","F"), active = T)))
Figure 9: 6个数据列表的Upset图,标记全部交集

UpSetR是一个很强大的包,功能很多,但有一个问题就是不能使用ggplot语句,导出的图也不是ggplot,这个我们可以用Y叔的ggplotify转换。当然其实还有一个叫ggupset的包可以使用ggplot语句的,但是需要数据处理,语法很复杂,有兴趣自己百度学习。


3 两组图的结合

Upset图有个缺点就是,如果按频数画图,右上角会空出很大一块,我们能不能把韦恩图和Upset图组合到一起,这个问题Y叔已经考虑并解决过了,可以参考转UpSet图为ggplot?

这里用到了Y叔叔的yyplot这个神包,然而这个包也是很难装上的,因为这是Y叔的私人包,需要补装很多包,而很多包都是只有Github上才能安装,比如gglayer这个包,你在CRAN和BiocManager上都上找不到的,只有Y叔的Github上的包才能安装,但是又经常访问困难,所以我Fork到了我的gitee上。而且用yyplot画韦恩图,需要用ggvenn这个函数,但是这个函数又要配置Java环境,所以个中曲折自己摸索吧。

remotes::install_git('https://gitee.com/swcyo/gglayer/'))
remotes::install_git('https://gitee.com/swcyo/yyplot/'))
library(yyplot)

我们可以不用Y叔的yyplot画韦恩图(因为很难安装),用ggVennDiagram的图即可,两种图的组合原理是两个ggplot图形的叠加,这里要用到Y叔的ggimge这个包(Y叔是无处不在的),最后效果见图 [10]所示。

library(ggplotify) #把别的图转为ggplot2
library(ggimage) # 组合图片
p1<-upset(fromList(list),nsets = 6,order.by = "freq",,queries = list(list(query = intersects, params = list("A","B","C","D","E","F"), active = T)))
g1<-as.ggplot(p1) # 转换为ggplot2图片
g2<-ggplot() +
  geom_sf(aes(fill = id), data = venn_region(data)) +
  geom_sf(aes(color = id), data = venn_setedge(data), show.legend = FALSE) +
  geom_sf_text(aes(label = name), data = venn_setlabel(data)) +
  geom_sf_label(aes(label = count), data = venn_region(data)) +
  theme_void()+
    theme(legend.position = 'none')
g3<-g1 + geom_subview(subview = g2 , x=.8, y=.7, w=.45, h=.45)
g3
Figure 10: 6个数据列表的Upset图,标记全部交集

但是这里有一个bug,就是左侧的bar不见了,不知道是不是我电脑的问题。。。

如果要用yyplot画韦恩图的画,可以这样,最后效果见 [11]

library(yyplot)
g4<-ggvenn(fromList(list))
g5<-g1 + geom_subview(subview = g4 + theme_void(), x=.8, y=.7, w=.5, h=.5)
g5
Figure 11: 6个数据列表的Upset图,标记全部交集

4 在线交互式作图

有很多网站可以在线交互式作图,有国产开发的,也有shiny,比如:

https://www.omicstudio.cn/tool/6

https://www.hiplot.com.cn

gehlenborglab.shinyapps.io/upsetr/

https://asntech.shinyapps.io/intervene/

更多的在线功能等着大家发现。

最后编辑于
©著作权归作者所有,转载或内容合作请联系作者
  • 序言:七十年代末,一起剥皮案震惊了整个滨河市,随后出现的几起案子,更是在滨河造成了极大的恐慌,老刑警刘岩,带你破解...
    沈念sama阅读 157,298评论 4 360
  • 序言:滨河连续发生了三起死亡事件,死亡现场离奇诡异,居然都是意外死亡,警方通过查阅死者的电脑和手机,发现死者居然都...
    沈念sama阅读 66,701评论 1 290
  • 文/潘晓璐 我一进店门,熙熙楼的掌柜王于贵愁眉苦脸地迎上来,“玉大人,你说我怎么就摊上这事。” “怎么了?”我有些...
    开封第一讲书人阅读 107,078评论 0 237
  • 文/不坏的土叔 我叫张陵,是天一观的道长。 经常有香客问我,道长,这世上最难降的妖魔是什么? 我笑而不...
    开封第一讲书人阅读 43,687评论 0 202
  • 正文 为了忘掉前任,我火速办了婚礼,结果婚礼上,老公的妹妹穿的比我还像新娘。我一直安慰自己,他们只是感情好,可当我...
    茶点故事阅读 52,018评论 3 286
  • 文/花漫 我一把揭开白布。 她就那样静静地躺着,像睡着了一般。 火红的嫁衣衬着肌肤如雪。 梳的纹丝不乱的头发上,一...
    开封第一讲书人阅读 40,410评论 1 211
  • 那天,我揣着相机与录音,去河边找鬼。 笑死,一个胖子当着我的面吹牛,可吹牛的内容都是我干的。 我是一名探鬼主播,决...
    沈念sama阅读 31,729评论 2 310
  • 文/苍兰香墨 我猛地睁开眼,长吁一口气:“原来是场噩梦啊……” “哼!你这毒妇竟也来了?” 一声冷哼从身侧响起,我...
    开封第一讲书人阅读 30,412评论 0 194
  • 序言:老挝万荣一对情侣失踪,失踪者是张志新(化名)和其女友刘颖,没想到半个月后,有当地人在树林里发现了一具尸体,经...
    沈念sama阅读 34,124评论 1 239
  • 正文 独居荒郊野岭守林人离奇死亡,尸身上长有42处带血的脓包…… 初始之章·张勋 以下内容为张勋视角 年9月15日...
    茶点故事阅读 30,379评论 2 242
  • 正文 我和宋清朗相恋三年,在试婚纱的时候发现自己被绿了。 大学时的朋友给我发了我未婚夫和他白月光在一起吃饭的照片。...
    茶点故事阅读 31,903评论 1 257
  • 序言:一个原本活蹦乱跳的男人离奇死亡,死状恐怖,灵堂内的尸体忽然破棺而出,到底是诈尸还是另有隐情,我是刑警宁泽,带...
    沈念sama阅读 28,268评论 2 251
  • 正文 年R本政府宣布,位于F岛的核电站,受9级特大地震影响,放射性物质发生泄漏。R本人自食恶果不足惜,却给世界环境...
    茶点故事阅读 32,894评论 3 233
  • 文/蒙蒙 一、第九天 我趴在偏房一处隐蔽的房顶上张望。 院中可真热闹,春花似锦、人声如沸。这庄子的主人今日做“春日...
    开封第一讲书人阅读 26,014评论 0 8
  • 文/苍兰香墨 我抬头看了看天上的太阳。三九已至,却和暖如春,着一层夹袄步出监牢的瞬间,已是汗流浃背。 一阵脚步声响...
    开封第一讲书人阅读 26,770评论 0 192
  • 我被黑心中介骗来泰国打工, 没想到刚下飞机就差点儿被人妖公主榨干…… 1. 我叫王不留,地道东北人。 一个月前我还...
    沈念sama阅读 35,435评论 2 269
  • 正文 我出身青楼,却偏偏与公主长得像,于是被迫代替她去往敌国和亲。 传闻我的和亲对象是个残疾皇子,可洞房花烛夜当晚...
    茶点故事阅读 35,312评论 2 260

推荐阅读更多精彩内容