GO的概念及用R做简单的富集分析






GO的概念及用R做简单的富集分析

小果  生信果  2022-12-18 19:00:17

收录于话题

#绘图#R#代码

GO(Gene Ontology)富集分析基因本体数据库是一个结构化的标准生物学模型,旨在建立基因及其产物知识的标准词汇体系。举个简单的例子,随机在中国挑出10000个人,这些人会有什么特征呢?如果按照男女划分,可以划分两大类,那如果按照省份划分则又可以划分34个大类,那如果按照高矮胖瘦划分呢?你又会划分到哪一类,会像小果一样人见人爱么。




文末有彩蛋


基因的划分也是如此,根据其不同的分类我们将他们划分为三大类:

  • 细胞组分(cellular component)

  • 分子功能(molecular function)

  • 生物学过程(biological process)

Cellular component解释的是基因存在在哪里,在细胞质还是在细胞核?如果存在细胞质那在哪个细胞器上?如果是在线粒体中那是存在线粒体膜上还是在线粒体的基质当中?这些信息都叫Cellular component。


Biological process是在说明该基因参与了哪些生物学过程,比如,它参与了rRNA的加工或参与了DNA的复制,这些信息都叫Biological process


Molecular function在讲该基因在分子层面的功能是什么?它是催化什么反应的?

差异基因GO分析

差异基因GO分析的关键是用统计学方法进行基因富集,分析这些基因参与了何种生物学功能、生物进程以及亚细胞定位,目前常用的基因富集分析法是基于超几何分布,用Fisher精确检验或卡方检验完成。

• 蛋白质或者基因可以通过ID对应或者序列注释的方法找到与之对应的GO号,而GO号可对应到Term,即功能类别或者细胞定位。

• 功能富集分析: 功能富集需要有一个参考数据集,通过该项分析可以找出在统计上显著富集的GO Term。

• GO功能分类是在某一功能层次上统计蛋白或者基因的数目或组成,往往是在GO的第二层次。此外也有研究都挑选一些Term,而后统计直接对应到该Term的基因或蛋白数。结果一般以柱状图或者饼图表示。

• 以差异基因作为前景基因,全部基因作为背景基因(参考基因),找出差异基因相关的GO分类,计算这些差异基因同GO 分类中某(几)个特定的分支的超几何分布关系,GO 分析会对每个有差异基因存在的GO 返回一个p-value,小的p 值表示差异基因在该GO 中出现了富集。

一般取n大于3,校正值(corrected p value)<0.05的条目作为显著性结果

• GO 分析对实验结果有提示的作用,通过差异基因的GO 分析,可以找到富集差异基因的GO分类条目,寻找不同样品的差异基因可能和哪些基因功能的改变有关。

下面的代码是小果请教大佬的,大佬给指了个方向,然后就做出来了,请大佬喝了奶茶。

也有人经常请小果喝奶茶,也管小果叫大佬,哈哈哈哈哈

输入数据

按照条目分类:


代码部分

###读取数据:r1=read.table("bp.txt",sep="t",header=F,as.is=T,quote="!")r2=read.table("cc.txt",sep="t",header=F,as.is=T,quote="!")r3=read.table("mf.txt",sep="t",header=F,as.is=T,quote="!")
##数据处理:x1=nrow(r1); x2=nrow(r2); x3=nrow(r3); x=x1+x2+x3m=c(r1[,2],r2[,2],r3[,2])l=c(r1[,1],r2[,1],r3[,1])y=ceiling(max(m)/10)*15
###画图:pdf(file="goBarplot.pdf",width=15)par(mar=c(20,4,3,3),mgp=c(0.8,0.3,0),cex.axis=.7)barplot(m,beside=T,col=c(rep(rgb(153/255,216/255,201/255),x1),rep(rgb(44/255,127/255,184/255),x2),rep(rgb(201/255,148/255,199/255),x3)),space=0,xaxs='i',yaxs='i',yaxt='n',las=2,names.arg=l,ylab="target genes")abline(h=0)par(xpd=T)lx=max(nchar(l))y1=1/4*y;y2=3/4*ysegments(0,-y1,0,-y2); segments(0,-y2,x,-y2); segments(x1,-y1,x1,-y2); segments(x1+x2,-y1,x1+x2,-y2); segments(x,-y1,x,-y2)text(x1/2,-(y2-1/10*y2),labels='biological_process',pos=1,cex=0.6,col=rgb(153/255,216/255,201/255))text(x1+x2/2,-(y2-1/10*y2),labels='cellular_component',pos=1,cex=0.6,col=rgb(44/255,127/255,184/255))text(x1+x2+x3/2,-(y2-1/10*y2),labels='molecular_function',pos=1,cex=0.6,col=rgb(201/255,148/255,199/255))axis(2)dev.off()

最后的结果图,是不是让人眼前一亮,下期将为你带来更多生信学习技巧。

本期的分享代码编号:f01

文末点击阅读原文或者扫描右侧二维码供邮箱。即可获得3天云服务器免费试用体验。


推荐阅读

扫码关注我们

shengxinguoer

生信果


生信硬核知识解答

和小果一起学生信

点击阅读原文,免费试用服务器