Mercurial > repos > pieterlukasse > prims_metabolomics

comparison 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
comparison
equal deleted inserted replaced
48:26b93438f30e 49:f772a5caa86a
1 ## read args:
2 args <- commandArgs(TRUE)
3 ## data files, e.g. "E:/Rworkspace/metaMS/data/data.zip" (with e.g. .CDF files) and unzip output dir, e.g. "E:/"
4 args.dataZip <- args[1]
5 args.zipExtrDir <- sub("\\.","_",paste(args[1],"dir", sep=""))
6 dir.create(file.path(args.zipExtrDir), showWarnings = FALSE, recursive = TRUE)
7 ## settings file, e.g. "E:/Rworkspace/metaMS/data/settings.r", should contain assignment to an object named "customMetaMSsettings"
8 args.settings <- args[2]
9
10 ## output file names, e.g. "E:/Rworkspace/metaMS/data/out.txt"
11 args.outPeakTable <- args[3]
12 args.xsetOut <- args[4]
13
14 ## report files
15 args.htmlReportFile <- args[5]
16 args.htmlReportFile.files_path <- args[6]
17
18
19 if (length(args) == 7)
20 {
21 args.outLogFile <- args[7]
22 # suppress messages:
23 # Send all STDERR to STDOUT using sink() see http://mazamascience.com/WorkingWithData/?p=888
24 msg <- file(args.outLogFile, open="wt")
25 sink(msg, type="message")
26 sink(msg, type="output")
27 }
28
29 cat("\nSettings used===============:\n")
30 cat(readChar(args.settings, 1e5))
31
32
33 tryCatch(
34 {
35 library(metaMS)
36
37 ## load the data files from a zip file
38 files <- unzip(args.dataZip, exdir=args.zipExtrDir)
39
40 ## load settings "script" into "customMetaMSsettings"
41 tempEnv <- new.env()
42 source(args.settings, local=tempEnv)
43 message(paste(" loaded : ", args.settings))
44 allSettings <- tempEnv[["customMetaMSsettings"]]
45
46 # trigger runLC:
47 LC <- runLC(files, settings = allSettings, nSlaves=20, returnXset = TRUE)
48
49 # write out runLC annotation results:
50 write.table(LC$PeakTable, args.outPeakTable, sep="\t", row.names=FALSE)
51
52 # save xset as rdata:
53 xsAnnotatePreparedData <- LC$xset
54 saveRDS(xsAnnotatePreparedData, file=args.xsetOut)
55
56 message("\nGenerating report.........")
57 # report
58 dir.create(file.path(args.htmlReportFile.files_path), showWarnings = FALSE, recursive = TRUE)
59 html <- "<html><body><h1>Info on alignment quality </h1>"
60 # TODO add (nr and mass error) and group size
61
62 message("\nPlotting figures... ")
63 figureName <- paste(args.htmlReportFile.files_path, "/figure_retcor.png", sep="")
64 html <- paste(html,"<img src='figure_retcor.png' /><br/>", sep="")
65 png( figureName, type="cairo", width=1100,height=600 )
66 retcor(LC$xset@xcmsSet, method="peakgroups", plottype = "mdevden")
67 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,
68 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
69 this option is added in metaMS API. </a><br/>", sep="")
70 devname = dev.off()
71
72
73 gt <- groups(LC$xset@xcmsSet)
74 groupidx1 <- which(gt[,"rtmed"] > 0 & gt[,"rtmed"] < 3000 & gt[,"npeaks"] > 3)
75
76 html <- paste(html,"</body><html>")
77 message("finished generating report")
78 write(html,file=args.htmlReportFile)
79 # unlink(args.htmlReportFile)
80 cat("\nWarnings================:\n")
81 str( warnings() )
82 },
83 error=function(cond) {
84 sink(NULL, type="message") # default setting
85 sink(stderr(), type="output")
86 message("\nERROR: ===========\n")
87 print(cond)
88 }
89 )