Entering edit mode
Tobias Straub
▴
430
@tobias-straub-2182
Last seen 9.6 years ago
hi
i was playing around with various background correction methods using
limma and discovered that whatever normexp.method i call for the
"normexp" method, the result is exactly the same. in fact the more i
digged into it the more i believe that normexp is not executed at all.
looking at the code for backgroundCorrect it more looks like some junk
of code has been lost. i also checked the latest dev version of limma
(3.0.3), looks the same. instead of doing normexp the execution jumps
into some "rma" kind of thing that just leaves me confused.
any clues someone? any correct code around somewhere?
best
T.
> backgroundCorrect
function (RG, method = "subtract", offset = 0, printer = RG$printer,
normexp.method = "saddle", verbose = TRUE)
{
if (!is.list(RG) && is.vector(RG))
RG <- as.matrix(RG)
if (is.matrix(RG)) {
method <- match.arg(method, c("none", "normexp", "saddle",
"neldermean", "bfgs", "rma", "mcgee"))
if (method == "normexp")
method = "saddle"
if (method != "none") {
for (j in 1:ncol(RG)) {
x <- RG[, j]
out <- normexp.fit(x, method = normexp.method)
RG[, j] <- normexp.signal(out$par, x)
if (verbose)
cat("Corrected array", j, "\n")
}
}
if (offset)
RG <- RG + offset
return(RG)
}
if (is.null(RG$Rb) != is.null(RG$Gb))
stop("Background values exist for one channel but not the
other")
method <- match.arg(method, c("none", "subtract", "half",
"minimum", "movingmin", "edwards", "normexp", "rma"))
if (method == "rma") {
method <- "normexp"
normexp.method <- "rma"
}
normexp.method <- match.arg(normexp.method, c("mle", "saddle",
"rma", "rma75", "mcgee"))
if (is.null(RG$Rb) && is.null(RG$Gb))
method <- "none"
switch(method, subtract = {
RG$R <- RG$R - RG$Rb
RG$G <- RG$G - RG$Gb
}, half = {
RG$R <- pmax(RG$R - RG$Rb, 0.5)
RG$G <- pmax(RG$G - RG$Gb, 0.5)
}, minimum = {
RG$R <- as.matrix(RG$R - RG$Rb)
RG$G <- as.matrix(RG$G - RG$Gb)
for (slide in 1:ncol(RG$R)) {
i <- RG$R[, slide] < 1e-18
if (any(i, na.rm = TRUE)) {
m <- min(RG$R[!i, slide], na.rm = TRUE)
RG$R[i, slide] <- m/2
}
i <- RG$G[, slide] < 1e-18
if (any(i, na.rm = TRUE)) {
m <- min(RG$G[!i, slide], na.rm = TRUE)
RG$G[i, slide] <- m/2
}
}
}, movingmin = {
RG$R <- RG$R - ma3x3.spottedarray(RG$Rb, printer = printer,
FUN = min, na.rm = TRUE)
RG$G <- RG$G - ma3x3.spottedarray(RG$Gb, printer = printer,
FUN = min, na.rm = TRUE)
}, edwards = {
one <- matrix(1, NROW(RG$R), 1)
delta.vec <- function(d, f = 0.1) {
quantile(d, mean(d < 1e-16, na.rm = TRUE) * (1 +
f), na.rm = TRUE)
}
sub <- as.matrix(RG$R - RG$Rb)
delta <- one %*% apply(sub, 2, delta.vec)
RG$R <- ifelse(sub < delta, delta * exp(1 - (RG$Rb +
delta)/RG$R), sub)
sub <- as.matrix(RG$G - RG$Gb)
delta <- one %*% apply(sub, 2, delta.vec)
RG$G <- ifelse(sub < delta, delta * exp(1 - (RG$Gb +
delta)/RG$G), sub)
}, normexp = , rma = {
if (verbose)
cat("Green channel\n")
RG$G <- backgroundCorrect(RG$G - RG$Gb, method = method,
verbose = verbose)
if (verbose)
cat("Red channel\n")
RG$R <- backgroundCorrect(RG$R - RG$Rb, method = method,
verbose = verbose)
})
RG$Rb <- NULL
RG$Gb <- NULL
if (offset) {
RG$R <- RG$R + offset
RG$G <- RG$G + offset
}
new("RGList", unclass(RG))
}
<environment: namespace:limma="">
> sessionInfo()
R version 2.9.2 RC (2009-08-17 r49303)
x86_64-apple-darwin9.7.0
locale:
en_US.UTF-8/en_US.UTF-8/C/C/en_US.UTF-8/en_US.UTF-8
attached base packages:
[1] stats graphics grDevices utils datasets methods base
other attached packages:
[1] vsn_3.12.0 Biobase_2.4.1 limma_2.18.3
loaded via a namespace (and not attached):
[1] affy_1.22.1 affyio_1.12.0 grid_2.9.2
lattice_0.17-25 preprocessCore_1.6.0
>
----------------------------------------------------------------------
Dr. Tobias Straub ++4989218075439 Adolf-Butenandt-Institute, M?nchen D
