[Robast-commits] r1046 - branches/robast-1.2/pkg/RobLoxBioC/R branches/robast-1.2/pkg/RobLoxBioC/man pkg/RobLoxBioC/R pkg/RobLoxBioC/man
noreply at r-forge.r-project.org
noreply at r-forge.r-project.org
Tue Jul 24 12:45:01 CEST 2018
Author: ruckdeschel
Date: 2018-07-24 12:45:01 +0200 (Tue, 24 Jul 2018)
New Revision: 1046
Modified:
branches/robast-1.2/pkg/RobLoxBioC/R/KolmogorovMinDist.R
branches/robast-1.2/pkg/RobLoxBioC/man/KolmogorovMinDist.Rd
pkg/RobLoxBioC/R/KolmogorovMinDist.R
pkg/RobLoxBioC/man/KolmogorovMinDist.Rd
Log:
[RobLoxBioC] cleaned leftovers from mergen in trunk and branch 1.2
Modified: branches/robast-1.2/pkg/RobLoxBioC/R/KolmogorovMinDist.R
===================================================================
--- branches/robast-1.2/pkg/RobLoxBioC/R/KolmogorovMinDist.R 2018-07-24 09:54:45 UTC (rev 1045)
+++ branches/robast-1.2/pkg/RobLoxBioC/R/KolmogorovMinDist.R 2018-07-24 10:45:01 UTC (rev 1046)
@@ -88,136 +88,59 @@
list("dist" = kmd.mat, "n" = ns.mat)
})
-setMethod("KolmogorovMinDist", signature(x = "BeadLevelList",
+setMethod("KolmogorovMinDist", signature(x = "beadLevelData",
D = "AbscontDistribution"),
- function(x, D, log = FALSE, imagesPerArray = 1, what = "G", probes = NULL, arrays = NULL){
+ function(x, D, log = FALSE, what = "Grn", probes = NULL, arrays = NULL){
BLData <- x
- arraynms <- arrayNames(BLData)
+ 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)
- 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")
- }
+
+ 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]
- 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) {
+ 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"]]
- 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]
- }
+
+ 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]
}
- 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)
- }
+ 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]
@@ -239,7 +162,7 @@
kmd1[index] <- res[["dist"]]
ns1[index] <- res[["n"]]
}
- }
+ }
len <- length(probes)
kmd <- numeric(len)
ns <- numeric(len)
@@ -251,3 +174,170 @@
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: branches/robast-1.2/pkg/RobLoxBioC/man/KolmogorovMinDist.Rd
===================================================================
--- branches/robast-1.2/pkg/RobLoxBioC/man/KolmogorovMinDist.Rd 2018-07-24 09:54:45 UTC (rev 1045)
+++ branches/robast-1.2/pkg/RobLoxBioC/man/KolmogorovMinDist.Rd 2018-07-24 10:45:01 UTC (rev 1046)
@@ -3,7 +3,7 @@
\alias{KolmogorovMinDist-methods}
\alias{KolmogorovMinDist,matrix,Norm-method}
\alias{KolmogorovMinDist,AffyBatch,AbscontDistribution-method}
-\alias{KolmogorovMinDist,BeadLevelList,AbscontDistribution-method}
+\alias{KolmogorovMinDist,beadLevelData,AbscontDistribution-method}
\title{Generic Function for Computing Minimum Kolmogorov Distance for Biological Data}
\description{
@@ -17,7 +17,7 @@
\S4method{KolmogorovMinDist}{AffyBatch,AbscontDistribution}(x, D, bg.correct = TRUE, pmcorrect = TRUE,
verbose = TRUE)
-\S4method{KolmogorovMinDist}{BeadLevelList,AbscontDistribution}(x, D, log = FALSE, imagesPerArray = 1, what = "G",
+\S4method{KolmogorovMinDist}{beadLevelData,AbscontDistribution}(x, D, log = FALSE, what = "Grn",
probes = NULL, arrays = NULL)
}
\arguments{
@@ -32,10 +32,8 @@
If \code{FALSE} only \code{log2(PM)} is used. }
\item{verbose}{ logical: if \code{TRUE}, some messages are printed. }
\item{log}{ if \code{TRUE}, then the log2 intensities for each bead-type are summarized. }
- \item{imagesPerArray}{ Specifies how many images (strips) there are per array.
- Normally 1 for a SAM and 1 or 2 for a BeadChip. The images (strips) from the same array
- will be combined so that each column in the output represents a sample. }
- \item{what}{ character string specifying which intensities/values to summarize. }
+ \item{what}{ character string specifying which intensities/values to summarize;
+ see \code{\link[beadarray]{getBeadData}}. }
\item{probes}{ Specify particular probes to summarize. If left \code{NULL} then all
the probes on the first array are used. }
\item{arrays}{ integer (scalar or vector) specifying the strips/arrays to summarize.
@@ -67,15 +65,16 @@
X <- matrix(rnorm(200, mean=ind*3, sd=(1-ind) + ind*9), nrow = 2)
KolmogorovMinDist(X, D = Norm())
-## using Affymetrix-Data
+## using Affymetrix data
data(SpikeIn)
probes <- log2(pm(SpikeIn))
(res <- KolmogorovMinDist(probes, Norm()))
boxplot(res$dist)
-\dontrun{
-## "Not run" just because of computation time
-require(affydata)
+\donttest{
+## \donttest because of check time
+## using Affymetrix data
+library(affydata)
data(Dilution)
res <- KolmogorovMinDist(Dilution[,1], Norm())
summary(res$dist)
@@ -86,14 +85,12 @@
uni.n <- min(res$n):max(res$n)
lines(uni.n, 1/(2*uni.n), col = "orange", lwd = 2)
legend("topright", legend = "minimal possible distance", fill = "orange")
-}
-## using Illumina-Data
-\dontrun{
-## "Not run" just because of computation time
-data(BLData)
-res <- KolmogorovMinDist(BLData, Norm(), arrays = 1)
-res1 <- KolmogorovMinDist(BLData, log = TRUE, Norm(), arrays = 1)
+## Illumina bead level data
+library(beadarrayExampleData)
+data(exampleBLData)
+res <- KolmogorovMinDist(exampleBLData, Norm(), arrays = 1)
+res1 <- KolmogorovMinDist(exampleBLData, Norm(), log = TRUE, arrays = 1)
summary(cbind(res$dist, res1$dist))
boxplot(list(res$dist, res1$dist), names = c("raw", "log-raw"))
sort(unique(res1$n))
Modified: pkg/RobLoxBioC/R/KolmogorovMinDist.R
===================================================================
--- pkg/RobLoxBioC/R/KolmogorovMinDist.R 2018-07-24 09:54:45 UTC (rev 1045)
+++ pkg/RobLoxBioC/R/KolmogorovMinDist.R 2018-07-24 10:45:01 UTC (rev 1046)
@@ -88,136 +88,59 @@
list("dist" = kmd.mat, "n" = ns.mat)
})
-setMethod("KolmogorovMinDist", signature(x = "BeadLevelList",
+setMethod("KolmogorovMinDist", signature(x = "beadLevelData",
D = "AbscontDistribution"),
- function(x, D, log = FALSE, imagesPerArray = 1, what = "G", probes = NULL, arrays = NULL){
+ function(x, D, log = FALSE, what = "Grn", probes = NULL, arrays = NULL){
BLData <- x
- arraynms <- arrayNames(BLData)
+ 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)
- 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")
- }
+
+ 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]
- 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) {
+ 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"]]
- 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]
- }
+
+ 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]
}
- 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)
- }
+ 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]
@@ -239,7 +162,7 @@
kmd1[index] <- res[["dist"]]
ns1[index] <- res[["n"]]
}
- }
+ }
len <- length(probes)
kmd <- numeric(len)
ns <- numeric(len)
@@ -251,3 +174,170 @@
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)
[TRUNCATED]
To get the complete diff run:
svnlook diff /svnroot/robast -r 1046
More information about the Robast-commits
mailing list