Mercurial > repos > pieterlukasse > prims_metabolomics

diff metaMS_cmd_pick_and_group.r @ 49:f772a5caa86a

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
Added more options and better documentation. Added MsClust support for parsing XCMS alignment results. Improved output reports for XCMS wrappers. New tools.
author pieter.lukasse@wur.nl
date Wed, 10 Dec 2014 22:03:27 +0100
parents
children 684d2341968c
line wrap: on
line diff
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/metaMS_cmd_pick_and_group.r	Wed Dec 10 22:03:27 2014 +0100
@@ -0,0 +1,89 @@
+## read args:
+args <- commandArgs(TRUE)
+## data files, e.g. "E:/Rworkspace/metaMS/data/data.zip" (with e.g. .CDF files) and unzip output dir, e.g. "E:/"
+args.dataZip <- args[1]
+args.zipExtrDir <- sub("\\.","_",paste(args[1],"dir", sep=""))
+dir.create(file.path(args.zipExtrDir), showWarnings = FALSE, recursive = TRUE)
+## settings file, e.g. "E:/Rworkspace/metaMS/data/settings.r", should contain assignment to an object named "customMetaMSsettings" 
+args.settings <- args[2]
+
+## output file names, e.g. "E:/Rworkspace/metaMS/data/out.txt"
+args.outPeakTable <- args[3]
+args.xsetOut <- args[4]
+
+## report files
+args.htmlReportFile <- args[5]
+args.htmlReportFile.files_path <- args[6]
+
+
+if (length(args) == 7)
+{
+	args.outLogFile <- args[7]
+	# suppress messages:
+	# Send all STDERR to STDOUT using sink() see http://mazamascience.com/WorkingWithData/?p=888
+	msg <- file(args.outLogFile, open="wt")
+	sink(msg, type="message") 
+	sink(msg, type="output")
+}
+
+cat("\nSettings used===============:\n")
+cat(readChar(args.settings, 1e5))
+
+
+tryCatch(
+        {
+	        library(metaMS)
+	
+			## load the data files from a zip file
+			files <- unzip(args.dataZip, exdir=args.zipExtrDir)
+			
+			## load settings "script" into "customMetaMSsettings"
+			tempEnv <- new.env() 
+			source(args.settings, local=tempEnv)
+			message(paste(" loaded : ", args.settings))
+			allSettings <- tempEnv[["customMetaMSsettings"]] 
+			
+			# trigger runLC: 
+			LC <- runLC(files, settings = allSettings, nSlaves=20, returnXset = TRUE)
+			
+			# write out runLC annotation results:
+			write.table(LC$PeakTable, args.outPeakTable, sep="\t", row.names=FALSE)
+			
+			# save xset as rdata:
+			xsAnnotatePreparedData <- LC$xset
+			saveRDS(xsAnnotatePreparedData, file=args.xsetOut)
+			
+			message("\nGenerating report.........")
+			# report
+			dir.create(file.path(args.htmlReportFile.files_path), showWarnings = FALSE, recursive = TRUE)
+			html <- "<html><body><h1>Info on alignment quality </h1>" 
+			# TODO add (nr and mass error) and group size
+			
+			message("\nPlotting figures... ")
+			figureName <- paste(args.htmlReportFile.files_path, "/figure_retcor.png", sep="")
+			html <- paste(html,"<img src='figure_retcor.png' /><br/>", sep="") 
+			png( figureName, type="cairo", width=1100,height=600 ) 
+			retcor(LC$xset@xcmsSet, method="peakgroups", plottype = "mdevden")
+			html <- paste(html,"<a>*NB: retention time correction plot based on 'peakgroups' option with default settings. This is not the plot matching the exact settings used in the run, 
+									but just intended to give a rough estimate of the retention time shifts present in the data. A more accurate plot will be available once
+                                    this option is added in metaMS API. </a><br/>", sep="")
+			devname = dev.off()
+			
+			
+			gt <- groups(LC$xset@xcmsSet)
+			groupidx1 <- which(gt[,"rtmed"] > 0 & gt[,"rtmed"] < 3000 & gt[,"npeaks"] > 3)
+			
+			html <- paste(html,"</body><html>")
+			message("finished generating report")
+			write(html,file=args.htmlReportFile)
+			# unlink(args.htmlReportFile)
+			cat("\nWarnings================:\n")
+			str( warnings() ) 
+		},
+        error=function(cond) {
+            sink(NULL, type="message") # default setting
+			sink(stderr(), type="output")
+            message("\nERROR: ===========\n")
+            print(cond)
+        }
+    )