Hi Chanchal,

I have not had a chance yet to look over your code. I will do so in the next
couple of days and get you a definitive answer if I can figure it out. From
the gut though, I suspect that I if remember Seth's implementation
correctly, he looks at the test set and universe set with respect to the
annotation package. If there is anything in the test or universe has yet to
be annotated in the package, these genes are removed before the analysis
takes place. I don't really know if *this* is the issue with which you are
concerned, but like I said, I will have a look at get back to you asap.

Cheers,
--tc


On Sat, Mar 1, 2008 at 2:53 PM, Chanchal Kumar <chanchal@biochem.mpg.de>
wrote:

> Hi James, Tony and All,
>
>   Thanks for your inputs and now I have assembled a small dataset which
> will help you reproduce the inconsistency in GO enrichment output using
> GOstats package.
>
> To summarize, my concern is when I compare any two "test" datasets against
> the common "universe" set then I get different "Category Size" output of
> "universe" for common GO terms across these two "test" datasets. The .txt
> data files (test1, test2, universe) are attached with the mail and the
> script is below. The ids are EntrezID for Human and may be redundant in each
> file. But then I am removing the redundancy before calculating enrichment.
> Also I am using a p-value cutoff of 1 just to get everything.
>
> After analysis I find discrepancies in these GO terms (there are more):
> GO:0016070 ,GO:0051179 ,GO:0006996, GO:0000074 , GO:0031325, GO:0007399,
> GO:0006464, GO:0030097, GO:0007586, GO:0006810, GO:0007243, GO:0048878
>
> I look forward to your inputs. In rest of the mail sessionInfo() is
> followed by script.
>
> >  sessionInfo()
> R version 2.6.2 (2008-02-08)
> i386-pc-mingw32
>
> locale:
> LC_COLLATE=English_United States.1252;LC_CTYPE=English_United States.1252;LC_MONETARY=English_United
> States.1252;LC_NUMERIC=C;LC_TIME=English_United States.1252
>
> attached base packages:
> [1] splines   tools     stats     graphics  grDevices utils     datasets
>  methods   base
>
> other attached packages:
>  [1] hgu95av2.db_2.0.2   GOstats_2.4.0       Category_2.4.0
>  genefilter_1.16.0   survival_2.34       RBGL_1.14.0
>  [7] annotate_1.16.1     xtable_1.5-2        GO.db_2.0.2
> AnnotationDbi_1.0.6 RSQLite_0.6-7       DBI_0.2-4
> [13] Biobase_1.16.3      graph_1.16.1
>
> loaded via a namespace (and not attached):
> [1] cluster_1.11.9
>
>
> ####################SCRIPT STARTS######################
> ######################################################
> # Load packages
> library("GOstats")
> library("hgu95av2.db")
>
>
> # read the universe set
> universe<-read.delim(file="universe.txt",sep="\t",header=T, colClasses =
> "character")
> universe<-as.vector(universe[,1])
> mode(universe)<-"character"
>
> ########## Analysing Test 1 ##################
>
> # read the test set 1
> test<-read.delim(file="test1.txt",sep="\t",header=T, colClasses =
> "character")
> test<-as.vector(test[,1])
> mode(test)<-"character"
>
> hgCutoff <- 1
>
> #testing for GO enrichment using conditional test
>
> params <- new("GOHyperGParams", geneIds = unique(test),
> universeGeneIds = unique(universe), annotation = "hgu95av2",
> ontology = "BP", pvalueCutoff = hgCutoff, conditional = TRUE,
> testDirection = "over")
>
> hgOver <- hyperGTest(params)
>
> # Choose only those categories which have at least 3 members in the
> universe set
> d<-summary(hgOver,categorySize =3)
> write.table(d,"Test1_GOBP.txt",row.names=F,col.names=T,sep="\t",quote=F)
>
>
> ########## Analysing Test 2 ##################
> # read the test set 2
> test<-read.delim(file="test2.txt",sep="\t",header=T, colClasses =
> "character")
> test<-as.vector(test[,1])
> mode(test)<-"character"
>
> hgCutoff <- 1
>
> #testing for GO enrichment using conditional test
>
> params <- new("GOHyperGParams", geneIds = unique(test),
> universeGeneIds = unique(universe), annotation = "hgu95av2",
> ontology = "BP", pvalueCutoff = hgCutoff, conditional = TRUE,
> testDirection = "over")
>
> hgOver <- hyperGTest(params)
>
> # Choose only those categories which have at least 3 members in the
> universe set
> d<-summary(hgOver,categorySize =3)
> write.table(d,"Test2_GOBP.txt",row.names=F,col.names=T,sep="\t",quote=F)
>
> ##########################SCRIPT ENDS##########################
> ###############################################################
>
>
> Best Regards,
> Chanchal
> ===============================
> Chanchal Kumar, Ph.D. Candidate
> Dept. of Proteomics and Signal Transduction
> Max Planck Institute of Biochemistry
> Am Klopferspitz 18
> 82152 D-Martinsried (near Munich)
> Germany
> e-mail: chanchal@biochem.mpg.de
> Phone: (Office) +49 (0) 89 8578 2296
> Fax:(Office) +49 (0) 89 8578 2219
> http://www.biochem.mpg.de/mann/
> ===============================
> ________________________________________
> From: tc.fhcrc@gmail.com [mailto:tc.fhcrc@gmail.com] On Behalf Of Tony
> Chiang
> Sent: Friday, February 29, 2008 3:01 PM
> To: James W. MacDonald
> Cc: Chanchal Kumar; bioconductor@stat.math.ethz.ch
> Subject: Re: [BioC] Inconsistency in GOstats results after GO enrichment
>
> Hi Chanchal,
>
> In addition to the small example, please provide the list with your output
> of sessionInfo(). This way, when we try to reproduce your example, we will
> try to use the same version of software packages that you are using.
>
> Cheers,
> --Tony
>
>
> On 2/29/08, James W. MacDonald <jmacdon@med.umich.edu> wrote:
> Hi Chanchal,
>
> Unless you can give us a small reproducible example that shows what you
> are talking about, there is really nothing anybody can do about this.
>
> Best,
>
> Jim
>
>
>
> Chanchal Kumar wrote:
> > Dear All,
> >
> >
> >
> >    I have been using the "GOstats" package for the GO enrichment on my
> > dataset. I have two cluster of genes(test1, test2) and a reference
> > set(universe set). After doing enrichment on the two sets for GO BP
> > w.r.t to the universe set I get a list of enriched terms. In both cases
> > I get "Cell Cycle" (GO:0007049) as enriched. Till there it's fine. But
> > strangely I see that the Category Size ( "Size" column of the output )
> > gives me different numbers for Cell cycle annotated genes in my universe
> > set( in test1 77, in test2 366). Now my concern is that given I use the
> > same universe set how I can get different numbers of Cell cycle
> > annotated universe genes? This is just one example and I see these
> > inconsistencies for many common terms.
> >
> >
> >
> > Any insights on this problem will be very helpful.
> >
> >
> >
> > Best Regards,
> >
> > Chanchal
> >
> > ===============================
> >
> > Chanchal Kumar, Ph.D. Candidate
> >
> > Dept. of Proteomics and Signal Transduction
> > Max Planck Institute of Biochemistry
> > Am Klopferspitz 18
> > 82152 D-Martinsried (near Munich)
> > Germany
> > e-mail: chanchal@biochem.mpg.de
> > Phone: (Office) +49 (0) 89 8578 2296
> >
> > Fax:(Office) +49 (0) 89 8578 2219
> > http://www.biochem.mpg.de/mann/
> > ===============================
> >
> >
> >
> >
> >       [[alternative HTML version deleted]]
> >
> > _______________________________________________
> > Bioconductor mailing list
> > Bioconductor@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
>
> _______________________________________________
> Bioconductor mailing list
> Bioconductor@stat.math.ethz.ch
> https://stat.ethz.ch/mailman/listinfo/bioconductor
> Search the archives:
> http://news.gmane.org/gmane.science.biology.informatics.conductor
>
>

	[[alternative HTML version deleted]]

