SegmentPlots.R

Code snipits and programs written in R, S or S-PLUS library(limma) library(NutriTools)

type <- list(GenePix = "GPR") platform <- "GenePix" typeplatform

sessionInfo() auxillaryDir <- "/Volumes/HD2/Data/Nutrigenomics/MultipleScans/Auxillary" dataDir <- file.path("/Volumes/HD2/Data/Nutrigenomics/MultipleScans/", typeplatform) targets <- readTargets(file.path(auxillaryDir, "targets.txt"))

1. 1) ------------------------ No background correction ------------------------ #

RG <- readGPR(dataDir=dataDir) colnames(RG) <- paste(targets$Filetype, targets$Scan, sep="") MA <- MA.RG(RG, bc.method="none")

slides <- colnames(RG)

plotDir <- "DiagnosticImages/ArrowPlots" dir.create(plotDir, recursive=TRUE)

slides <- colnames(RG) for(i in seq(1,32, by=3)){

 for(j in 0:1) {
   print(paste("I:",i, " J:",j, sep=""))
   scanType <- c("Low","Medium","High")

   # Calculate MA scale slope coefficients
   coefs <- rep(NA, length=nrow(MA))
   for(k in seq(nrow(MA))) {
     x <- cbind(1, MA$A[k, (i+j):(i+j+1)])
     y <- MA$M[k, (i+j):(i+j+1)]
     if(!is.na(c(x[,2]))  %*% is.na(y)) {
       coefs[k] <- lm.fit(x,y)$coef[2]
     }
   }
   # Create colour lookup table
   cols <- rev(topo.colors(length(coefs)))
   
   pdf(file.path(plotDir, paste(platform, "-RGplot-", targets[i, "Filetype"], "-", scanType[j+2],"vs",scanType[j+1], ".pdf", sep="")), version = "1.4")

   plot(5:16, 5:16, type="n", xlab=expression(log[2](G)), ylab=expression(log[2](R)))
   counter <- 1
   for( k in rev(order(abs(coefs)))) {
     rgbCol <- col2rgb(cols[counter])
     counter <- counter+1
     arrows(log2(RG$G[k,slides[i+j]]), log2(RG$R[k, slides[i+j]]),  log2(RG$G[k,slides[i+j+1]]), log2(RG$R[k, slides[i+j+1]]), col=rgb(rgbCol[1], rgbCol[2], rgbCol[3], 0.5*255, maxColorValue=255), length=0.05)
   }
    # Overlay poor quality spots based on weight from either scan
   for(k in seq(nrow(RG))) {
     if(RG$weights[k, slides[i]]==0 | RG$weights[k, slides[i+j]]==0) {
       arrows(log2(RG$G[k,slides[i+j]]), log2(RG$R[k, slides[i+j]]),  log2(RG$G[k,slides[i+j+1]]), log2(RG$R[k, slides[i+j+1]]), col=rgb(0,0,0, 0.5), length=0.05)
     }
   }

   dev.off()

   pdf(file.path(plotDir, paste(platform,"-MAplot-", targets[i, "Filetype"], "-", scanType[j+2],"vs",scanType[j+1], ".pdf", sep="")), version = "1.4")

   plot(5:16, -5:6, type="n", xlab=expression(A==over(log[2](R)-log[2](G),2)), ylab=expression(M==log[2](R)-log[2](G)))
   counter <- 1
   for( k in rev(order(abs(coefs)))) {
     rgbCol <- col2rgb(cols[counter])
     counter <- counter+1
     arrows((MA$A[k,slides[i+j]]),(MA$M[k, slides[i+j]]),  (MA$A[k,slides[i+j+1]]), (MA$M[k, slides[i+j+1]]), col=rgb(rgbCol[1], rgbCol[2], rgbCol[3], 0.5*255, maxColorValue=255), length=0.05)
   }
   # Overlay poor quality spots based on weight from either scan
   for(k in seq(nrow(MA))) {
     if(MA$weights[k, slides[i]]==0 | MA$weights[k, slides[i+j]]==0) {
       arrows((MA$A[k,slides[i+j]]),(MA$M[k, slides[i+j]]),  (MA$A[k,slides[i+j+1]]), (MA$M[k, slides[i+j+1]]), col=rgb(0,0,0, 0.2), length=0.05)
     }
   }
   
   dev.off()
 }

}

1. 2) ------------------------- Background correction -------------------------- #

RGbc <- backgroundCorrect(RG, method="subtract") MAbc <- MA.RG(RGbc)

slides <- colnames(RGbc) for(i in seq(1,32, by=3)){

 for(j in 0:1) {
   print(paste("I:",i, " J:",j, sep=""))
   scanType <- c("Low","Medium","High")

       # Calculate MA scale slope coefficients
   coefs <- rep(NA, length=nrow(MA))
   for(k in seq(nrow(MA))) {
     x <- cbind(1, MA$A[k, (i+j):(i+j+1)])
     y <- MA$M[k, (i+j):(i+j+1)]
     if(!is.na(c(x[,2]))  %*% is.na(y)) {
       coefs[k] <- lm.fit(x,y)$coef[2]
     }
   }
   # Create colour lookup table
   cols <- rev(topo.colors(length(coefs)))
   
   pdf(file.path(plotDir, paste(platform,"-RGplot-bc-", targets[i, "Filetype"], "-", scanType[j+2],"vs",scanType[j+1], ".pdf", sep="")), version = "1.4")

   plot(c(0,16), c(0,16), type="n", xlab=expression(log[2](G)), ylab=expression(log[2](R)))
   counter <- 1
   for(k in seq(nrow(MAbc))) {
     arrows(log2(RGbc$G[k,slides[i+j]]), log2(RGbc$R[k, slides[i+j]]),  log2(RGbc$G[k,slides[i+j+1]]), log2(RGbc$R[k, slides[i+j+1]]), col=rgb(rgbCol[1], rgbCol[2], rgbCol[3], 0.5*255, maxColorValue=255), length=0.05)
   }
   # Overlay poor quality spots based on weight from either scan
   for(k in seq(nrow(MA))) {
     if(RG$weights[k, slides[i]]==0 | RG$weights[k, slides[i+j]]==0) {
       arrows(log2(RG$G[k,slides[i+j]]), log2(RG$R[k, slides[i+j]]),  log2(RG$G[k,slides[i+j+1]]), log2(RG$R[k, slides[i+j+1]]), col=rgb(0,0,0,0.5), length=0.05)
     }
   }

   dev.off()

   pdf(file.path(plotDir, paste(platform,"-MAplot-bc-", targets[i, "Filetype"], "-", scanType[j+2],"vs",scanType[j+1], ".pdf", sep="")), version = "1.4")

   plot(c(0,16), c(-8,8), type="n", xlab=expression(A==over(log[2](R)-log[2](G),2)), ylab=expression(M==log[2](R)-log[2](G)))
   coefs <- rep(0, length=nrow(MAbc)) # NA's will not evaluate a slope so plot these first by assigning 0
   for(k in seq(nrow(MAbc))) {
     x <- cbind(1, MAbc$A[k, (i+j):(i+j+1)])
     y <- MAbc$M[k, (i+j):(i+j+1)]
     if(!is.na(c(x[,2]))  %*% is.na(y)) {
       coefs[k] <- lm.fit(x,y)$coef[2]
     }
   }
   counter <- 1
   for( k in rev(order(abs(coefs)))) {  # Reverse or forward ordering is informative here
     rgbCol <- col2rgb(cols[counter])
     counter <- counter+1
     arrows((MAbc$A[k,slides[i+j]]),(MAbc$M[k, slides[i+j]]),  (MAbc$A[k,slides[i+j+1]]), (MAbc$M[k, slides[i+j+1]]), col=rgb(rgbCol[1], rgbCol[2], rgbCol[3], maxColorValue=255), length=0.05) #  0.5*255,
   }
   # Overlay poor quality spots based on weight from either scan
   for(k in seq(nrow(MA))) {
     if(MA$weights[k, slides[i]]==0 | MA$weights[k, slides[i+j]]==0) {
       arrows((MA$A[k,slides[i+j]]),(MA$M[k, slides[i+j]]),  (MA$A[k,slides[i+j+1]]), (MA$M[k, slides[i+j+1]]), col=rgb(0,0,0,0.5), length=0.05)
     }
   }
   dev.off()
 }

}