Entering edit mode
Alexander C Cambon
▴
40
@alexander-c-cambon-1251
Last seen 10.2 years ago
I am trying to use the "Custom Processing Methods" available in
package "affy" to create new expression methods. But for some reason I
am getting "NA"'s for the expression values. I start out using 4 cel
files with rae230a data.
> x
AffyBatch object
size of arrays=602x602 features (11330 kb)
cdf=RAE230A (15923 affyids)
number of samples=4
number of genes=15923
annotation=rae230a
## And I proceed by using the example in the "Custom Processing
Methods (How To)" vignette on page 4:
> computeExprVal.huber <- function(probes) {
+ res <- apply(probes, 2, huber)
+ mu <- unlist(lapply(res, function(x) x$mu))
+ s <- unlist(lapply(res, function(x) x$s))
+ return(list(exprs = mu, se.exprs = s))
+ }
>generateExprSet.methods <- c(generateExprSet.methods, "huber")
###Then I use the expresso function... and get an error message (see
below)
>eset.bg.huber<-expresso(x, bgcorrect.method="mas",
pmcorrect.method="pmonly",
+ summary.method="huber")
background correction: mas
PM/MM correction : pmonly
expression values: huber
background correcting...done.
normalizing...done.
Error in match.arg(summary.method, express.summary.stat.methods) :
'arg' should be one of avgdiff, liwong, mas, medianpolish,
playerout
## So, after looking at the error, I decide to add (perhaps wrongly)
the line
>
express.summary.stat.methods<-c(express.summary.stat.methods,"huber")
## and I rerun...
>eset.bg.trm<-expresso(x, bgcorrect.method="mas",
pmcorrect.method="pmonly",
+ summary.method="huber")
## Now I get:
background correction: mas
normalization:
PM/MM correction : pmonly
expression values: huber
background correcting...done.
normalizing...done.
15923 ids to be processed
| |
|####################|
## This looks reassuring, and I start thinking I am over the hump.
However when I ask for the first few lines, I get all ##"NA"'s for
expression ..
> exprs(eset.bg.huber)[1:10,1:4]
AK_0A.CEL AK_0B.CEL AK_4A.CEL AK_4B.CEL
1367452_at NA NA NA NA
1367453_at NA NA NA NA
1367454_at NA NA NA NA
1367455_at NA NA NA NA
1367456_at NA NA NA NA
1367457_at NA NA NA NA
1367458_at NA NA NA NA
1367459_at NA NA NA NA
###But when I do the same thing, using a "built-in" expresso method
(using "x"), I get
> exprs(eset44)[1:10,1:4]
AK_0A.CEL AK_0B.CEL AK_4A.CEL AK_4B.CEL
1367452_at 564.68667 400.93697 476.22935 396.00641
1367453_at 330.41054 281.98687 279.02867 305.24909
1367454_at 158.99028 165.78677 134.08273 151.77046
1367455_at 302.39252 223.56433 192.13621 204.74223
1367456_at 881.07080 931.81819 792.18972 583.09672
1367457_at 139.75233 142.67949 140.01944 153.02465
1367458_at 144.34063 154.16194 176.33189 190.89218
1367459_at 1010.72613 709.47988 886.92512 957.64834
1367460_at 48.19831 51.70533 59.43114 54.13907
1367461_at 213.03203 194.16510 210.99492 182.22481
Any ideas what I am doing wrong?
Alex
Alexander Cambon
Biostatistician
School of Public Health
Dept of Biostatistics and Bioinformatics
University of Louisville
R 2.1.0
Bioconductor 1.6
Windows XP
Dell Optiplex
1 GB RAM