git-svn-id: file:///home/git/hedgehog.fhcrc.org/bioconductor/trunk/madman/Rpacks/crlmm@49125 bc3139a8-67e5-0310-9ffc-ced21a209358
Rob Scharp authored on 30/08/2010 19:40:16 | ... | ... |
@@ -56,6 +56,7 @@ importFrom(mvtnorm, dmvnorm) |
| 56 | 56 |
|
| 57 | 57 |
importFrom(ellipse, ellipse) |
| 58 | 58 |
|
| 59 |
+<<<<<<< HEAD |
|
| 59 | 60 |
importFrom(ff, ffdf, physical.ff, physical.ffdf) |
| 60 | 61 |
importClassesFrom(oligoClasses, ffdf, ff_matrix) |
| 61 | 62 |
|
| ... | ... |
@@ -67,6 +68,19 @@ export(crlmm, |
| 67 | 68 |
ellipseCenters, |
| 68 | 69 |
genotype, |
| 69 | 70 |
readIdatFiles, |
| 71 |
+======= |
|
| 72 |
+importFrom(ff, ffdf) |
|
| 73 |
+ |
|
| 74 |
+exportClasses(CNSetLM, ffdf, list) |
|
| 75 |
+exportMethods(open, "[", show, lM) |
|
| 76 |
+export(crlmm, |
|
| 77 |
+ crlmmCopynumber, |
|
| 78 |
+ crlmmIllumina, |
|
| 79 |
+## ellipse, |
|
| 80 |
+ genotype, |
|
| 81 |
+ readIdatFiles, |
|
| 82 |
+ readIdatFiles2, |
|
| 83 |
+>>>>>>> Exporting list and ffdf classes so that lM<- assignment works |
|
| 70 | 84 |
snprma, |
| 71 | 85 |
snprma2, |
| 72 | 86 |
crlmm2, |
| ... | ... |
@@ -347,8 +347,143 @@ ACN.X <- function(object, allele, i, j){
|
| 347 | 347 |
if(allele=="B") acn[[k]] <- cb |
| 348 | 348 |
} |
| 349 | 349 |
} |
| 350 |
+<<<<<<< HEAD |
|
| 350 | 351 |
return(acn) |
| 351 | 352 |
} |
| 353 |
+======= |
|
| 354 |
+ x |
|
| 355 |
+}) |
|
| 356 |
+ |
|
| 357 |
+ |
|
| 358 |
+setMethod("lM", "CNSetLM", function(object) object@lM)
|
|
| 359 |
+setReplaceMethod("lM", c("CNSetLM", "list_or_ffdf"), function(object, value){
|
|
| 360 |
+ object@lM <- value |
|
| 361 |
+ object |
|
| 362 |
+}) |
|
| 363 |
+ |
|
| 364 |
+ |
|
| 365 |
+ |
|
| 366 |
+setMethod("open", "CNSetLM", function(con,...){
|
|
| 367 |
+ callNextMethod(con,...) |
|
| 368 |
+ lapply(physical(lM(con)), open) |
|
| 369 |
+}) |
|
| 370 |
+ |
|
| 371 |
+setAs("SnpSuperSet", "CNSetLM", function(from, to){
|
|
| 372 |
+ stopifnot("batch" %in% varLabels(protocolData(from)))
|
|
| 373 |
+ cnSet <- new("CNSetLM",
|
|
| 374 |
+ alleleA=A(from), |
|
| 375 |
+ alleleB=B(from), |
|
| 376 |
+ call=snpCall(from), |
|
| 377 |
+ callProbability=snpCallProbability(from), |
|
| 378 |
+ CA=initializeBigMatrix("CA", nrow(from), ncol(from)),
|
|
| 379 |
+ CB=initializeBigMatrix("CB", nrow(from), ncol(from)),
|
|
| 380 |
+ annotation=annotation(from), |
|
| 381 |
+ featureData=featureData(from), |
|
| 382 |
+ experimentData=experimentData(from), |
|
| 383 |
+ protocolData=protocolData(from), |
|
| 384 |
+ phenoData=phenoData(from)) |
|
| 385 |
+ lM(cnSet) <- initializeParamObject(list(featureNames(cnSet), unique(protocolData(from)$batch))) |
|
| 386 |
+ return(cnSet) |
|
| 387 |
+}) |
|
| 388 |
+ |
|
| 389 |
+setMethod("computeCopynumber", "CNSet",
|
|
| 390 |
+ function(object, |
|
| 391 |
+ MIN.OBS, |
|
| 392 |
+ DF.PRIOR, |
|
| 393 |
+ bias.adj, |
|
| 394 |
+ prior.prob, |
|
| 395 |
+ seed, |
|
| 396 |
+ verbose, |
|
| 397 |
+ GT.CONF.THR, |
|
| 398 |
+ PHI.THR, |
|
| 399 |
+ nHOM.THR, |
|
| 400 |
+ MIN.NU, |
|
| 401 |
+ MIN.PHI, |
|
| 402 |
+ THR.NU.PHI, |
|
| 403 |
+ thresholdCopynumber){
|
|
| 404 |
+ ## to do the bias adjustment, initial estimates of the parameters are needed |
|
| 405 |
+ ## The initial estimates are gotten by running computeCopynumber with cnOptions[["bias.adj"]]=FALSE |
|
| 406 |
+ cnOptions <- list( |
|
| 407 |
+ MIN.OBS=MIN.OBS, |
|
| 408 |
+ DF.PRIOR=DF.PRIOR, |
|
| 409 |
+ bias.adj=bias.adj, |
|
| 410 |
+ prior.prob=prior.prob, |
|
| 411 |
+ seed=seed, |
|
| 412 |
+ verbose=verbose, |
|
| 413 |
+ GT.CONF.THR=GT.CONF.THR, |
|
| 414 |
+ PHI.THR=PHI.THR, |
|
| 415 |
+ nHOM.THR=nHOM.THR, |
|
| 416 |
+ MIN.NU=MIN.NU, |
|
| 417 |
+ MIN.PHI=MIN.PHI, |
|
| 418 |
+ THR.NU.PHI=THR.NU.PHI, |
|
| 419 |
+ thresholdCopynumber=thresholdCopynumber) |
|
| 420 |
+ bias.adj <- cnOptions[["bias.adj"]] |
|
| 421 |
+ if(bias.adj & all(is.na(CA(object)))){
|
|
| 422 |
+ cnOptions[["bias.adj"]] <- FALSE |
|
| 423 |
+ } |
|
| 424 |
+ object <- computeCopynumber.CNSet(object, cnOptions) |
|
| 425 |
+ if(bias.adj & !cnOptions[["bias.adj"]]){
|
|
| 426 |
+ ## Do a second iteration with bias adjustment |
|
| 427 |
+ cnOptions[["bias.adj"]] <- TRUE |
|
| 428 |
+ object <- computeCopynumber.CNSet(object, cnOptions) |
|
| 429 |
+ } |
|
| 430 |
+ object |
|
| 431 |
+}) |
|
| 432 |
+setMethod("copyNumber", "CNSet", function(object){
|
|
| 433 |
+ I <- isSnp(object) |
|
| 434 |
+ ffIsLoaded <- class(calls(object))[[1]]=="ff" |
|
| 435 |
+ CA <- CA(object) |
|
| 436 |
+ CB <- CB(object) |
|
| 437 |
+ if(ffIsLoaded){
|
|
| 438 |
+ open(CA) |
|
| 439 |
+ open(CB) |
|
| 440 |
+ CA <- as.matrix(CA[,]) |
|
| 441 |
+ CB <- as.matrix(CB[,]) |
|
| 442 |
+ } |
|
| 443 |
+ CN <- CA + CB |
|
| 444 |
+ ##For nonpolymorphic probes, CA is the total copy number |
|
| 445 |
+ CN[!I, ] <- CA(object)[!I, ] |
|
| 446 |
+ CN <- CN/100 |
|
| 447 |
+ CN |
|
| 448 |
+}) |
|
| 449 |
+ |
|
| 450 |
+setMethod("totalCopyNumber", "CNSet", function(object, i, j){
|
|
| 451 |
+ if(missing(i) & missing(j)){
|
|
| 452 |
+ if(inherits(CA(object), "ff") | inherits(CA(object), "ffdf")) stop("Must specify i and/or j for ff objects")
|
|
| 453 |
+ } |
|
| 454 |
+ if(missing(i) & !missing(j)){
|
|
| 455 |
+ snp.index <- which(isSnp(object)) |
|
| 456 |
+ cn.total <- as.matrix(CA(cnSet)[, j]) |
|
| 457 |
+ cb <- as.matrix(CB(cnSet)[snp.index, j] ) |
|
| 458 |
+ cn.total[snp.index, ] <- cn.total[snp.index, ] + cb |
|
| 459 |
+ } |
|
| 460 |
+ if(!missing(i) & missing(j)){
|
|
| 461 |
+ snp.index <- intersect(which(isSnp(object)), i) |
|
| 462 |
+ cn.total <- as.matrix(CA(cnSet)[i, ]) |
|
| 463 |
+ cb <- as.matrix(CB(cnSet)[snp.index, ]) |
|
| 464 |
+ cn.total[snp.index, ] <- cn.total[snp.index, ] + cb |
|
| 465 |
+ } |
|
| 466 |
+ if(!missing(i) & !missing(j)){
|
|
| 467 |
+ snp.index <- intersect(which(isSnp(object)), i) |
|
| 468 |
+ cn.total <- as.matrix(CA(cnSet)[i, j]) |
|
| 469 |
+ cb <- as.matrix(CB(cnSet)[snp.index, j]) |
|
| 470 |
+ cn.total[snp.index, ] <- cn.total[snp.index, ] + cb |
|
| 471 |
+ } |
|
| 472 |
+ cn.total <- cn.total/100 |
|
| 473 |
+ dimnames(cn.total) <- NULL |
|
| 474 |
+ return(cn.total) |
|
| 475 |
+}) |
|
| 476 |
+ |
|
| 477 |
+##setMethod("copyNumber", "CNSet", function(object){
|
|
| 478 |
+## I <- isSnp(object) |
|
| 479 |
+## CA <- CA(object) |
|
| 480 |
+## CB <- CB(object) |
|
| 481 |
+## CN <- CA + CB |
|
| 482 |
+## ##For nonpolymorphic probes, CA is the total copy number |
|
| 483 |
+## CN[!I, ] <- CA(object)[!I, ] |
|
| 484 |
+## CN |
|
| 485 |
+##}) |
|
| 486 |
+>>>>>>> Exporting list and ffdf classes so that lM<- assignment works |
|
| 352 | 487 |
|
| 353 | 488 |
|
| 354 | 489 |
ACN <- function(object, allele, i , j){
|