[BioC] limma normexp background correction bug?

[Tobias Straub]

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