| fd653332 |
setMethod("show", "CNSetLM", function(object){
|
| e40ceb44 |
callNextMethod(object)
cat("lM: ", length(lM(object)), " elements \n")
print(names(lM(object)))
})
|
| d5bc779e |
|
| fd653332 |
setMethod("[", "CNSetLM", function(x, i, j, ..., drop=FALSE){
|
| e40ceb44 |
x <- callNextMethod(x, i, j, ..., drop=drop)
if(!missing(i)){
if(class(lM(x)) == "ffdf"){
lM(x) <- lapply(physical(lM(x)), function(x, i){open(x); x[i, ]}, i=i)
} else {
lM(x) <- lapply(lM(x), function(x, i) x[i, , drop=FALSE], i=i)
}
}
x
})
|
| d5bc779e |
|
| ea66649b |
|
| fd653332 |
setMethod("lM", "CNSetLM", function(object) object@lM)
|
| 4c48a3d8 |
setReplaceMethod("lM", c("CNSetLM", "ffdf"), function(object, value){
object@lM <- value
object
})
setReplaceMethod("lM", c("CNSetLM", "list"), function(object, value){
|
| e40ceb44 |
object@lM <- value
object
})
|
| d5bc779e |
setMethod("open", "CNSetLM", function(con,...){
callNextMethod(con,...)
lapply(physical(lM(con)), open)
})
|
| c00f715a |
setAs("SnpSuperSet", "CNSetLM", function(from, to){
|
| d5bc779e |
stopifnot("batch" %in% varLabels(protocolData(from)))
|
| c00f715a |
cnSet <- new("CNSetLM",
alleleA=A(from),
alleleB=B(from),
call=snpCall(from),
callProbability=snpCallProbability(from),
CA=initializeBigMatrix("CA", nrow(from), ncol(from)),
CB=initializeBigMatrix("CB", nrow(from), ncol(from)),
annotation=annotation(from),
featureData=featureData(from),
experimentData=experimentData(from),
|
| d5bc779e |
protocolData=protocolData(from),
|
| c00f715a |
phenoData=phenoData(from))
|
| d5bc779e |
lM(cnSet) <- initializeParamObject(list(featureNames(cnSet), unique(protocolData(from)$batch)))
|
| c00f715a |
return(cnSet)
})
|
| 5fcfed7d |
|
| 1c2a3f47 |
setMethod("computeCopynumber", "CNSet",
function(object,
MIN.OBS,
DF.PRIOR,
bias.adj,
prior.prob,
seed,
verbose,
GT.CONF.THR,
PHI.THR,
nHOM.THR,
MIN.NU,
MIN.PHI,
|
| c00f715a |
THR.NU.PHI,
|
| 1c2a3f47 |
thresholdCopynumber){
|
| 5fcfed7d |
## to do the bias adjustment, initial estimates of the parameters are needed
## The initial estimates are gotten by running computeCopynumber with cnOptions[["bias.adj"]]=FALSE
|
| 1c2a3f47 |
cnOptions <- list(
MIN.OBS=MIN.OBS,
|
| c00f715a |
DF.PRIOR=DF.PRIOR,
|
| 1c2a3f47 |
bias.adj=bias.adj,
prior.prob=prior.prob,
seed=seed,
verbose=verbose,
|
| c00f715a |
GT.CONF.THR=GT.CONF.THR,
|
| 1c2a3f47 |
PHI.THR=PHI.THR,
nHOM.THR=nHOM.THR,
MIN.NU=MIN.NU,
MIN.PHI=MIN.PHI,
THR.NU.PHI=THR.NU.PHI,
|
| c00f715a |
thresholdCopynumber=thresholdCopynumber)
|
| 5fcfed7d |
bias.adj <- cnOptions[["bias.adj"]]
if(bias.adj & all(is.na(CA(object)))){
cnOptions[["bias.adj"]] <- FALSE
}
|
| 7f14b4e4 |
object <- computeCopynumber.CNSet(object, cnOptions)
|
| 5fcfed7d |
if(bias.adj & !cnOptions[["bias.adj"]]){
## Do a second iteration with bias adjustment
cnOptions[["bias.adj"]] <- TRUE
object <- computeCopynumber.CNSet(object, cnOptions)
}
object
|
| 96d10c6e |
})
|
| 7e1d14dd |
|
| 5fcfed7d |
setMethod("copyNumber", "CNSet", function(object){
|
| 74345385 |
message("This accessor will be deprecated.")
|
| 5fcfed7d |
I <- isSnp(object)
|
| c87212d7 |
ffIsLoaded <- inherits(CA(object), "ff")
|
| 5fcfed7d |
CA <- CA(object)
CB <- CB(object)
|
| c2b6e3a8 |
if(ffIsLoaded){
open(CA)
open(CB)
|
| 86556c4b |
CA <- as.matrix(CA[,])
CB <- as.matrix(CB[,])
|
| c2b6e3a8 |
}
|
| 5fcfed7d |
CN <- CA + CB
##For nonpolymorphic probes, CA is the total copy number
CN[!I, ] <- CA(object)[!I, ]
CN
})
|
| 96d10c6e |
|
| 5fcfed7d |
setMethod("ellipse", "CNSet", function(x, copynumber, batch, ...){
ellipse.CNSet(x, copynumber, batch, ...)
|
| 96d10c6e |
})
|
| 5fcfed7d |
##setMethod("ellipse", "CNSet", function(x, copynumber, ...){
ellipse.CNSet <- function(x, copynumber, batch, ...){
if(nrow(x) > 1) stop("only 1 snp at a time")
##batch <- unique(x$batch)
if(missing(batch)){
stop("must specify batch")
|
| 96d10c6e |
}
|
| 5fcfed7d |
if(length(batch) > 1) stop("batch variable not unique")
nuA <- getParam(x, "nuA", batch)
nuB <- getParam(x, "nuB", batch)
phiA <- getParam(x, "phiA", batch)
phiB <- getParam(x, "phiB", batch)
tau2A <- getParam(x, "tau2A", batch)
tau2B <- getParam(x, "tau2B", batch)
sig2A <- getParam(x, "sig2A", batch)
sig2B <- getParam(x, "sig2B", batch)
corrA.BB <- getParam(x, "corrA.BB", batch)
corrB.AA <- getParam(x, "corrB.AA", batch)
corr <- getParam(x, "corr", batch)
for(CN in copynumber){
for(CA in 0:CN){
CB <- CN-CA
A.scale <- sqrt(tau2A*(CA==0) + sig2A*(CA > 0))
B.scale <- sqrt(tau2B*(CB==0) + sig2B*(CB > 0))
scale <- c(A.scale, B.scale)
if(CA == 0 & CB > 0) rho <- corrA.BB
if(CA > 0 & CB == 0) rho <- corrB.AA
if(CA > 0 & CB > 0) rho <- corr
if(CA == 0 & CB == 0) rho <- 0
lines(ellipse(x=rho, centre=c(log2(nuA+CA*phiA),
log2(nuB+CB*phiB)),
scale=scale), ...)
|
| 96d10c6e |
}
}
}
|
| 7f14b4e4 |
setMethod("totalCopyNumber",
signature=signature(object="CNSet", ...),
function(object, ...){
is.ff <- is(CA(object), "ff") | is(CA(object), "ffdf")
dotArgs <- names(list(...))
missing.i <- "i" %in% dotArgs
missing.j <- "j" %in% dotArgs
if(missing.i & missing.j){
if(is.ff) stop("Must specify i and/or j for ff objects")
}
if(!missing.i) {
i <- dotArgs[["i"]]
snp.index <- intersect(i, which(isSnp(object)))
##which rows in the return matrix are snps
snp.index2 <- which(isSnp(object)[i])
} else {
i <- 1:nrow(object)
snp.index2 <- snp.index <- which(isSnp(object))
}
if(!missing.j){
j <- dotArgs[["j"]]
} else j <- 1:ncol(object)
cn.total <- as.matrix(CA(object)[i, j])
if(length(snp.index) > 0){
cb <- as.matrix(CB(object)[snp.index, j])
## snps <- (1:nrow(cn.total))[i %in% snp.index]
cn.total[snp.index2, ] <- cn.total[snp.index, ] + cb
}
cn.total <- cn.total/100
return(cn.total)
## if(missing.i){
#### if(missing.i & !missing.j){
## snp.index <- which(isSnp(object))
## cn.total <- as.matrix(CA(object)[, j])
## if(length(snp.index) > 0){
## cb <- as.matrix(CB(object)[snp.index, j])
#### snps <- (1:nrow(cn.total))[i %in% snp.index]
## cn.total[snp.index, ] <- cn.total[snp.index, ] + cb
## }
## } else{
## snp.index <- intersect(which(isSnp(object)), i)
## cn.total <- as.matrix(CA(object)[i, ])
## if(length(snp.index) > 0){
## cb <- as.matrix(CB(object)[snp.index, ])
## snps <- (1:nrow(cn.total))[i %in% snp.index]
## cn.total[snps, ] <- cn.total[snps, ] + cb
## }
## }
## if(!missing(i) & missing(j)){
## snp.index <- intersect(which(isSnp(object)), i)
## cn.total <- as.matrix(CA(object)[i, ])
## if(length(snp.index) > 0){
## cb <- as.matrix(CB(object)[snp.index, ])
## snps <- (1:nrow(cn.total))[i %in% snp.index]
## cn.total[snps, ] <- cn.total[snps, ] + cb
## }
## }
## if(!missing(i) & !missing(j)){
## snp.index <- intersect(which(isSnp(object)), i)
## cn.total <- as.matrix(CA(object)[i, j])
## if(length(snp.index) > 0){
## cb <- as.matrix(CB(object)[snp.index, j])
## snps <- (1:nrow(cn.total))[i %in% snp.index]
## cn.total[snps, ] <- cn.total[snps, ] + cb
## }
## }
## cn.total <- cn.total/100
## dimnames(cn.total) <- NULL
## return(cn.total)
})
|
