goseq

goseq是一个R包,用于寻找GO terms,即基因富集分析。

GO terms是标准化描述基因或基因产物的词汇,包括三方面,cellular component,molecular funciton,biological process。

每个GO term都有一个GO ID,比如 GO:006260,每个GO term背后都有一系列的相关基因。

GO分析的目的:在差异性基因分析后,我们可能得到很多差异基因,这些基因里的一部分可能跟某个生物过程相关,或几个生物过程相关。经过GO分析后,我们就能将差异性基因具体的生物功能展示出来,为下一步研究做准备。

GOseq需要输入的文件:

1.所有有count的genes。

2.差异性表达的genes。

3.genome信息,基因长度信息。#对于许多模式基因组来说,这些内容都被做成了独立的R包。

4.GO terms包。

>source("http://bioconductor.org/biocLite.R")
>biocLite("goseq")
>biocLite("geneLenDataBase")#genome,genes信息
>biocLite("org.Dm.eg.db")#果蝇的GO terms

>library("goseq")
>library("geneLenDataBase")
>library("org.Dm.eg.db")

>DEG<-read.table("DEG",header=FALSE)
>ALL<-read.table("ALL",header=FALSE)
#DEG:差异性基因表 ALL:所有基因表(数据框格式)

>DEG.vector<-c(t(DEG))
>ALL.vector<-c(t(ALL))
#把数据格式转化为vector,便于下步操作

>gene.vector=as.integer(ALL.vector%in%DEG.vector)
#生成二进制的gene vector(1代表差异性基因,0代表非差异性基因)
>names(gene.vector)<-ALL.vector

>pwf=nullp(gene.vector,"dm3","ensGene")
#生成probability weighting function.

>GO.wall=goseq(pwf,"dm3","ensGene")
#生成GO terms ID 。这边的疑问:genes 没有mapping 到GO categories。 goseq函数有一个选项:gene2cat,如果gene2cat=NULL,则goseq会自动调用getgo函数实现mapping功能,并将输出值gene2cat。

>enriched.GO=GO.wall$category[GO.wall$over_represented_pvalue<.05]
#生成差异性 GO terms ID

>library(GO.db)
>capture.output(for(go in enriched.GO[1:length(enriched.GO)]){
print(GOTERM[go])
cat("___________")
}
,file="SigGo.txt")
#生成具体的GO TERM详解

  

时间: 2024-11-08 19:51:15

goseq的相关文章

GO分析-GOseq的使用教程

GOseq的介绍 GOseq是一个R包,用于寻找GO terms,即基因富集分析.此方法基于 Wallenius non-central hyper-geometric distribution.相对于普通的超几何分布(Hyper-geometric distribution),此分布的特点是从某个类别中抽取个体的概率与从某个类别之外抽取一个个体的概率是不同的,这种概率的不同是通过对基因长度的偏好性进行估计得到的,从而能更为准确地计算出 GO term 被差异基因富集的概率. 1.GOseq的安