New to R/Bioconductor and all that so forgive my lack of R chops or understanding of objects like ExpressionSet. Just finished reading the manuals...
I am looking to create a function that takes a GDS id and various control GSM ids and find the most up-regulated gene and most down-regulated gene. Looking at GDS1731 with controls of GSM84549, GSM84550, GSM84551, I can setup the following:
library('GEOquery') library('samr') gds <- getGEO('GDS1731') eset <- GDS2eSet(gds) x = exprs(eset) y = c(rep(1,3),rep(2,9)) data <- list(x = x, y = y, logged2 = FALSE) samr.out <- samr(data = data, resp.type = "Two class unpaired", assay.type = "array") delta.table <- samr.compute.delta.table(samr.out) delta=.5 siggenes.table <- samr.compute.siggenes.table(samr.out, delta, data, delta.table) up_genes <- siggenes.table$genes.up lo_genes <- siggenes.table$genes.lo
And then take the first row of upgenes or logenes.
What I can't figure out is how to get the geneid and genenames into the mix. They look like parameters to the SAM method but not samr. Nor are they parameters to the samr.compute.siggenes.table - that I can see. Also its not clear to me how to get geneid and genenames out of the ExpressionSet - all depending on what gene identifier I want to use and what is in the GDS/GPL.
Any clues thoughts?