inst/doc/ADaCGH2.R

In ADaCGH2: Analysis of big data from aCGH experiments using parallel computing and ff objects

### R code from vignette source 'ADaCGH2.Rnw'

###################################################
### code chunk number 1: ADaCGH2.Rnw:329-334
###################################################
library(ADaCGH2)

data(inputEx)
summary(inputEx)
head(inputEx)


###################################################
### code chunk number 2: ADaCGH2.Rnw:425-430
###################################################
fnameRdata <- list.files(path = system.file("data", package = "ADaCGH2"),
                     full.names = TRUE, pattern = "inputEx.RData")

inputToADaCGH(ff.or.RAM = "RAM",
                      RDatafilename = fnameRdata)


###################################################
### code chunk number 3: ADaCGH2.Rnw:458-461
###################################################
data(inputEx) ## make inputEx available as a data frame with that name
inputToADaCGH(ff.or.RAM = "RAM",
                      dataframe = inputEx)


###################################################
### code chunk number 4: ADaCGH2.Rnw:472-476
###################################################
data(inputEx)
cgh.dat <- inputEx[, -c(1, 2, 3)]
chrom.dat <- as.integer(inputEx[, 2])
pos.dat <- inputEx[, 3]


###################################################
### code chunk number 5: ADaCGH2.Rnw:492-499
###################################################
fnametxt <- list.files(path = system.file("data", package = "ADaCGH2"),
                         full.names = TRUE, pattern = "inputEx.txt")

##    You might want to adapt mc.cores to your hardware
tmp <- inputToADaCGH(ff.or.RAM = "RAM",
                     textfilename = fnametxt,
                     mc.cores = 2)


###################################################
### code chunk number 6: ADaCGH2.Rnw:560-561 (eval = FALSE)
###################################################
## help(pSegment)


###################################################
### code chunk number 7: ADaCGH2.Rnw:564-569
###################################################
##    You might want to adapt mc.cores to your hardware
haar.RAM.fork <- pSegmentHaarSeg(cgh.dat, chrom.dat,
                                   merging = "MAD",
                                 typeParall = "fork",
                                 mc.cores = 2)


###################################################
### code chunk number 8: ADaCGH2.Rnw:575-577
###################################################
lapply(haar.RAM.fork, head)
summary(haar.RAM.fork[[1]])


###################################################
### code chunk number 9: ADaCGH2.Rnw:592-601
###################################################
##    You might want to adapt mc.cores to your hardware
pChromPlot(haar.RAM.fork,
           cghRDataName = cgh.dat,
           chromRDataName = chrom.dat,
           posRDataName = pos.dat,
           probenamesRDataName = probenames.dat,
           imgheight = 350,
           typeParall = "fork",
           mc.cores = 2)


###################################################
### code chunk number 10: ADaCGH2.Rnw:659-665
###################################################
if(.Platform$OS.type != "windows") {

originalDir <- getwd()
## make it explicit where we are
print(originalDir)
}


###################################################
### code chunk number 11: ADaCGH2.Rnw:668-673
###################################################
if(.Platform$OS.type != "windows") {
if(!file.exists("ADaCGH2_vignette_tmp_dir"))
  dir.create("ADaCGH2_vignette_tmp_dir")
setwd("ADaCGH2_vignette_tmp_dir")
}


###################################################
### code chunk number 12: ADaCGH2.Rnw:715-721
###################################################
if(.Platform$OS.type != "windows") {
fnameRdata <- list.files(path = system.file("data", package = "ADaCGH2"),
                     full.names = TRUE, pattern = "inputEx.RData")
inputToADaCGH(ff.or.RAM = "ff",
                      RDatafilename = fnameRdata)
}


###################################################
### code chunk number 13: ADaCGH2.Rnw:771-778 (eval = FALSE)
###################################################
## mcparallel(inputToADaCGH(ff.or.RAM = "ff",
##                                  RDatafilename = fnameRData),
##                                  silent = FALSE)
##   tableChromArray <- mccollect()
##   if(inherits(tableChromArray, "try-error")) {
##     stop("ERROR in input data conversion")
##   }


###################################################
### code chunk number 14: ADaCGH2.Rnw:792-797
###################################################
if(.Platform$OS.type != "windows") {
data(inputEx) ## make inputEx available as a data frame with that name
inputToADaCGH(ff.or.RAM = "ff",
                      dataframe = inputEx)
}


###################################################
### code chunk number 15: ADaCGH2.Rnw:813-822
###################################################
if(.Platform$OS.type != "windows") {
fnametxt <- list.files(path = system.file("data", package = "ADaCGH2"),
                         full.names = TRUE, pattern = "inputEx.txt")

##    You might want to adapt mc.cores to your hardware
inputToADaCGH(ff.or.RAM = "ff",
              textfilename = fnametxt,
              mc.cores = 2)
}


###################################################
### code chunk number 16: ADaCGH2.Rnw:873-885
###################################################
if(.Platform$OS.type != "windows") {
## Adapt number of nodes to your hardware    
number.of.nodes <- 2 ##detectCores()
cl2 <- parallel::makeCluster(number.of.nodes, "PSOCK")
parallel::clusterSetRNGStream(cl2)
parallel::setDefaultCluster(cl2) 
parallel::clusterEvalQ(NULL, library("ADaCGH2"))

wdir <- getwd()
parallel::clusterExport(NULL, "wdir")
parallel::clusterEvalQ(NULL, setwd(wdir))
}


###################################################
### code chunk number 17: ADaCGH2.Rnw:915-916 (eval = FALSE)
###################################################
## help(pSegment)


###################################################
### code chunk number 18: ADaCGH2.Rnw:918-924
###################################################
if(.Platform$OS.type != "windows") {
    haar.ff.cluster <- pSegmentHaarSeg("cghData.RData",
                                   "chromData.RData", 
                                   merging = "MAD",
                                   typeParall = "cluster")
}


###################################################
### code chunk number 19: ADaCGH2.Rnw:932-936
###################################################
if(.Platform$OS.type != "windows") {
lapply(haar.ff.cluster, open)
summary(haar.ff.cluster[[1]][,])
}


###################################################
### code chunk number 20: ADaCGH2.Rnw:949-960
###################################################
if(.Platform$OS.type != "windows") {
save(haar.ff.cluster, file = "hs_mad.out.RData", compress = FALSE)

pChromPlot(outRDataName = "hs_mad.out.RData",
           cghRDataName = "cghData.RData",
           chromRDataName = "chromData.RData",
           posRDataName = "posData.RData",
           probenamesRDataName = "probeNames.RData",
           imgheight = 350,
           typeParall = "cluster")
}


###################################################
### code chunk number 21: ADaCGH2.Rnw:965-968
###################################################
if(.Platform$OS.type != "windows") {
parallel::stopCluster(cl2)
}


###################################################
### code chunk number 22: ADaCGH2.Rnw:1028-1034
###################################################
if(.Platform$OS.type != "windows") {
fnameRdata <- list.files(path = system.file("data", package = "ADaCGH2"),
                     full.names = TRUE, pattern = "inputEx.RData")
inputToADaCGH(ff.or.RAM = "ff",
                      RDatafilename = fnameRdata)
}


###################################################
### code chunk number 23: ADaCGH2.Rnw:1048-1055 (eval = FALSE)
###################################################
## mcparallel(inputToADaCGH(ff.or.RAM = "ff",
##                                  RDatafilename = fnameRdata), 
##            silent = FALSE)
## tableChromArray <- collect()
## if(inherits(tableChromArray, "try-error")) {
##   stop("ERROR in input data conversion")
## }


###################################################
### code chunk number 24: ADaCGH2.Rnw:1068-1077
###################################################
if(.Platform$OS.type != "windows") {
fnametxt <- list.files(path = system.file("data", package = "ADaCGH2"),
                         full.names = TRUE, pattern = "inputEx.txt")

##    You might want to adapt mc.cores to your hardware
inputToADaCGH(ff.or.RAM = "ff",
              textfilename = fnametxt,
              mc.cores = 2)
}


###################################################
### code chunk number 25: ADaCGH2.Rnw:1138-1156
###################################################
if( (.Platform$OS.type == "unix") && (Sys.info()['sysname'] != "Darwin") ) {
  fnametxt <- list.files(path = system.file("data", package = "ADaCGH2"),
                          full.names = TRUE, pattern = "inputEx.txt")
  if(file.exists("cuttedFile")) {
    stop("The cuttedFile directory already exists. ",
         "Did you run this vignette from this directory before? ",
         "You will not want to do that, unless you modify the arguments ",
         "to inputToADaCGH below")
  } else dir.create("cuttedFile")
  setwd("cuttedFile")
  ## You might want to adapt mc.cores to your hardware
  cutFile(fnametxt, 1, 2, 3, sep = "\t", mc.cores = 2)
  cuttedFile.dir <- getwd()
  setwd("../")
} else {
  cuttedFile.dir <- system.file("example-datadir", 
                                 package = "ADaCGH2")
}


###################################################
### code chunk number 26: chunk0001
###################################################
if(.Platform$OS.type != "windows") {
## You might want to adapt mc.cores to your hardware    
inputToADaCGH(ff.or.RAM = "ff",
              path = cuttedFile.dir,
              verbose = TRUE,
              mc.cores = 2)
}


###################################################
### code chunk number 27: chunk0002
###################################################
if(.Platform$OS.type != "windows") {
## You might want to adapt mc.cores to your hardware
haar.ff.fork <- pSegmentHaarSeg("cghData.RData",
                                "chromData.RData",
                                 merging = "MAD",
                                typeParall = "fork",
                                mc.core = 2)
}


###################################################
### code chunk number 28: chunk0003
###################################################
if(.Platform$OS.type != "windows") {
save(haar.ff.fork, file = "haar.ff.fork.RData", compress = FALSE)

##    You might want to adapt mc.cores to your hardware
pChromPlot(outRDataName = "haar.ff.fork.RData",
           cghRDataName = "cghData.RData",
           chromRDataName = "chromData.RData",
           posRDataName = "posData.RData",
           probenamesRDataName = "probeNames.RData",
           imgheight = 350,
           typeParall = "fork",
           mc.cores = 2)
}


###################################################
### code chunk number 29: chunk0004
###################################################
if(.Platform$OS.type != "windows") {
require("limma")
require("snapCGH")
datadir <- system.file("testdata", package = "snapCGH")
targets <- readTargets("targets.txt", path = datadir)
RG1 <- read.maimages(targets$FileName, path = datadir, source = "genepix")


## This is snapCGH-specific
RG1 <- read.clonesinfo("cloneinfo.txt", RG1, path = datadir) 
RG1$printer <- getLayout(RG1$genes)
types <- readSpotTypes("SpotTypes.txt", path = datadir)
RG1$genes$Status <- controlStatus(types, RG1)

RG1$design <- c(-1, -1)


RG2 <- backgroundCorrect(RG1, method = "minimum") ## class RGList
MA <- normalizeWithinArrays(RG2, method = "median") ## class MAList
class(MA)
## now obtain an object of class SegList
MA2 <- processCGH(MA, method.of.averaging = mean, ID = "ID") 
class(MA2)
}


###################################################
### code chunk number 30: chunk0005
###################################################
if(.Platform$OS.type != "windows") {
tmp <- inputToADaCGH(MAList = MA, 
                             robjnames = c("cgh-ma.dat", "chrom-ma.dat",
                                          "pos-ma.dat", "probenames-ma.dat"))

tmp <- inputToADaCGH(MAList = MA2, 
                             robjnames = c("cgh-ma.dat", "chrom-ma.dat",
                                          "pos-ma.dat", "probenames-ma.dat"),
                             minNumPerChrom = 4)
}


###################################################
### code chunk number 31: chunk0006
###################################################
if(.Platform$OS.type != "windows") {
targets.limma <- readTargets("targets.txt", path = datadir)
RG.limma <- read.maimages(targets.limma, path = datadir, 
                          source="genepix")
RG.limma <- backgroundCorrect(RG.limma, method="normexp", 
                              offset=50)
MA.limma <- normalizeWithinArrays(RG.limma)
}


###################################################
### code chunk number 32: chunk0007
###################################################
if(.Platform$OS.type != "windows") {
fclone <- list.files(path = system.file("testdata", package = "snapCGH"),
                     full.names = TRUE, pattern = "cloneinfo.txt")
fclone
tmp <- inputToADaCGH(MAList = MA.limma, 
                             cloneinfo = fclone,
                             robjnames = c("cgh-ma.dat", "chrom-ma.dat",
                                          "pos-ma.dat", "probenames-ma.dat"))
}


###################################################
### code chunk number 33: chunk0008
###################################################
if(.Platform$OS.type != "windows") {
acloneinfo <- MA$genes
tmp <- inputToADaCGH(MAList = MA.limma, 
                             cloneinfo = acloneinfo,
                             robjnames = c("cgh-ma.dat", "chrom-ma.dat",
                                          "pos-ma.dat", "probenames-ma.dat"))
}


###################################################
### code chunk number 34: chunk0009
###################################################
if(.Platform$OS.type != "windows") {
forcghr <- outputToCGHregions(haar.ff.cluster)
if(require(CGHregions)) {
  regions1 <- CGHregions(na.omit(forcghr))
  regions1
}
}


###################################################
### code chunk number 35: chunk00010
###################################################
if(.Platform$OS.type != "windows") {
## We are done with the executable code in the vignette.
## Restore the directory
setwd(originalDir)
print(getwd())
}


###################################################
### code chunk number 36: chunk00011
###################################################
if(.Platform$OS.type != "windows") {
## Remove the tmp dir. Sys.sleep to prevent Windoze problems.
## Sys.sleep(1)
## What is in that dir?
dir("ADaCGH2_vignette_tmp_dir")
unlink("ADaCGH2_vignette_tmp_dir", recursive = TRUE)
## Sys.sleep(1)
}


###################################################
### code chunk number 37: ADaCGH2.Rnw:2263-2264
###################################################
sessionInfo()


###################################################
### code chunk number 38: ADaCGH2.Rnw:2271-2274
###################################################
## try(closeAllConnections()) ## nope, breaks things
## Remove this sys.sleep; is this sometimes breaking windoze build in BioC?
## Sys.sleep(1) ## I hope we do not need this