搜索
查看: 3746|回复: 3

[annotation] GO.db/topGO/goseq R包

[复制链接]

29

主题

131

帖子

1178

积分

金牌会员

Rank: 6Rank: 6

积分
1178
发表于 2017-8-19 22:57:46 | 显示全部楼层 |阅读模式
本帖最后由 anlan 于 2017-8-19 23:07 编辑

GO.db
GO.db是一个用注释maps来描述Gene Ontology的一个R包,其在很多GO注释及富集的R包中被调用,应用广泛,每半年更新一次。主要用来描述各个Term之间的父子节点联系以及Term的信息。

至于使用的话,可以用names(GOMAPCOUNTS)来查看主要有几个函数,以BP为例:
所有父节点信息 GOBPANCESTOR
上一级父节点信息 GOBPPARENTS
所有子节点信息 GOBPCHILDREN
下一级子节点信息 GOBPOFFSPRING
每个Term所有描述信息 GOTERM

使用方式个人总结下有以下几种:
转化为list
     ancestor <- as.list(GOBPANCESTOR)
转化为data.frame
     ancestor <- toTable(GOBPANCESTOR)
直接查看某个GO id的信息
     ancestor <- GOBPANCESTOR$"GO:0001580"
用处:一些R包会调用GO.db获得某个GO id的描述信息,还有些会用GO.db(例如topGO)的每个GO之间的关系来绘制 directed acyclic graph (DAG 有向无环图),而我使用这个R包是想将几个子节点的父节点也参与GO富集(通常的GO富集R包都不会有这个功能)。

topGO
topGO这个R包可以用来做GO富集,而我接触到这个R包则是因为它能将富集结果以有向无环图形式展示,这一点是其他大部分做GO富集的R包所有没有的,就算是有这个功能的R包一般也调用了topGO。

但从富集这个功能来说,topGO比其他的R包没有较大的优势,以及除了拥有不止一种的统计学方法以及算法来计算富集的P值,一般的步骤如下(一般用于有模式物种或者的有bioconductor注释包的物种):
[AppleScript] 纯文本查看 复制代码
sampleGOdata <- new("topGOdata",     description = "Simple session",,     allGenes = geneList, geneSel = topDiffGenes,     nodeSize = 10,     annot = annFUN.db, affyLib = affyLib)

annFUN.db表示从bioconductor注释包中提取注释信息的函数(其实还有其他函数可以使用),affyLib则是对应注释包

然后就是使用经典的fisher检验进行富集分析
[AppleScript] 纯文本查看 复制代码
resultFisher <- runTest(sampleGOdata, algorithm = "classic", statistic = "fisher")
当然topGO这个R包并不是只能做有注释R包的物种,它提供了一种自定义的注释,其实就是输入你自己的GO注释,然后构建topGOdata对象,接着就跟上述模式物种富集方法类似了,从而可以使用其他功能。我们需要准备一个文件go_bp.txt(以BP为例),第一列是gene id,第二列则是gene对应的GO id,每个GO id以逗号分隔,然后使用readMappings读入,然后再输入差异基因列表geneid.txt
[AppleScript] 纯文本查看 复制代码
geneid2GO <- readMappings(file = "go_bp.txt")
        genenames <- names(geneid2GO)
        gene_data <- read.table(file = "geneid.txt")
        gene_id <- gene_data[,1]
        genelist <- factor(as.integer(genenames %in% gene_id))
        names(genelist) <- genenames
        
        GOdata <- new("topGOdata", allGenes = genelist,, 
              annot = annFUN.gene2GO, gene2GO = geneid2GO)
这个也算一个不错的功能,但是这个功能有一些GO富集的R包也有。。。。。。

最后则是算是这个R包特有的将GO富集结果以DAG(有向无环图)展示
[AppleScript] 纯文本查看 复制代码
showSigOfNodes(GOdata, pvalue, firstSigNodes = 10, useInfo = 'all')

goseq
goseq也是一个做GO富集的R包(也做KEGG富集,但是KEGG的数据包已经在R里面好久没更新了,因为KEGG数据库收费了)。

goseq相比其他富集R包,有个特别的地方:其富集算法使用了Wallenius non-central hyper-geometric distribution,相对于普通的 Hyper-geometric distribution,此分布的特点是从某个类别中抽取个体的概率与从某个类别之外抽取一个个体的概率是不同的,这种概率的不同是通过对基因长度的偏好性进行估计得到的,从而其认为能更为准确地计算出 GOterm被差异基因富集的概率。

因此在做富集前,需要对gene长度进行计算,如果有参考基因组的话,直接选择这个R包所提供的基因组注释信息来提取gene长度
[AppleScript] 纯文本查看 复制代码
pwf=nullp(genes,"hg19","ensGene")


然后进行富集
[AppleScript] 纯文本查看 复制代码
GO.wall=goseq(pwf,"hg19","ensGene")

当然goseq也能像topGO一样接受自定义的注释信息,唯一的不同都是需要自行统计gene的长度,然后再导入,比如差异基因gene_id,自行提供的基因对应GO信息gene2GO
[AppleScript] 纯文本查看 复制代码
pwf=nullp(genes,bias.data=gene_lengths)

gene_lengths是每个基因的长度

然后进行自定义GO注释的富集分析
[AppleScript] 纯文本查看 复制代码
go <- goseq(pwf = pwf, id = pro_id, gene2cat = uniprotid2GO, method = "Hypergeometric", use_genes_without_cat = TRUE)
Summary
上述的3个R包,各有特点吧,如果没有特别要求的话,后两者都能满足GO富集的需求。GO.db包是很多GO分析包的基础,了解下也不错。

最近看了下clusterProfiler包代码,发现其也是调用了GO.db包进行一些GO信息的提取,并也调用了topGO包来展现富集分析的结果,而且也支持自定义的GO注释信息的输入来进行GO富集分析,所以说也蛮全面的,而且是使用起来也比较友好。

其实最好的说明还是查看三者的说明文档,里面都有详细的教程,看一下一般就懂了
[size=0em]​



回复

使用道具 举报

2

主题

4

帖子

63

积分

注册会员

Rank: 2

积分
63
发表于 2018-6-4 11:15:35 | 显示全部楼层
你好,我想问下,我是自行下载的基因及其对应的长度,还有注释信息,但是不知道怎么加载到命令行里
,您能帮我解答一下吗
回复 支持 反对

使用道具 举报

2

主题

4

帖子

63

积分

注册会员

Rank: 2

积分
63
发表于 2018-6-4 11:20:59 | 显示全部楼层
我的第一个是csv文件,第一列是基因名,第二列是对应的基因长度-bp
第二个文件是直接下载好的注释信息
回复 支持 反对

使用道具 举报

29

主题

131

帖子

1178

积分

金牌会员

Rank: 6Rank: 6

积分
1178
 楼主| 发表于 2018-6-7 22:32:06 | 显示全部楼层
zhang1012 发表于 2018-6-4 11:20
我的第一个是csv文件,第一列是基因名,第二列是对应的基因长度-bp
第二个文件是直接下载好的注释信息
...

长度的话,直接写入我上面帖子的pwf里就行,GO注释信息的话,则需要自行处理下格式,具体格式可以看这个包的说明文档,其有列出示例的
回复 支持 反对

使用道具 举报

您需要登录后才可以回帖 登录 | 立即注册

本版积分规则

QQ|手机版|小黑屋|生信技能树    

GMT+8, 2018-6-25 15:31 , Processed in 0.116538 second(s), 23 queries .

Powered by Discuz! X3.2

© 2001-2013 Comsenz Inc.