我想做路径浓缩分析。我有21个显着基因的列表,以及我想要检查的多种路径类型。检查KEGG途径、GOterms、配合物等的富集情况。
我在一个旧的BioC帖子上找到了这个代码示例。然而,我很难自己适应它。
首先,1-这是什么意思?我不知道这个多重冒号语法。
hyperg <- Category:::.doHyperGInternal我不明白这句话是怎么回事。hyperg.test是一个需要传递给它的3个变量的函数,对吗?这一行是否以某种方式将"genes.by.pathways、significant.genes和all.geneIDs“传递给thr hyperg.test?
pVals.by.pathway<-t(sapply(genes.by.pathway, hyperg.test, significant.genes, all.geneIDs))我想要修改的代码
library(KEGGREST)
library(org.Hs.eg.db)
# created named list, length 449, eg:
# path:hsa00010: "Glycolysis / Gluconeogenesis"
pathways <- keggList("pathway", "hsa")
# make them into KEGG-style human pathway identifiers
human.pathways <- sub("path:", "", names(pathways))
# for demonstration, just use the first ten pathways
demo.pathway.ids <- head(human.pathways, 10)
demo.pathways <- setNames(keggGet(demo.pathway.ids), demo.pathway.ids)
genes.by.pathway <- lapply(demo.pathways, function(demo.pathway) {
demo.pathway$GENE[c(TRUE, FALSE)]
})
all.geneIDs <- keys(org.Hs.eg.db)
# chose one of these for demonstration. the first (a whole genome random
# set of 100 genes) has very little enrichment, the second, a random set
# from the pathways themselves, has very good enrichment in some pathways
set.seed(123)
significant.genes <- sample(all.geneIDs, size=100)
#significant.genes <- sample(unique(unlist(genes.by.pathway)), size=10)
# the hypergeometric distribution is traditionally explained in terms of
# drawing a sample of balls from an urn containing black and white balls.
# to keep the arguments straight (in my mind at least), I use these terms
# here also
hyperg <- Category:::.doHyperGInternal
hyperg.test <-
function(pathway.genes, significant.genes, all.genes, over=TRUE)
{
white.balls.drawn <- length(intersect(significant.genes, pathway.genes))
white.balls.in.urn <- length(pathway.genes)
total.balls.in.urn <- length(all.genes)
black.balls.in.urn <- total.balls.in.urn - white.balls.in.urn
balls.pulled.from.urn <- length(significant.genes)
hyperg(white.balls.in.urn, black.balls.in.urn,
balls.pulled.from.urn, white.balls.drawn, over)
}
pVals.by.pathway <-
t(sapply(genes.by.pathway, hyperg.test, significant.genes, all.geneIDs))
print(pVals.by.pathway)发布于 2014-11-24 14:16:40
您得到错误的原因是,您似乎没有从生物导体安装Category包。我怀疑这是因为三冒号运算符:::。这个运算符非常类似于双冒号操作符::。使用::,您可以在不加载包的情况下访问导出的对象,而:::则允许访问非导出的对象(在本例中,hyperg函数来自Category)。如果安装Category包,则代码运行时不会出错。
关于sapply语句:
pVals.by.pathway<-t(sapply(genes.by.pathway, hyperg.test, significant.genes, all.geneIDs))你可以把它分解成不同的部分来理解它。首先,sapply迭代gene.by.pathway的元素,并将它们传递给hyperg.test的第一个参数。以下参数是两个加法参数。这是有点不清楚,我个人建议人们明确识别参数,以避免意外的惊喜,并避免需要完全相同的顺序。在这种情况下,这有点重复,但这是避免愚蠢的bug的好方法(例如,将significant.genes放在all.geneIds之后)
重写:
pVals.by.pathway <-
t(sapply(genes.by.pathway, hyperg.test, significant.genes=significant.genes, all.genes=all.geneIDs))一旦这个循环完成,sapply函数将把输出简化为一个矩阵。但是,通过采用转置式t,输出更加方便用户。
一般来说,当我试图理解复杂的apply语句时,我发现最好将它们分解成更小的部分,看看对象本身的样子。
https://stackoverflow.com/questions/27106208
复制相似问题