diff script.R @ 0:bd8fd161908b draft default tip

planemo upload for repository https://github.com/TGAC/earlham-galaxytools/tree/master/tools/plotheatmap commit 7e9d6b283879b45d5fe673767668c06104d31609-dirty
author earlhaminst
date Wed, 12 Jul 2017 14:29:51 -0400
parents
children
line wrap: on
line diff
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/script.R	Wed Jul 12 14:29:51 2017 -0400
@@ -0,0 +1,91 @@
+#!/usr/bin/env Rscript
+
+args = commandArgs(TRUE)
+countDataPath = args[1]
+statsDataPath = args[2]
+logFC = args[3]
+logCPM = args[4]
+pValue = args[5]
+fdr = args[6]
+
+clusterRow = args[7]
+clusterCol = args[8]
+hclustMethod = args[9]
+
+mgColumnNm = as.numeric(args[10])
+mgRowNm = as.numeric(args[11])
+
+pdfWidth = as.numeric(args[12])
+pdfHeight = as.numeric(args[13])
+
+
+if(clusterRow == "Yes"){
+  clusterRow = TRUE
+} else {
+  clusterRow = NA
+}
+
+if(clusterCol == "Yes"){
+  clusterCol = TRUE
+} else {
+  clusterCol = NA
+}
+
+require(preprocessCore)
+require(gplots)
+
+#prepare counts data --------------------------------------------------------
+countData = read.table(countDataPath, comment = "",
+                       sep = "\t")
+
+groups = sapply(as.character(countData[1, -1]), strsplit, ":")
+groups = as.vector(t(countData[1, -1]))
+
+names = as.vector(t(countData[2, -1]))
+
+countData = countData[-c(1,2), ]
+rownames(countData) = countData[, 1]
+countData = countData[, -1]
+colnames(countData) = names
+
+countData = countData[, order(groups, names)]
+
+# prepare stats data ------------------------------------------------------
+
+statsData = read.table(statsDataPath, sep = "\t",
+                       header = T)
+
+colnames(statsData)[-1] = sapply(colnames(statsData)[-1], function(x){
+  unlist(strsplit(x, ".", fixed = T))[3]
+})
+
+wh = which(abs(statsData$logFC) >= logFC & statsData$logCPM >= logCPM & statsData$PValue <= pValue & statsData$FDR <= fdr)
+
+for(i in 1:ncol(countData)){
+  countData[,i] = as.numeric(countData[,i])
+}
+
+countDataNorm = normalize.quantiles(as.matrix(countData), copy = T)
+countDataNormLog = log(countDataNorm + 1, 2)
+
+colnames(countDataNormLog) = colnames(countData)
+rownames(countDataNormLog) = rownames(countData)
+
+#svg("heatmap.svg", width = 3+length(names), height = 1/2*length(wh))
+pdf("heatmap.pdf", width = pdfWidth, height = pdfHeight)
+
+heatmap.2(
+  countDataNormLog[wh, ],
+  density.info=c("none"),
+  hclustfun = function(x) hclust(x, method = hclustMethod),
+  distfun = function(x) as.dist(1-cor(t(x))),  
+  col = bluered(50), 
+  scale = "row", 
+  trace = "none", 
+  Rowv = clusterRow, 
+  Colv = clusterCol,
+  margins = c(mgColumnNm, mgRowNm)
+)
+
+dev.off()
+