[BioC] help with multtest and annaffy

James W. MacDonald jmacdon at med.umich.edu
Tue May 18 21:45:27 CEST 2010


Hi Katie,

Taylor, Katie wrote:
> Hi,
> 
> I have followed the annaffy vignette online and successfully created the HTML documents for annotating the aafExpr dataset. When I try to do the same with my data (GDS2609) which I downloaded from GEO I can only get so far and then R becomes unresponsive and is forced to close. I have copied the session info below:
> 
> 'help.start()' for an HTML browser interface to help.
> Type 'q()' to quit R.
> 
>> library(affy); library(limma)
> Loading required package: Biobase
> 
> Welcome to Bioconductor
> 
>   Vignettes contain introductory material. To view, type
>   'openVignette()'. To cite Bioconductor, see
>   'citation("Biobase")' and for packages 'citation(pkgname)'.
> 
>> Data <- ReadAffy()
>> Gds2609 <- Data
>> eset <- rma(Gds2609)
> Background correcting
> Normalizing
> Calculating Expression
>> pData(eset)
>              sample
> GSM93789.CEL      1
> GSM93920.CEL      2
> GSM93952.CEL      3
> GSM93954.CEL      4
>> strain <- c("n", "n", "c", "c")
>> design <- model.matrix(~factor(strain))
>> colnames(design) <- c("n", "c")
>> design
>   n c
> 1 1 1
> 2 1 1
> 3 1 0
> 4 1 0
> attr(,"assign")
> [1] 0 1
> attr(,"contrasts")
> attr(,"contrasts")$`factor(strain)`
> [1] "contr.treatment"

Unless you are using limma, the above (after pData(eset)) is superfluous.

>  library(multtest)
>> class <- as.integer(pData(Gds2609)$covar1) - 1

As Sean Davis mentioned the last time you asked, class isn't a very good 
name, as it is also a function name. Using something like classlabel 
instead is a more reasonable approach. In addition, you might show what 
is in this vector.

In addition, it looks like you are getting mixed up in all your 
gyrations above (why call the AffyBatch 'Data' and then immediately copy 
it into another AffyBatch called Gds2609? You are just chewing up RAM 
unnecessarily.). The 'Gds2609' object is an AffyBatch, not an 
ExpressionSet. Although exprs() will work on an AffyBatch, the matrix 
you are outputting is over ten times the size of the matrix that I think 
you want (exprs(eset)).


Best,

Jim




>>  teststat <- mt.teststat(exprs(Gds2609), class) ##at this point R is forced to close each time I try it. 
>> sessionInfo() ## had to do this command before teststat otherwise R closes before I can do it
> R version 2.11.0 (2010-04-22) 
> i386-pc-mingw32 
> 
> locale:
> [1] LC_COLLATE=English_United Kingdom.1252 
> [2] LC_CTYPE=English_United Kingdom.1252   
> [3] LC_MONETARY=English_United Kingdom.1252
> [4] LC_NUMERIC=C                           
> [5] LC_TIME=English_United Kingdom.1252    
> 
> attached base packages:
> [1] stats     graphics  grDevices utils     datasets  methods   base     
> 
> other attached packages:
> [1] multtest_2.4.0       hgu133plus2cdf_2.6.0 limma_3.4.0         
> [4] affy_1.26.0          Biobase_2.8.0       
> 
> loaded via a namespace (and not attached):
> [1] affyio_1.16.0         MASS_7.3-5            preprocessCore_1.10.0
> [4] splines_2.11.0        survival_2.35-8       tools_2.11.0         
> 
> 
> The commands that I want to use after this are as follows:
> index <- order(abs(teststat), decreasing = TRUE)
> probeids <- featureNames(Gds2609)[index]
> library(annaffy)
> aaf.handler()
> anncols <- aaf.handler()[c(1:3, 8:9, 11:13)]
> anntable <- aafTableAnn(probeids[1:250], "hgu133plus2.db", anncols)
> saveHTML(anntable, "example1.html", title = "Example Table without Data")
> testtable <- aafTable(`t-statistic` = teststat[index[1:250]],
> signed = TRUE)
> table <- merge(anntable, testtable)
> exprtable <- aafTableInt(Gds2609, probeids = probeids[1:250])
> table <- merge(table, exprtable)
> saveHTML(table, "example2.html", title = "Example Table with Data")
> saveText(table, "example2.txt")
> 
> 
> If anyone has any ideas why R keeps closing on me like this I would really appreciate your help. I've tried running R on different computers as I didn't know if it was a memory problem but the example in the annaffy pdf works fine. I have no idea what is happening and could really do with some help.
> 
> Best Wishes,
> 
> Katie
> 
> _______________________________________________
> 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
Douglas Lab
University of Michigan
Department of Human Genetics
5912 Buhl
1241 E. Catherine St.
Ann Arbor MI 48109-5618
734-615-7826
**********************************************************
Electronic Mail is not secure, may not be read every day, and should not be used for urgent or sensitive issues 



More information about the Bioconductor mailing list