[BioC] Gostats question

[James W. MacDonald]

Hi Dongxiao,

Zhu, Dongxiao wrote:
> Hi, I am trying to use Gostats to calculate enrichment or depletion of
> GO terms. The below is the code adapted from vignette. Assume we start
> from a mouse 4302 expressionset object "exprSet1000". It always returns
> empty table no matter what data I try.
> 
> Can somebody help me out?

Yes. If you look at the code below, you will notice that your geneIds 
are identical to your universeGeneIds. I don't think this makes any 
sense - you want the geneIds to be a subset of your universeGeneIds 
(which technically I guess they are, so you want a subset where geneIds 
< universeGeneIds).

Best,

Jim


> 
> Thanks.
> 
> ## General Filtering, remove probe sets without Entrez ID
> entrezIds <- mget(featureNames(exprSet1000), envir = mouse4302ENTREZID)
> haveEntrezId <- names(entrezIds)[sapply(entrezIds, function(x)
> !is.na(x))]
> numNoEntrezId <- length(featureNames(exprSet1000)) -
> length(haveEntrezId)
> bcrAblOrNeg <- exprSet1000[haveEntrezId, ]
> 
> ## Remove probe sets for which we have no GO annotation.
> haveGo <- sapply(mget(featureNames(bcrAblOrNeg), mouse4302GO),
> function(x) {
> if (length(x) == 1 && is.na(x))
> FALSE
> else TRUE
> })
> numNoGO <- sum(!haveGo)
> bcrAblOrNeg <- bcrAblOrNeg[haveGo, ]
> 
> ## remove redundant probe set ids
> bcrAblOrNegIqr <- apply(exprs(bcrAblOrNeg), 1, IQR)
> numNsWithDups <- length(featureNames(bcrAblOrNeg))
> uniqGenes <- findLargest(featureNames(bcrAblOrNeg), bcrAblOrNegIqr,
> "mouse4302")
> bcrAblOrNeg <- bcrAblOrNeg[uniqGenes, ]
> 
> entrezUniverse <- unlist(mget(featureNames(bcrAblOrNeg),
> mouse4302ENTREZID))
> if (any(duplicated(entrezUniverse))) stop("error in gene universe: can't
> have duplicate Entrez")  
> selectedEntrezIds <- unlist(mget(featureNames(bcrAblOrNeg),
> mouse4302ENTREZID))
> 
> hgCutoff <- 1
> params <- new("GOHyperGParams", geneIds = selectedEntrezIds, 
> universeGeneIds = entrezUniverse, annotation = "mouse4302",
> ontology = "BP", pvalueCutoff = hgCutoff, conditional = FALSE,
> testDirection = "over")
> 
> paramsCond <- params
> conditional(paramsCond) <- TRUE
> hgOver <- hyperGTest(params)
> hgCondOver <- hyperGTest(paramsCond)
> 
> 
> 
> -------------------------------------------------
> Dongxiao Zhu
> Stowers Institute for Medical Research
> Phone: 816-926-4456
> FAX: 816-926-4674
> 
> 
> 	[[alternative HTML version deleted]]
> 
> _______________________________________________
> Bioconductor mailing list
> Bioconductor at stat.math.ethz.ch
> https://stat.ethz.ch/mailman/listinfo/bioconductor
> Search the archives: http://news.gmane.org/gmane.science.biology.informatics.conductor


-- 
James W. MacDonald, M.S.
Biostatistician
Affymetrix and cDNA Microarray Core
University of Michigan Cancer Center
1500 E. Medical Center Drive
7410 CCGC
Ann Arbor MI 48109
734-647-5623


**********************************************************
Electronic Mail is not secure, may not be read every day, and should not be used for urgent or sensitive issues.