comparison limma_voom.R @ 13:d5a940112511 draft

planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/limma_voom commit 8560b34a261fde200bd77dc2e817d55d386ac811
author iuc
date Sun, 30 Sep 2018 10:51:29 -0400
parents 81796eb60bd0
children 3133e833b3ce
comparison
equal deleted inserted replaced
12:81796eb60bd0 13:d5a940112511
369 stripOutPdf[i] <- makeOut(paste0("stripcharts_", con, ".pdf")) 369 stripOutPdf[i] <- makeOut(paste0("stripcharts_", con, ".pdf"))
370 mdvolOutPng[i] <- makeOut(paste0("mdvolplot_", con, ".png")) 370 mdvolOutPng[i] <- makeOut(paste0("mdvolplot_", con, ".png"))
371 topOut[i] <- makeOut(paste0(deMethod, "_", con, ".tsv")) 371 topOut[i] <- makeOut(paste0(deMethod, "_", con, ".tsv"))
372 glimmaOut[i] <- makeOut(paste0("glimma_", con, "/MD-Plot.html")) 372 glimmaOut[i] <- makeOut(paste0("glimma_", con, "/MD-Plot.html"))
373 } 373 }
374 filtOut <- makeOut(paste0(deMethod, "_filtcounts.tsv")) 374 filtOut <- makeOut(paste0(deMethod, "_", "filtcounts"))
375 normOut <- makeOut(paste0(deMethod, "_normcounts.tsv")) 375 normOut <- makeOut(paste0(deMethod, "_", "normcounts"))
376 rdaOut <- makeOut(paste0(deMethod, "_analysis.RData")) 376 rdaOut <- makeOut(paste0(deMethod, "_analysis.RData"))
377 sessionOut <- makeOut("session_info.txt") 377 sessionOut <- makeOut("session_info.txt")
378 378
379 # Initialise data for html links and images, data frame with columns Label and 379 # Initialise data for html links and images, data frame with columns Label and
380 # Link 380 # Link
449 449
450 if (wantFilt) { 450 if (wantFilt) {
451 print("Outputting filtered counts") 451 print("Outputting filtered counts")
452 filt_counts <- data.frame(data$genes, data$counts) 452 filt_counts <- data.frame(data$genes, data$counts)
453 write.table(filt_counts, file=filtOut, row.names=FALSE, sep="\t", quote=FALSE) 453 write.table(filt_counts, file=filtOut, row.names=FALSE, sep="\t", quote=FALSE)
454 linkData <- rbind(linkData, data.frame(Label=paste0(deMethod, "_", "filtcounts.tsv"), Link=paste0(deMethod, "_", "filtcounts.tsv"), stringsAsFactors=FALSE)) 454 linkData <- rbind(linkData, data.frame(Label=paste0(deMethod, "_", "filtcounts.tsv"), Link=paste0(deMethod, "_", "filtcounts"), stringsAsFactors=FALSE))
455 } 455 }
456 456
457 # Plot Density 457 # Plot Density
458 if ("d" %in% plots) { 458 if ("d" %in% plots) {
459 # PNG 459 # PNG
721 plotData <- logCPM 721 plotData <- logCPM
722 722
723 # Save normalised counts (log2cpm) 723 # Save normalised counts (log2cpm)
724 if (wantNorm) { 724 if (wantNorm) {
725 write.table(logCPM, file=normOut, row.names=TRUE, sep="\t", quote=FALSE) 725 write.table(logCPM, file=normOut, row.names=TRUE, sep="\t", quote=FALSE)
726 linkData <- rbind(linkData, c((paste0(deMethod, "_", "normcounts.tsv")), (paste0(deMethod, "_", "normcounts.tsv")))) 726 linkData <- rbind(linkData, c((paste0(deMethod, "_", "normcounts.tsv")), (paste0(deMethod, "_", "normcounts"))))
727 } 727 }
728 } else { 728 } else {
729 # limma-voom approach 729 # limma-voom approach
730 voomOutPdf <- makeOut("voomplot.pdf") 730 voomOutPdf <- makeOut("voomplot.pdf")
731 voomOutPng <- makeOut("voomplot.png") 731 voomOutPng <- makeOut("voomplot.png")
772 772
773 # Save normalised counts (log2cpm) 773 # Save normalised counts (log2cpm)
774 if (wantNorm) { 774 if (wantNorm) {
775 norm_counts <- data.frame(vData$genes, vData$E) 775 norm_counts <- data.frame(vData$genes, vData$E)
776 write.table(norm_counts, file=normOut, row.names=FALSE, sep="\t", quote=FALSE) 776 write.table(norm_counts, file=normOut, row.names=FALSE, sep="\t", quote=FALSE)
777 linkData <- rbind(linkData, c((paste0(deMethod, "_", "normcounts.tsv")), (paste0(deMethod, "_", "normcounts.tsv")))) 777 linkData <- rbind(linkData, c((paste0(deMethod, "_", "normcounts.tsv")), (paste0(deMethod, "_", "normcounts"))))
778 } 778 }
779 779
780 # Fit linear model and estimate dispersion with eBayes 780 # Fit linear model and estimate dispersion with eBayes
781 voomFit <- contrasts.fit(voomFit, contrasts) 781 voomFit <- contrasts.fit(voomFit, contrasts)
782 if (wantRobust) { 782 if (wantRobust) {
1050 } 1050 }
1051 } 1051 }
1052 1052
1053 cata("<h4>Tables:</h4>\n") 1053 cata("<h4>Tables:</h4>\n")
1054 for (i in 1:nrow(linkData)) { 1054 for (i in 1:nrow(linkData)) {
1055 if (grepl(".tsv", linkData$Link[i])) { 1055 if (grepl("counts$", linkData$Link[i])) {
1056 HtmlLink(linkData$Link[i], linkData$Label[i])
1057 } else if (grepl(".tsv", linkData$Link[i])) {
1056 HtmlLink(linkData$Link[i], linkData$Label[i]) 1058 HtmlLink(linkData$Link[i], linkData$Label[i])
1057 } 1059 }
1058 } 1060 }
1059 1061
1060 if (wantRda) { 1062 if (wantRda) {