富集分析后的结果如何画chord弦状图

Chord Diagram弦状图-用来画富集分析的结果

前面看的文献中,有一个很漂亮的圈图,https://mp.weixin.qq.com/s/NIqSTvJ916ZGk_dwsRlH5w,就是kegg和go分析后的基因和通路之间的关系图,用cnetplot函数画出来的图也是一样的含义,不过这个圈图看起来就真的很好看,就想画,主要其实颜色是可以选择的,用https://www.sioe.cn/yingyong/yanse-rgb-16/中的颜色就可以选择,所以可以画出来更好看的图。不过后面在数据变换的时候,就还是没转过来可以直接用最后面画图用的数据,就还在思考如何一步步构建包中的中间数据从而得到足后的数据,后来老大告诉我直接构建后面的数据就可以了。

弦状图这种类型的图表可视化了实体之间的相互关系。实体之间的连接用于显示它们共享某些相同的东西。节点沿圆排列,通过使用圆弧将点之间的关系连接在一起。值被分配给每个连接,它由每个弧的大小成比例地表示。颜色可以用来把数据分组到不同的类别,这有助于进行比较和区分组。不过当有太多的连接显示时,可能就会出现混乱。

这使得弦图非常适合于比较数据集内或不同数据组之间的相似性。

使用R语言中的绘图包,比如 circlize, RCircos, CIRCUS, and OmicCircos等来绘制。

一张显示人口流动的图如下,来自http://download.gsb.bund.de/BIB/global_flow/

image-20200617112544546

GOplot包可用于生物数据的可视化。但是要注意该包不能用于执行这些分析,只能把分析结果进行可视化

这个包的用法参考:http://wencke.github.io/。主要是可以把富集分析后的结果进行圈图展示。链接中的讲解代码很长,提炼后的代码如下

david <- EC$david
genelist <- EC$genelist
circ <- circle_dat(EC$david, EC$genelist)
chord <- chord_dat(data = circ, genes = EC$genes, process = EC$process)
GOChord(chord, space = 0.02, gene.order = 'logFC', gene.space = 0.25, gene.size = 5)

画出的图如下

image-20200617120210980

EC是这个GOplot包中的实例数据,为了能用这个包中的作图函数,就先了解这个包中的数据。那么david和genelist的数据长下面这样。

image-20200617114924282

生成的中间数据chord长下面这样。

image-20200617115516739

最后用GOChord函数画图的数据chord的数据格式如下图所示。这是我们要把自己的数据变成这样的数据的样子。就是下面的chord数据格式。变成下面这样就可以画圈图了。

image-20200617114629848

上面最后画图的输入数据函数就是chord,最后的logFC是从前面的差异分析的结果得来的,也就是genelist中的logFC。然后chord数据除了logFC以外,前面的数据就是如果某个基因在通路中,就是1,否则就是0。这中数据就是最前面代码中的chord_dat函数实现的,但是问题是这个chord 数据是从circ得到的,而circ是从david和genelist得到的,david数据是富集分析的结果,genelist是差异分析的结果,现在我没有差异分析的结果,我只是有前面比如通过cox回归筛选得到的基因名,怎么办呢?

问题解决的办法就是,直接制作数据格式到chord这样的数据,把logFC变为0就可以了。

  • 下面的代码是从前面cox回归筛选到的基因开始的
load(file = 'sur_logKM_cox.Rdata')
cox=rownames(cox_results[cox_results[,4]<0.01,])
log=names(log_rank_p[log_rank_p<0.01])
surGenes=intersect(cox,log)

上面得到的surGenes如下:就是一堆字符串啦

image-20200617121812423

但是用老大出品的这个AnnoProbe进行ID转换,得到了如下图所示非常方便后续进行操作的表格


library(AnnoProbe) #用老大出品的这个AnnoProbe进行ID转换,
surGenes=annoGene(surGenes,ID_type = 'SYMBOL','human') 


tail(sort(table(surGenes$biotypes)))
head(surGenes)

image-20200617122118975

接下来富集分析

library(ggplot2)
library(clusterProfiler)
library(org.Hs.eg.db)
df <- bitr(unique(surGenes$SYMBOL), fromType = "SYMBOL",
           toType = c( "ENTREZID"),
           OrgDb = org.Hs.eg.db)#富集分析要求数据数据的是ENTREZID,要从前面得到的SYMBOL用bitr函数得到ENTREZID
head(df)   

gene_up=df$ENTREZID

#下面是为了画图做的插曲
#将bp的结果保存
bp_go <- enrichGO(gene_up, 
                  OrgDb = "org.Hs.eg.db", 
                  ont="BP",
                  pvalueCutoff  = 0.001,
                  qvalueCutoff  = 0.001) 

enrichgo=DOSE::setReadable(bp_go, OrgDb='org.Hs.eg.db',keyType='ENTREZID')#setReadable函数直接进行ID转换了,从ENTREZID再得到symbol

mm <- enrichgo@result[1:5,c(2,7,8)]

得到的mm数据如下(mm是我名字的缩写,哈哈)

image-20200617122407370
tmp=do.call(rbind,
        apply(mm, 1,function(x){
          data.frame(go=x[1],
                     gene=strsplit(x[3],'/')[[1]])
        })
)

得到的tmp数据格式如下

image-20200617123212419

接下来就是将这个长形数据用reshape包中的dcast函数变为宽形数据,如下

library(reshape2)
tmp2=dcast(tmp,go~gene)
tmp2[is.na(tmp2)]=0
rownames(tmp2)=tmp2[,1]
tmp2=tmp2[,-1]
tmp2=t(tmp2)
tmp2[tmp2!=0]=1
tmp2=as.data.frame(tmp2)
tmp2$logFC=0
cg=rownames(tmp2)
tmp2=apply(tmp2,2,as.numeric)
rownames(tmp2)=cg

得到的tmp2数据如下,这就是已经做好了和前面的实例数据chord长得一摸一样的数据了。

image-20200617123041839

接下来画图

GOChord(tmp2, space = 0.02, gene.order = 'logFC', gene.space = 0.25, gene.size = 5)
image-20200617123440677

区别是我们没有用变化倍数进行聚类,所以没有logFC的渐变的变化,不过已经超级好看了!

如果想有倍数变化的结果,那么就是一开始就是有差异分析的结果,就可以下面的代码,得到的就是本文第二张图片了。

# Create the plot
chord <- chord_dat(data = circ, genes = EC$genes, process = EC$process)
GOChord(chord, space = 0.02, gene.order = 'logFC', gene.space = 0.25, gene.size = 5)

参考的代码全套在http://wencke.github.io/,也可以搜索到翻译过来的中文版的教程。

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