[Robast-commits] r1333 - in pkg/RobLoxBioC: . R man

noreply at r-forge.r-project.org noreply at r-forge.r-project.org
Sun Jan 19 17:58:15 CET 2025


Author: stamats
Date: 2025-01-19 17:58:14 +0100 (Sun, 19 Jan 2025)
New Revision: 1333

Modified:
   pkg/RobLoxBioC/DESCRIPTION
   pkg/RobLoxBioC/R/KolmogorovMinDist.R
   pkg/RobLoxBioC/man/robloxbioc.Rd
Log:
removed .Internal

Modified: pkg/RobLoxBioC/DESCRIPTION
===================================================================
--- pkg/RobLoxBioC/DESCRIPTION	2025-01-12 14:32:59 UTC (rev 1332)
+++ pkg/RobLoxBioC/DESCRIPTION	2025-01-19 16:58:14 UTC (rev 1333)
@@ -1,6 +1,6 @@
 Package: RobLoxBioC
 Version: 1.2.3
-Date: 2025-01-12
+Date: 2025-01-19
 Title: Infinitesimally Robust Estimators for Preprocessing -Omics Data
 Description: Functions for the determination of optimally robust influence curves and
             estimators for preprocessing omics data, in particular gene expression data (Kohl

Modified: pkg/RobLoxBioC/R/KolmogorovMinDist.R
===================================================================
--- pkg/RobLoxBioC/R/KolmogorovMinDist.R	2025-01-12 14:32:59 UTC (rev 1332)
+++ pkg/RobLoxBioC/R/KolmogorovMinDist.R	2025-01-19 16:58:14 UTC (rev 1333)
@@ -1,343 +1,343 @@
-###############################################################################
-## Methods for KolmogorovMinDist: computation of minimum Kolmogorov distance
-###############################################################################
-setMethod("KolmogorovMinDist", signature(x = "matrix",
-                                         D = "Norm"),
-    function(x, D, mad0 = 1e-4){
-        stopifnot(is.numeric(x))
-        ksdist <- function(pars, x){
-            if(any(ina <- is.na(x))) x <- x[!ina]
-            y <- .Internal(sort(x, FALSE))
-            p.y <- pnorm(y, mean = pars[1], sd = pars[2])
-            n <- length(y)
-            max(p.y - (seq_len(n)-1)/n, seq_len(n)/n - p.y)
-        }
-        Med <- rowMedians(x, na.rm = TRUE)
-        Mad <- pmax(rowMedians(abs(x-Med), na.rm = TRUE)/qnorm(0.75), mad0)
-        startPars <- cbind(Med, Mad)
-        M <- nrow(x)
-        n <- ncol(x)
-        res <- matrix(NA, nrow = M, ncol = 2)
-        for(i in seq_len(M)){
-            temp <- try(optim(par = startPars[i,], fn = ksdist, x = x[i,], 
-                              method = "L-BFGS-B", lower = c(-Inf, 1e-15), 
-                              upper = c(Inf, Inf))$value, silent = TRUE)
-            if(!inherits(temp, "try-error")){
-                res[i,1] <- temp
-                res[i,2] <- n
-            }
-        }
-        list("dist" = res[,1], "n" = res[,2])
-    })
-
-setMethod("KolmogorovMinDist", signature(x = "AffyBatch",
-                                         D = "AbscontDistribution"),
-    function(x, D, bg.correct = TRUE, pmcorrect = TRUE, verbose = TRUE){
-        if(bg.correct){
-            if(verbose) cat("Background correcting ...")
-            x <- affy::bg.correct.mas(x, griddim = 16)
-            if(verbose) cat(" done.\n")
-        }
-        n <- length(x)
-        ids <- featureNames(x)
-        m <- length(ids)
-
-        CDFINFO <- getCdfInfo(x)
-        INDEX <- sapply(ids, get, envir = CDFINFO)
-        NROW <- unlist(lapply(INDEX, nrow))
-        nr <- as.integer(names(table(NROW)))
-
-        intensData <- intensity(x)
-        kmd <- numeric(m*n)
-        ns <- numeric(m*n)
-        if(pmcorrect){
-            if(verbose) cat("Calculating PM/MM ...")
-            div <- function(INDEX, x){
-                l.pm <- INDEX[,1]
-                if(ncol(INDEX) == 2)
-                    l.mm <- INDEX[,2]
-                else
-                    l.mm <- integer()
-                x[l.pm, , drop = FALSE]/x[l.mm, , drop = FALSE]
-            }
-            res <- lapply(INDEX, div, x = intensData)
-            if(verbose) cat(" done.\n")
-        }else{
-            if(verbose) cat("Extract PM data ...")
-            pm.only <- function(INDEX, x){
-                l.pm <- INDEX[,1]
-                x[l.pm, , drop = FALSE]
-            }
-            res <- lapply(INDEX, pm.only, x = intensData)
-            if(verbose) cat(" done.\n")
-        }
-        if(verbose) cat("Computing minimum Kolmogorov distance ...")
-        for(k in nr){
-            ind <- which(NROW == k)
-            temp <- matrix(do.call(rbind, res[ind]), nrow = k)
-            ind1 <-  as.vector(sapply(seq_len(n)-1, function(x, ind, m){ ind + x*m }, ind = ind, m = m))
-            temp.res <- KolmogorovMinDist(log2(t(temp)), D = D)
-            kmd[ind1] <- temp.res[["dist"]]
-            ns[ind1] <- temp.res[["n"]]
-        }
-        if(verbose) cat(" done.\n")
-        kmd.mat <- matrix(kmd, nrow = m)
-        ns.mat <- matrix(ns, nrow = m)
-        dimnames(kmd.mat) <- list(ids, sampleNames(x))
-        dimnames(ns.mat) <- list(ids, sampleNames(x))
-        list("dist" = kmd.mat, "n" = ns.mat)
-    })
-
-setMethod("KolmogorovMinDist", signature(x = "beadLevelData",
-                                         D = "AbscontDistribution"),
-    function(x, D, log = FALSE, what = "Grn", probes = NULL, arrays = NULL){
-        BLData <- x
-        arraynms <- sectionNames(BLData)
-        if(!is.null(arrays) && !is.character(arrays)) arraynms <- arraynms[arrays]
-        if(is.character(arrays)) arraynms <- which(arraynms %in% arrays)
-        len <- length(arraynms)
-
-        tmp <- BLData[[1]]
-        what <- match.arg(what, colnames(tmp))
-        if(is.na(what))
-          stop("Could not find bead data of type ", what, "\n The available choices are: ",
-               paste(colnames(tmp)), "\n")
-
-        pr <- getBeadData(BLData, what = "ProbeID", array = 1)
-        finten <- getBeadData(BLData, what = what, array = 1)
-        if(log) finten <- log2(finten)
-        nasinf <- !is.finite(finten) | is.na(finten)
-        finten <- finten[!nasinf]
-
-        if(is.null(probes)) probes <- sort(unique(pr))
-        probes <- probes[probes > 0 & !is.na(probes)]
-        noprobes <- length(probes)
-        pr <- pr[!nasinf]
-        G.kmd <- matrix(0, nrow = noprobes, ncol = len)
-        G.ns <- matrix(0, nrow = noprobes, ncol = len)
-        colnames(G.kmd) <- colnames(G.ns) <- arraynms
-
-        i <- 1
-        while (i <= len) {
-            probeIDs <- as.integer(pr)
-            start <- 0
-            G.res <- kmdBeadLevel(x = finten, D = D, probeIDs = probeIDs, probes = probes)
-            G.kmd[,i] <- G.res[["dist"]]
-            G.ns[,i] <- G.res[["n"]]
-
-            i <- i + 1
-            rm(probeIDs)
-            gc()
-
-            pr <- getBeadData(BLData, what = "ProbeID", array = i)
-            finten <- getBeadData(BLData, what = what, array = i)
-            if(log) finten <- log2(finten)
-            nasinf <- !is.finite(finten) | is.na(finten)
-            pr <- pr[!nasinf]
-            finten <- finten[!nasinf]
-        }
-        rownames(G.kmd) <- rownames(G.ns) <- probes
-        res <- list("dist" = G.kmd, "n" = G.ns)
-        return(res)
-    })
-
-kmdBeadLevel <- function(x, D, probeIDs, probes){
-    comIDs <- intersect(probeIDs, probes)
-    x <- x[probeIDs %in% comIDs]
-    probeIDs <- probeIDs[probeIDs %in% comIDs]
-    noBeads <- as.vector(table(probeIDs))
-    noBeads.uni <- as.integer(names(table(noBeads)))
-    probes1 <- comIDs
-    len1 <- length(probes1)
-    kmd1 <- numeric(len1)
-    ns1 <- numeric(len1)
-    for(i in seq(along = noBeads.uni)){
-        index <- noBeads == noBeads.uni[i]
-        IDs <- probes1[index]
-        if(noBeads.uni[i] == 1){
-            kmd1[index] <- 0.5
-        }else{
-            temp <- matrix(x[probeIDs %in% IDs], ncol = noBeads.uni[i], byrow = TRUE)
-            res <- KolmogorovMinDist(temp, D = D)
-            kmd1[index] <- res[["dist"]]
-            ns1[index] <- res[["n"]]
-        }
-    }
-    len <- length(probes)
-    kmd <- numeric(len)
-    ns <- numeric(len)
-    nas <- !(probes %in% comIDs)
-    kmd[nas] <- NA
-    kmd[!nas] <- kmd1
-    ns[nas] <- NA
-    ns[!nas] <- ns1
-
-    list("dist" = kmd, "n" = ns)
-}
-
-## leftover in trunk from revision 824 ....
-
-# setMethod("KolmogorovMinDist", signature(x = "BeadLevelList",
-#                                          D = "AbscontDistribution"),
-#     function(x, D, log = FALSE, imagesPerArray = 1, what = "G", probes = NULL, arrays = NULL){
-#         BLData <- x
-#         arraynms <- arrayNames(BLData)
-#         if(!is.null(arrays) && !is.character(arrays)) arraynms <- arraynms[arrays]
-#         if(is.character(arrays)) arraynms <- which(arraynms %in% arrays)
-#         len <- length(arraynms)
-#         what <- match.arg(what, c("G", "R", "RG", "M", "A", "beta"))
-#         whatelse <- ""
-#         if(what == "RG"){
-#             if(BLData at arrayInfo$channels == "two"){
-#                 what <- "G"
-#                 whatelse <- "R"
-#             }else{
-#                 stop("Need two-channel data to calculate summary R and G values")
-#             }
-#         }
-#         if(imagesPerArray == 1){
-#             sel <- getArrayData(BLData, what = "ProbeID", array = arraynms[1]) != 0
-#             pr <- getArrayData(BLData, what = "ProbeID", array = arraynms[1])[sel]
-#             finten <- getArrayData(BLData, what = what, log = log, array = arraynms[1])[sel]
-#             nasinf <- !is.finite(finten) | is.na(finten)
-#             finten <- finten[!nasinf]
-#        }
-#        else if(imagesPerArray == 2){
-#            if(length(arraynms)%%2 != 0)
-#                stop("Need an even number of arrays when 'imagesPerArray=2'")
-#            arrayord <- order(arraynms)
-#            arraynms <- arraynms[arrayord]
-#            tmp <- unlist(strsplit(arraynms, "_"))
-#            chipnums <- tmp[seq(1, length(tmp), by = 3)]
-#            pos <- tmp[seq(2, length(tmp), by = 3)]
-#            stripnum <- as.numeric(tmp[seq(3, length(tmp), by = 3)])
-#            check <- ((chipnums[seq(1, len, by = 2)] == chipnums[seq(2, len, by = 2)])
-#                      & (pos[seq(1, len, by = 2)] == pos[seq(2, len, by = 2)])
-#                      & (stripnum[seq(1, len, by = 2)] == stripnum[seq(2, len, by = 2)] - 1))
-#            if(sum(check) != length(check)) stop("Missing arrays")
-#            sel1 <- getArrayData(BLData, what = "ProbeID", array = arraynms[1]) != 0
-#            sel2 <- getArrayData(BLData, what = "ProbeID", array = arraynms[2]) != 0
-#            pr <- append(getArrayData(BLData, what = "ProbeID", array = arraynms[1])[sel1],
-#                        getArrayData(BLData, what = "ProbeID", array = arraynms[2])[sel2])
-#            finten <- append(getArrayData(BLData, what = what, log = log, array = arraynms[1])[sel1],
-#                            getArrayData(BLData, what = what, log = log, array = arraynms[2])[sel2])
-#            nasinf <- !is.finite(finten) | is.na(finten)
-#            finten <- finten[!nasinf]
-#        }else{
-#            stop("You can only specify 1 or 2 images per array")
-#        }
-#        if(is.null(probes)) probes <- sort(unique(pr))
-#        probes <- probes[probes > 0 & !is.na(probes)]
-#        noprobes <- length(probes)
-#        pr <- pr[!nasinf]
-#        if (imagesPerArray == 1) {
-#            G.kmd <- matrix(0, nrow = noprobes, ncol = len)
-#            G.ns <- matrix(0, nrow = noprobes, ncol = len)
-#            colnames(G.kmd) <- colnames(G.ns) <- arraynms
-#            if (BLData at arrayInfo$channels == "two" && !is.null(BLData[[arraynms[1]]]$R) && whatelse == "R")
-#                R.kmd <- R.ns <- G.kmd
-#            else R.kmd <- R.ns <- NULL
-#        }
-#        else if (imagesPerArray == 2) {
-#            G.kmd <- matrix(0, nrow = noprobes, ncol = (len/2))
-#            G.ns <- matrix(0, nrow = noprobes, ncol = (len/2))
-#            colnames(G.kmd) <- colnames(G.ns) <- arraynms[seq(1, len, by = 2)]
-#            if (BLData at arrayInfo$channels == "two" && !is.null(BLData[[arraynms[1]]]$R) && whatelse == "R")
-#                R.kmd <- R.ns <- G.kmd
-#            else R.kmd <- R.ns <- NULL
-#        }
-#        i <- j <- 1
-#        while (j <= len) {
-#            probeIDs <- as.integer(pr)
-#            start <- 0
-#            G.res <- kmdBeadLevel(x = finten, D = D, probeIDs = probeIDs, probes = probes)
-#            G.kmd[,i] <- G.res[["dist"]]
-#            G.ns[,i] <- G.res[["n"]]
-#            if (BLData at arrayInfo$channels == "two" && !is.null(BLData[[arraynms[i]]]$R) && whatelse == "R") {
-#                if (imagesPerArray == 1) {
-#                    finten <- getArrayData(BLData, what = whatelse, log = log, array = arraynms[i])[sel]
-#                    nasinf <- !is.finite(finten) | is.na(finten)
-#                    finten <- finten[!nasinf]
-#                    binten <- rep(0, length(finten))
-#                }
-#                else if (imagesPerArray == 2) {
-#                    finten <- append(getArrayData(BLData, what = whatelse, log = log, array = arraynms[j])[sel1],
-#                                    getArrayData(BLData, what = whatelse, log = log, array = arraynms[j + 1])[sel2])
-#                    nasinf <- !is.finite(finten) | is.na(finten)
-#                    finten <- finten[!nasinf]
-#                    binten <- rep(0, length(finten))
-#                }
-#                start <- 0
-#                R.res[,i] <- kmdBeadLevel(x = finten, D = D, probeIDs = probeIDs, probes = probes)
-#                R.kmd[,i] <- R.res[["dist"]]
-#                R.ns[,i] <- R.res[["n"]]
-#            }
-#            j <- j + imagesPerArray
-#            i <- i + 1
-#            rm(probeIDs)
-#            gc()
-#            if ((imagesPerArray == 1) && (i <= len)) {
-#                sel <- getArrayData(BLData, what = "ProbeID", array = arraynms[i]) != 0
-#                pr <- getArrayData(BLData, what = "ProbeID", array = arraynms[i])[sel]
-#                finten <- getArrayData(BLData, what = what, log = log, array = arraynms[i])[sel]
-#                nasinf <- !is.finite(finten) | is.na(finten)
-#                pr <- pr[!nasinf]
-#                finten <- finten[!nasinf]
-#            }
-#            else if ((imagesPerArray == 2) && (j < len)) {
-#                sel1 <- getArrayData(BLData, what = "ProbeID", array = arraynms[j]) != 0
-#                sel2 <- getArrayData(BLData, what = "ProbeID", array = arraynms[j + 1]) != 0
-#                pr <- append(getArrayData(BLData, what = "ProbeID", array = arraynms[j])[sel1],
-#                             getArrayData(BLData, what = "ProbeID", array = arraynms[j + 1])[sel2])
-#                finten <- append(getArrayData(BLData, what = what, log = log, array = arraynms[j])[sel1],
-#                                 getArrayData(BLData, what = what, log = log, array = arraynms[j + 1])[sel2])
-#                nasinf <- !is.finite(finten) | is.na(finten)
-#                pr <- pr[!nasinf]
-#                finten <- finten[!nasinf]
-#            }
-#        }
-#        if (whatelse == "R") {
-#            rownames(G.kmd) <- rownames(G.ns) <- rownames(R.kmd) <- rownames(R.ns) <- probes
-#            res <- list(G = list("dist" = G.kmd, "n" = G.ns),
-#                        R = list("dist" = R.kmd, "n" = R.ns))
-#        }
-#        else {
-#            rownames(G.kmd) <- rownames(G.ns) <- probes
-#            res <- list("dist" = G.kmd, "n" = G.ns)
-#        }
-#        return(res)
-#    })
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
+###############################################################################
+## Methods for KolmogorovMinDist: computation of minimum Kolmogorov distance
+###############################################################################
+setMethod("KolmogorovMinDist", signature(x = "matrix",
+                                         D = "Norm"),
+    function(x, D, mad0 = 1e-4){
+        stopifnot(is.numeric(x))
+        ksdist <- function(pars, x){
+            if(any(ina <- is.na(x))) x <- x[!ina]
+            y <- sort(x, FALSE)
+            p.y <- pnorm(y, mean = pars[1], sd = pars[2])
+            n <- length(y)
+            max(p.y - (seq_len(n)-1)/n, seq_len(n)/n - p.y)
+        }
+        Med <- rowMedians(x, na.rm = TRUE)
+        Mad <- pmax(rowMedians(abs(x-Med), na.rm = TRUE)/qnorm(0.75), mad0)
+        startPars <- cbind(Med, Mad)
+        M <- nrow(x)
+        n <- ncol(x)
+        res <- matrix(NA, nrow = M, ncol = 2)
+        for(i in seq_len(M)){
+            temp <- try(optim(par = startPars[i,], fn = ksdist, x = x[i,], 
+                              method = "L-BFGS-B", lower = c(-Inf, 1e-15), 
+                              upper = c(Inf, Inf))$value, silent = TRUE)
+            if(!inherits(temp, "try-error")){
+                res[i,1] <- temp
+                res[i,2] <- n
+            }
+        }
+        list("dist" = res[,1], "n" = res[,2])
+    })
+
+setMethod("KolmogorovMinDist", signature(x = "AffyBatch",
+                                         D = "AbscontDistribution"),
+    function(x, D, bg.correct = TRUE, pmcorrect = TRUE, verbose = TRUE){
+        if(bg.correct){
+            if(verbose) cat("Background correcting ...")
+            x <- affy::bg.correct.mas(x, griddim = 16)
+            if(verbose) cat(" done.\n")
+        }
+        n <- length(x)
+        ids <- featureNames(x)
+        m <- length(ids)
+
+        CDFINFO <- getCdfInfo(x)
+        INDEX <- sapply(ids, get, envir = CDFINFO)
+        NROW <- unlist(lapply(INDEX, nrow))
+        nr <- as.integer(names(table(NROW)))
+
+        intensData <- intensity(x)
+        kmd <- numeric(m*n)
+        ns <- numeric(m*n)
+        if(pmcorrect){
+            if(verbose) cat("Calculating PM/MM ...")
+            div <- function(INDEX, x){
+                l.pm <- INDEX[,1]
+                if(ncol(INDEX) == 2)
+                    l.mm <- INDEX[,2]
+                else
+                    l.mm <- integer()
+                x[l.pm, , drop = FALSE]/x[l.mm, , drop = FALSE]
+            }
+            res <- lapply(INDEX, div, x = intensData)
+            if(verbose) cat(" done.\n")
+        }else{
+            if(verbose) cat("Extract PM data ...")
+            pm.only <- function(INDEX, x){
+                l.pm <- INDEX[,1]
+                x[l.pm, , drop = FALSE]
+            }
+            res <- lapply(INDEX, pm.only, x = intensData)
+            if(verbose) cat(" done.\n")
+        }
+        if(verbose) cat("Computing minimum Kolmogorov distance ...")
+        for(k in nr){
+            ind <- which(NROW == k)
+            temp <- matrix(do.call(rbind, res[ind]), nrow = k)
+            ind1 <-  as.vector(sapply(seq_len(n)-1, function(x, ind, m){ ind + x*m }, ind = ind, m = m))
+            temp.res <- KolmogorovMinDist(log2(t(temp)), D = D)
+            kmd[ind1] <- temp.res[["dist"]]
+            ns[ind1] <- temp.res[["n"]]
+        }
+        if(verbose) cat(" done.\n")
+        kmd.mat <- matrix(kmd, nrow = m)
+        ns.mat <- matrix(ns, nrow = m)
+        dimnames(kmd.mat) <- list(ids, sampleNames(x))
+        dimnames(ns.mat) <- list(ids, sampleNames(x))
+        list("dist" = kmd.mat, "n" = ns.mat)
+    })
+
+setMethod("KolmogorovMinDist", signature(x = "beadLevelData",
+                                         D = "AbscontDistribution"),
+    function(x, D, log = FALSE, what = "Grn", probes = NULL, arrays = NULL){
+        BLData <- x
+        arraynms <- sectionNames(BLData)
+        if(!is.null(arrays) && !is.character(arrays)) arraynms <- arraynms[arrays]
+        if(is.character(arrays)) arraynms <- which(arraynms %in% arrays)
+        len <- length(arraynms)
+
+        tmp <- BLData[[1]]
+        what <- match.arg(what, colnames(tmp))
+        if(is.na(what))
+          stop("Could not find bead data of type ", what, "\n The available choices are: ",
+               paste(colnames(tmp)), "\n")
+
+        pr <- getBeadData(BLData, what = "ProbeID", array = 1)
+        finten <- getBeadData(BLData, what = what, array = 1)
+        if(log) finten <- log2(finten)
+        nasinf <- !is.finite(finten) | is.na(finten)
+        finten <- finten[!nasinf]
+
+        if(is.null(probes)) probes <- sort(unique(pr))
+        probes <- probes[probes > 0 & !is.na(probes)]
+        noprobes <- length(probes)
+        pr <- pr[!nasinf]
+        G.kmd <- matrix(0, nrow = noprobes, ncol = len)
+        G.ns <- matrix(0, nrow = noprobes, ncol = len)
+        colnames(G.kmd) <- colnames(G.ns) <- arraynms
+
+        i <- 1
+        while (i <= len) {
+            probeIDs <- as.integer(pr)
+            start <- 0
+            G.res <- kmdBeadLevel(x = finten, D = D, probeIDs = probeIDs, probes = probes)
+            G.kmd[,i] <- G.res[["dist"]]
+            G.ns[,i] <- G.res[["n"]]
+
+            i <- i + 1
+            rm(probeIDs)
+            gc()
+
+            pr <- getBeadData(BLData, what = "ProbeID", array = i)
+            finten <- getBeadData(BLData, what = what, array = i)
+            if(log) finten <- log2(finten)
+            nasinf <- !is.finite(finten) | is.na(finten)
+            pr <- pr[!nasinf]
+            finten <- finten[!nasinf]
+        }
+        rownames(G.kmd) <- rownames(G.ns) <- probes
+        res <- list("dist" = G.kmd, "n" = G.ns)
+        return(res)
+    })
+
+kmdBeadLevel <- function(x, D, probeIDs, probes){
+    comIDs <- intersect(probeIDs, probes)
+    x <- x[probeIDs %in% comIDs]
+    probeIDs <- probeIDs[probeIDs %in% comIDs]
+    noBeads <- as.vector(table(probeIDs))
+    noBeads.uni <- as.integer(names(table(noBeads)))
+    probes1 <- comIDs
+    len1 <- length(probes1)
+    kmd1 <- numeric(len1)
+    ns1 <- numeric(len1)
+    for(i in seq(along = noBeads.uni)){
+        index <- noBeads == noBeads.uni[i]
+        IDs <- probes1[index]
+        if(noBeads.uni[i] == 1){
+            kmd1[index] <- 0.5
+        }else{
+            temp <- matrix(x[probeIDs %in% IDs], ncol = noBeads.uni[i], byrow = TRUE)
+            res <- KolmogorovMinDist(temp, D = D)
+            kmd1[index] <- res[["dist"]]
+            ns1[index] <- res[["n"]]
+        }
+    }
+    len <- length(probes)
+    kmd <- numeric(len)
+    ns <- numeric(len)
+    nas <- !(probes %in% comIDs)
+    kmd[nas] <- NA
+    kmd[!nas] <- kmd1
+    ns[nas] <- NA
+    ns[!nas] <- ns1
+
+    list("dist" = kmd, "n" = ns)
+}
+
+## leftover in trunk from revision 824 ....
+
+# setMethod("KolmogorovMinDist", signature(x = "BeadLevelList",
+#                                          D = "AbscontDistribution"),
+#     function(x, D, log = FALSE, imagesPerArray = 1, what = "G", probes = NULL, arrays = NULL){
+#         BLData <- x
+#         arraynms <- arrayNames(BLData)
+#         if(!is.null(arrays) && !is.character(arrays)) arraynms <- arraynms[arrays]
+#         if(is.character(arrays)) arraynms <- which(arraynms %in% arrays)
+#         len <- length(arraynms)
+#         what <- match.arg(what, c("G", "R", "RG", "M", "A", "beta"))
+#         whatelse <- ""
+#         if(what == "RG"){
+#             if(BLData at arrayInfo$channels == "two"){
+#                 what <- "G"
+#                 whatelse <- "R"
+#             }else{
+#                 stop("Need two-channel data to calculate summary R and G values")
+#             }
+#         }
+#         if(imagesPerArray == 1){
+#             sel <- getArrayData(BLData, what = "ProbeID", array = arraynms[1]) != 0
+#             pr <- getArrayData(BLData, what = "ProbeID", array = arraynms[1])[sel]
+#             finten <- getArrayData(BLData, what = what, log = log, array = arraynms[1])[sel]
+#             nasinf <- !is.finite(finten) | is.na(finten)
+#             finten <- finten[!nasinf]
+#        }
+#        else if(imagesPerArray == 2){
+#            if(length(arraynms)%%2 != 0)
+#                stop("Need an even number of arrays when 'imagesPerArray=2'")
+#            arrayord <- order(arraynms)
+#            arraynms <- arraynms[arrayord]
+#            tmp <- unlist(strsplit(arraynms, "_"))
+#            chipnums <- tmp[seq(1, length(tmp), by = 3)]
+#            pos <- tmp[seq(2, length(tmp), by = 3)]
+#            stripnum <- as.numeric(tmp[seq(3, length(tmp), by = 3)])
+#            check <- ((chipnums[seq(1, len, by = 2)] == chipnums[seq(2, len, by = 2)])
+#                      & (pos[seq(1, len, by = 2)] == pos[seq(2, len, by = 2)])
+#                      & (stripnum[seq(1, len, by = 2)] == stripnum[seq(2, len, by = 2)] - 1))
+#            if(sum(check) != length(check)) stop("Missing arrays")
+#            sel1 <- getArrayData(BLData, what = "ProbeID", array = arraynms[1]) != 0
+#            sel2 <- getArrayData(BLData, what = "ProbeID", array = arraynms[2]) != 0
+#            pr <- append(getArrayData(BLData, what = "ProbeID", array = arraynms[1])[sel1],
+#                        getArrayData(BLData, what = "ProbeID", array = arraynms[2])[sel2])
+#            finten <- append(getArrayData(BLData, what = what, log = log, array = arraynms[1])[sel1],
+#                            getArrayData(BLData, what = what, log = log, array = arraynms[2])[sel2])
+#            nasinf <- !is.finite(finten) | is.na(finten)
+#            finten <- finten[!nasinf]
+#        }else{
+#            stop("You can only specify 1 or 2 images per array")
+#        }
+#        if(is.null(probes)) probes <- sort(unique(pr))
+#        probes <- probes[probes > 0 & !is.na(probes)]
+#        noprobes <- length(probes)
+#        pr <- pr[!nasinf]
+#        if (imagesPerArray == 1) {
+#            G.kmd <- matrix(0, nrow = noprobes, ncol = len)
+#            G.ns <- matrix(0, nrow = noprobes, ncol = len)
+#            colnames(G.kmd) <- colnames(G.ns) <- arraynms
+#            if (BLData at arrayInfo$channels == "two" && !is.null(BLData[[arraynms[1]]]$R) && whatelse == "R")
+#                R.kmd <- R.ns <- G.kmd
+#            else R.kmd <- R.ns <- NULL
+#        }
+#        else if (imagesPerArray == 2) {
+#            G.kmd <- matrix(0, nrow = noprobes, ncol = (len/2))
+#            G.ns <- matrix(0, nrow = noprobes, ncol = (len/2))
+#            colnames(G.kmd) <- colnames(G.ns) <- arraynms[seq(1, len, by = 2)]
+#            if (BLData at arrayInfo$channels == "two" && !is.null(BLData[[arraynms[1]]]$R) && whatelse == "R")
+#                R.kmd <- R.ns <- G.kmd
+#            else R.kmd <- R.ns <- NULL
+#        }
+#        i <- j <- 1
+#        while (j <= len) {
+#            probeIDs <- as.integer(pr)
+#            start <- 0
+#            G.res <- kmdBeadLevel(x = finten, D = D, probeIDs = probeIDs, probes = probes)
+#            G.kmd[,i] <- G.res[["dist"]]
+#            G.ns[,i] <- G.res[["n"]]
+#            if (BLData at arrayInfo$channels == "two" && !is.null(BLData[[arraynms[i]]]$R) && whatelse == "R") {
+#                if (imagesPerArray == 1) {
+#                    finten <- getArrayData(BLData, what = whatelse, log = log, array = arraynms[i])[sel]
+#                    nasinf <- !is.finite(finten) | is.na(finten)
+#                    finten <- finten[!nasinf]
+#                    binten <- rep(0, length(finten))
+#                }
+#                else if (imagesPerArray == 2) {
+#                    finten <- append(getArrayData(BLData, what = whatelse, log = log, array = arraynms[j])[sel1],
+#                                    getArrayData(BLData, what = whatelse, log = log, array = arraynms[j + 1])[sel2])
+#                    nasinf <- !is.finite(finten) | is.na(finten)
+#                    finten <- finten[!nasinf]
+#                    binten <- rep(0, length(finten))
+#                }
+#                start <- 0
+#                R.res[,i] <- kmdBeadLevel(x = finten, D = D, probeIDs = probeIDs, probes = probes)
+#                R.kmd[,i] <- R.res[["dist"]]
+#                R.ns[,i] <- R.res[["n"]]
+#            }
+#            j <- j + imagesPerArray
+#            i <- i + 1
+#            rm(probeIDs)
+#            gc()
+#            if ((imagesPerArray == 1) && (i <= len)) {
+#                sel <- getArrayData(BLData, what = "ProbeID", array = arraynms[i]) != 0
+#                pr <- getArrayData(BLData, what = "ProbeID", array = arraynms[i])[sel]
+#                finten <- getArrayData(BLData, what = what, log = log, array = arraynms[i])[sel]
+#                nasinf <- !is.finite(finten) | is.na(finten)
+#                pr <- pr[!nasinf]
+#                finten <- finten[!nasinf]
+#            }
+#            else if ((imagesPerArray == 2) && (j < len)) {
+#                sel1 <- getArrayData(BLData, what = "ProbeID", array = arraynms[j]) != 0
+#                sel2 <- getArrayData(BLData, what = "ProbeID", array = arraynms[j + 1]) != 0
+#                pr <- append(getArrayData(BLData, what = "ProbeID", array = arraynms[j])[sel1],
+#                             getArrayData(BLData, what = "ProbeID", array = arraynms[j + 1])[sel2])
+#                finten <- append(getArrayData(BLData, what = what, log = log, array = arraynms[j])[sel1],
+#                                 getArrayData(BLData, what = what, log = log, array = arraynms[j + 1])[sel2])
+#                nasinf <- !is.finite(finten) | is.na(finten)
+#                pr <- pr[!nasinf]
+#                finten <- finten[!nasinf]
+#            }
+#        }
+#        if (whatelse == "R") {
+#            rownames(G.kmd) <- rownames(G.ns) <- rownames(R.kmd) <- rownames(R.ns) <- probes
+#            res <- list(G = list("dist" = G.kmd, "n" = G.ns),
+#                        R = list("dist" = R.kmd, "n" = R.ns))
+#        }
+#        else {
+#            rownames(G.kmd) <- rownames(G.ns) <- probes
+#            res <- list("dist" = G.kmd, "n" = G.ns)
+#        }
+#        return(res)
+#    })
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+

Modified: pkg/RobLoxBioC/man/robloxbioc.Rd
===================================================================
--- pkg/RobLoxBioC/man/robloxbioc.Rd	2025-01-12 14:32:59 UTC (rev 1332)
+++ pkg/RobLoxBioC/man/robloxbioc.Rd	2025-01-19 16:58:14 UTC (rev 1333)
@@ -1,187 +1,186 @@
-\name{robloxbioc}
-\alias{robloxbioc}
-\alias{robloxbioc-methods}
-\alias{robloxbioc,matrix-method}
-\alias{robloxbioc,AffyBatch-method}
-\alias{robloxbioc,beadLevelData-method}
-
-\title{Generic Function for Preprocessing Biological Data}
-\description{
-  Generic function for preprocessing biological data using optimally robust
-  (rmx) estimators; confer Rieder (1994), Kohl (2005), Rieder et al (2008).
-}
-\usage{
-robloxbioc(x, ...)
-
-\S4method{robloxbioc}{matrix}(x, eps = NULL, eps.lower = 0, eps.upper = 0.05, steps = 3L, 
-           fsCor = TRUE, mad0 = 1e-4)
-
-\S4method{robloxbioc}{AffyBatch}(x, bg.correct = TRUE, pmcorrect = TRUE, normalize = FALSE, 
-           add.constant = 32, verbose = TRUE, eps = NULL, 
-           eps.lower = 0, eps.upper = 0.05, steps = 3L, fsCor = TRUE, 
-           mad0 = 1e-4, contrast.tau = 0.03, scale.tau = 10, 
-           delta = 2^(-20), sc = 500)
-
-\S4method{robloxbioc}{beadLevelData}(x, channelList = list(greenChannel), probeIDs = NULL, 
-           useSampleFac = FALSE, sampleFac = NULL, weightNames = "wts", 
-           removeUnMappedProbes = TRUE, eps = NULL, eps.lower = 0, 
-           eps.upper = 0.05, steps = 3L, fsCor = TRUE, mad0 = 1e-4)
-}
-\arguments{
-  \item{x}{ biological data. }
-  \item{\dots}{ additional parameters. }
-  \item{eps}{ positive real (0 < \code{eps} <= 0.5): amount of gross errors. 
-        See details below. }
-  \item{eps.lower}{ positive real (0 <= \code{eps.lower} <= \code{eps.upper}): 
-        lower bound for the amount of gross errors. See details below. }
-  \item{eps.upper}{ positive real (\code{eps.lower} <= \code{eps.upper} <= 0.5): 
-        upper bound for the amount of gross errors. See details below. }
-  \item{steps}{ positive integer. k-step is used to compute the optimally robust estimator. }
-  \item{fsCor}{ logical: perform finite-sample correction. See function \code{\link[RobLox]{finiteSampleCorrection}}. }
-  \item{mad0}{ scale estimate used if computed MAD is equal to zero}
-  \item{bg.correct}{ if \code{TRUE} MAS 5.0 background correction is performed;
-    confer \code{\link[affy:bgc]{bg.correct.mas}}. }
-  \item{pmcorrect}{ method used for PM correction; \code{TRUE} calls an algorithm which is 
-    comparable to the algorithm of MAS 5.0; confer \code{\link[affy:pmcorrect]{pmcorrect.mas}}. 
-    If \code{FALSE} only the PM intensities are used. }
-  \item{normalize}{ logical: if \code{TRUE}, Affymetrix scale normalization is performed. }
-  \item{add.constant}{ constant added to the MAS 5.0 expression values before the normalization
-    step. Improves the variance of the measure one no longer devides by numbers close to 0
-    when computing fold-changes. }
-  \item{verbose}{ logical: if \code{TRUE}, some messages are printed. }
-  \item{contrast.tau}{ a number denoting the contrast tau parameter; confer the MAS 5.0 
-    PM correction algorithm. }
-  \item{scale.tau}{ a number denoting the scale tau parameter; confer the MAS 5.0 
-    PM correction algorithm. }
-  \item{delta}{ a number denoting the delta parameter; confer the MAS 5.0 
-    PM correction algorithm. }
-  \item{sc}{ value at which all arrays will be scaled to. }
-  \item{channelList}{ List of objects of class illuminaChannel that defines the
-    summarisation to be performed where in our setup only the slots \code{transFun}
-    and \code{name} have an effect on the computations. Setting the slots 
-    \code{outlierFun}, \code{exprFun}, and \code{varFun} has no effect. In any 
-    case rmx estimators are applied. }
-  \item{probeIDs}{ Vector of ArrayAddressIDs to be included in the summarized object. 
-    The default is to summarize all probes. }
-  \item{useSampleFac}{if \code{TRUE} sections belonging to the same biological sample
-    will be combined. The default is to summarize each section separately.}
[TRUNCATED]

To get the complete diff run:
    svnlook diff /svnroot/robast -r 1333


More information about the Robast-commits mailing list