[BioC] Strange FDR values using GSA package with large numbers of permutations
Wolfgang Huber
whuber at embl.de
Fri Sep 18 11:53:29 CEST 2009
Hi Alexander
it would be great if someone subscribing to this list has an insight on
this question. However, note that GSA is not a Bioconductor package, so
you may need to explore further avenue such as contacting its authors or
examining the source code of the function.
Also, please note that the probability of getting a useful answer is
increased if you provide a readily reproducible example.
Best wishes
Wolfgang
PS It is curious that the numeric values in your example are reported as
character strings.
Alexander C Cambon wrote:
> I am using the GSA package for gene set analysis. I have normalize data from a microarray experiment with four groups and four replicates (16 total arrays).
>
> I downloaded the set of curated gene sets ("c2.all.v2.5.symbols.gmt") from the Broad Institute web page MSigDB (I did register) and read the data into R using the GSA package as follows:
>
> geneset.obj<- GSA.read.gmt("c2.all.v2.5.symbols.gmt")
>
> GSA.obj3<-GSA(x,y, genenames=gn, genesets=geneset.obj$genesets, resp.type="Multiclass", nperms=10000)
>
> (I also tried a smaller number of permutations).
>
> Then, I got the gene set list using
>
> GSA.3<-GSA.listsets(GSA.obj3, geneset.names=geneset.obj$geneset.names,FDRcut=.5)
>
>> dim(GSA.3$positive)
> [1] 721 5
>
>
> I noticed that when I used 10000 permutations, the p-values for the gene sets varied all the way from 0.0013 for the top gene set to 0.9434 for the gene sets at the bottom, but the FDR values for all but the last of the 721 gene sets stayed at 0.0749. This did not happen, at least not this severely, with a smaller number of permutations.
>
> Does anyone have an explanation?
>
> Alexander Cambon
> Biostatistician
> Dept of Bioinformatics and Biostatistics
> School of Public Health and Information Sciences
> University of Louisville
> Louisville, KY
>
>
>
> Here is a sample of some of the first ones (I x'ed out the gene set numbers and names)
>
> Gene_set Gene_set_name Score p-value FDR
>
> "xxx" "xxx" "0.0994" "0.0013" "0.0749"
> "xxx" "xxx" "0.1622" "0.0013" "0.0749"
> "xx" "xxxx" "0.3668" "0.0016" "0.0749"
> "xxx" "xxxx" "0.2887" "0.0016" "0.0749"
>
>
> Here are the last ones
> Gene_set Gene_set_name Score p-value FDR
>
>
> "xxx" "xxx" "0.0285" "0.8882" "0.0749"
> "xxx" " xxx" "0.024" "0.901" "0.0749"
> "xxx" "xxx" "0.0035" "0.9434" "0.0749"
>
>
> I am using Widows XP
>> sessionInfo()
> R version 2.9.1 (2009-06-26)
> 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] stats graphics grDevices utils datasets methods base
>
> other attached packages:
> [1] GSA_1.0
>
> _______________________________________________
> 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
More information about the Bioconductor
mailing list