DiagnosticImages.R

From Organic Design wiki
Revision as of 22:17, 8 January 2006 by Sven (talk | contribs) (document)
(diff) ← Older revision | Latest revision (diff) | Newer revision → (diff)

library(limma)

source("readData.R")

  1. 1) Create diagnostics directories

dir.exists <- function(x, curdir = ".", rootdir) {

 if( missing(rootdir)) {
   rootdir <- getwd()
 }
 setwd(curdir)
 
 x <- unlist(strsplit(x,"/"))
 if( !file.exists(x[1]) ){
   dir.create(x[1])
   print(paste("Creating", file.path(getwd(), x[1])))
 }
 curdir <- x[1]
 if(length(x)==1) {
   setwd(rootdir)
 } else {
   return(dir.exists(x[-1], curdir = curdir,  rootdir = rootdir))
 }

}

bc.correct <- c("none","subtract","minimum", "movingmin") method.within <- c("none","loess") method.between <- c("none","scale", "quantile")

  1. 3) Making Diagnostic directories

plotDir <- "tmp" dir.exists(file.path(plotDir, "Densityplots")) for(i in c("ImagePlots","MAplots")) {

 dir.exists(file.path(plotDir, i))
 if(i == "ImagePlots") {
   for(x in bc.correct) {
     for(y in method.within) {
       for(z in method.between) {
         for( j in c("R","G", "M","A")) {
           dir.exists(file.path(plotDir, i, paste(x,"-",y,"-",z, sep=""),j))
           if(x=="none" & y =="none" & z =="none") {
             for(k in c("Rb", "Gb")) {
               dir.exists(file.path(plotDir, i, paste(x,"-",y,"-",z, sep=""),k))
             }
           }
         }
       }
     }
   }
 }
 if(i == "MAplots") {
   for(x in bc.correct) {
     for(y in method.within) {
       for(z in method.between) {
         dir.exists(file.path(plotDir, i, paste(x,"-",y,"-",z, sep="")))
       }
     }
   }
 }

}


  1. 2) Diagnostic plots for raw and loess normalized data (with background subtraction)
  1. Density plots

for(i in bc.correct) {

 for( j in method.within) {
   MA <- normalizeWithinArrays(RG, method=j, bc.method=i)
   for (k in method.between) {
     MA.Between <- normalizeBetweenArrays(MA, method=k)
     print(paste("[ Graphing:",i,j,k,"]", sep=" "))
     png(file.path(plotDir, "Densityplots", paste(i,"-",j,"-",k,".png", sep="")))
     plotDensities(MA.Between)
     dev.off()
   }
 }

}

  1. MAplots

for(i in bc.correct) {

 for( j in method.within) {
   MA <- normalizeWithinArrays(RG, method=j, bc.method=i)
   for (k in method.between) {
     MA.Between <- normalizeBetweenArrays(MA, method=k)
     print(paste("[ Graphing:",i,j,k,"]", sep=" "))
     plotMA3by2(MA.Between, path=file.path(plotDir, "MAplots", paste(i,"-",j,"-",k, sep="")))
   }
 }

}

  1. Imageplots

for(i in bc.correct) {

 for( j in method.within) { # handles background correction and within normalization without overwriting RG
   MA <- normalizeWithinArrays(RG, method=j, bc.method=i) 
   for (k in method.between) {
     MA.Between <- normalizeBetweenArrays(MA, method=k)
     print(paste("[ Graphing:",i,j,k,"]", sep=" "))
     imageplot3by2(MA.Between, z="M", path = file.path(plotDir,"Imageplots",  paste(i,"-",j,"-",k, sep=""), "M"))
     imageplot3by2(MA.Between, z="A", path = file.path(plotDir,"Imageplots",  paste(i,"-",j,"-",k, sep=""), "A"))
     imageplot3by2(RG.MA(MA.Between), z="G", path = file.path(plotDir,"Imageplots", paste(i,"-",j,"-",k, sep=""),"G"),
                   low="black", high="green")
     imageplot3by2(RG.MA(MA.Between), z="R", path = file.path(plotDir,"Imageplots", paste(i,"-",j,"-",k, sep=""),"R"),
                   low="black", high="red")
     if(i=="none" & j=="none" & k=="none") {
       imageplot3by2(RG, z="Gb", path = file.path(plotDir,"Imageplots", paste(i,"-",j,"-",k, sep=""),"Gb"),
                     low="black", high="green")
       imageplot3by2(RG, z="Rb", path = file.path(plotDir,"Imageplots", paste(i,"-",j,"-",k, sep=""),"Rb"),
                     low="black", high="red")
     }
   }
 }

}


  1. Splom pairs plots

library(lattice) tmp <- as.data.frame(rbind(RG$R, RG$G)) tmp$channel <- rep(c("log2(Red)","log2(Green)"), each=nrow(tmp)/2)

dir.exists(file.path(plotDir,"PairsPlots"))

for( i in seq(1, ncol(RG), by=3)) {

 png(file.path(plotDir, "PairsPlots", paste(targets[,3][i],".png",sep="")))
 print(splom(~ log2(tmp[,i:(i+2)])| tmp$channel, xlab=targets[,3][i]))
 dev.off()

}

  1. Comparing Pool 3 Cy5 to Pool 11 Cy5 (same treatment)

P3ind <- grep("P3", targets[,"Filetype"]) P9ind <- grep("P9", targets[,"Filetype"])

tmp <- RG tmp$R <- RG$R[,P3ind] tmp$Rb <- RG$Rb[,P3ind] tmp$G <- RG$R[,P9ind] tmp$Gb <- RG$Rb[,P9ind] tmp$weights <- NULL tmp

  1. Shows when scanning is clipping

par(mfrow=c(2,2)) plotMA(tmp, array=1) plotMA(tmp, array=2) plotMA(tmp, array=3)