富集分析后的结果如何画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-20200617112544546GOplot包可用于生物数据的可视化。但是要注意该包不能用于执行这些分析,只能把分析结果进行可视化。
这个包的用法参考: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-20200617120210980EC是这个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-20200617122407370tmp=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/,也可以搜索到翻译过来的中文版的教程。