comparison nmr_annotation2d/annotationRmn2DWrapper.R @ 0:8035235e46c7 draft

Uploaded
author marie-tremblay-metatoul
date Mon, 23 Dec 2019 09:26:20 -0500
parents
children
comparison
equal deleted inserted replaced
-1:000000000000 0:8035235e46c7
1 #!/usr/local/public/bin/Rscript --vanilla --slave --no-site-file
2
3 ## 201919016 2DNmrAnnotation_1.0.0.R
4 ## Marie Tremblay-Franco
5 ## MetaboHUB: The French Infrastructure for Metabolomics and Fluxomics
6 ## www.metabohub.fr/en
7 ## marie.tremblay-franco@toulouse.inra.fr
8
9 runExampleL <- FALSE
10
11 if(runExampleL) {
12 ##------------------------------
13 ## Example of arguments
14 ##------------------------------
15 }
16
17
18 ##------------------------------
19 ## Options
20 ##------------------------------
21 strAsFacL <- options()$stringsAsFactors
22 options(stringsAsFactors = FALSE)
23
24 ##------------------------------
25 ## Constants
26 ##------------------------------
27 topEnvC <- environment()
28 flagC <- "\n"
29
30
31 ##-------------------------
32 ## Input parameters reading
33 ##-------------------------
34
35 ##------------------------------
36 ## R libraries laoding
37 ##------------------------------
38 library(batch)
39 library(dplyr)
40 library(ggplot2)
41 library(openxlsx)
42 library(stringr)
43 library(tidyr)
44
45 if(!runExampleL)
46 argLs <- parseCommandArgs(evaluate=FALSE)
47 logFile <- argLs[["logOut"]]
48 sink(logFile)
49
50 cat("\tPACKAGE INFO\n")
51 sessionInfo()
52
53 ##------------------------------
54 ## Functions
55 ##------------------------------
56 source_local <- function(fname)
57 {
58 argv <- commandArgs(trailingOnly = FALSE)
59 base_dir <- dirname(substring(argv[grep("--file=", argv)], 8))
60 source(paste(base_dir, fname, sep="/"))
61 }
62 #Import the different functions
63 source_local("annotationRmn2D.R")
64 source_local("annotationRmn2DGlobale.R")
65 source_local("viridis.R")
66
67 ## Input parameter values
68 fileToAnnotate <- argLs[[1]]
69 # Chosen sequence(s)
70 cosy <- 0
71 hmbc <- 0
72 hsqc <- 0
73 jres <- 0
74 tocsy <- 0
75 ## sequences <- str_split(argLs[[2]], ",")[[1]]
76 ## for (s in 1:length(sequences))
77 ## {
78 ## argv <- commandArgs(trailingOnly = FALSE)
79 ## currentDir <- dirname(substring(argv[grep("--file=", argv)], 8))
80 ## if (sequences[s]=="cosy"){
81 ## cosy <- 1
82 ## load(paste(currentDir, "BdDReference_COSY.RData", sep="/"))
83 ## }else if(sequences[s]=="hmbc"){
84 ## hmbc <- 1
85 ## load(paste(currentDir, "BdDReference_HMBC.RData", sep="/"))
86 ## }else if(sequences[s]=="hsqc"){
87 ## hsqc <- 1
88 ## load(paste(currentDir, "BdDReference_HSQC.RData", sep="/"))
89 ## }else if(sequences[s]=="jres"){
90 ## jres <- 1
91 ## load(paste(currentDir, "BdDReference_JRES.RData", sep="/"))
92 ## }else if(sequences[s]=="tocsy"){
93 ## tocsy <- 1
94 ## load(paste(currentDir, "BdDReference_TOCSY.RData", sep="/"))
95 ## }else
96 ## stop("No chosen sequence", call.=FALSE)
97 ## }
98
99 if (argLs[[2]]=='yes')
100 {
101 argv <- commandArgs(trailingOnly = FALSE)
102 currentDir <- dirname(substring(argv[grep("--file=", argv)], 8))
103 cosy <- 1
104 load(paste(currentDir, "BdDReference_COSY.RData", sep="/"))
105 }
106
107 if (argLs[[3]]=='yes')
108 {
109 argv <- commandArgs(trailingOnly = FALSE)
110 currentDir <- dirname(substring(argv[grep("--file=", argv)], 8))
111 jres <- 1
112 load(paste(currentDir, "BdDReference_JRES.RData", sep="/"))
113 }
114
115 if (argLs[[4]]=='yes')
116 {
117 argv <- commandArgs(trailingOnly = FALSE)
118 currentDir <- dirname(substring(argv[grep("--file=", argv)], 8))
119 hmbc <- 1
120 load(paste(currentDir, "BdDReference_HMBC.RData", sep="/"))
121 }
122
123 if (argLs[[5]]=='yes')
124 {
125 argv <- commandArgs(trailingOnly = FALSE)
126 currentDir <- dirname(substring(argv[grep("--file=", argv)], 8))
127 hsqc <- 1
128 load(paste(currentDir, "BdDReference_HSQC.RData", sep="/"))
129 }
130
131 if (argLs[[6]]=='yes')
132 {
133 argv <- commandArgs(trailingOnly = FALSE)
134 currentDir <- dirname(substring(argv[grep("--file=", argv)], 8))
135 tocsy <- 1
136 load(paste(currentDir, "BdDReference_TOCSY.RData", sep="/"))
137 }
138
139 if (argLs[[2]]=='no' & argLs[[3]]=='no' & argLs[[4]]=='no' & argLs[[5]]=='no' & argLs[[6]]=='no')
140 stop("No chosen sequence", call.=FALSE)
141
142
143 # User database
144
145
146 # Allowed chemical shifts
147 tolPpm1 <- argLs$tolppm1
148 tolPpm2HJRes <- argLs$tolppmJRES
149 tolPpm2C <- argLs$tolppm2
150 # Threshold to remove metabolites (probability score < threshold)
151 seuil <- argLs$threshold
152 # Remove metabolites when multiple assignations?
153 unicite <- str_to_upper(argLs$unicity)
154
155 ## Output paramater values
156 AnnotationGraph <- argLs[["AnnotationGraph"]]
157
158 print(argLs)
159
160 ## ANNOTATION
161 st0=Sys.time()
162 pdf(AnnotationGraph,onefile=TRUE)
163 annotationMelange <- annotationRmn2DGlobale(fileToAnnotate, tolPpm1=tolPpm1, tolPpm2HJRes=tolPpm2HJRes,
164 tolPpm2C=tolPpm2C, cosy=cosy, hmbc=hmbc, hsqc=hsqc,
165 jres=jres, tocsy=tocsy, seuil=seuil, unicite=unicite)
166 dev.off()
167
168 if (cosy==1)
169 {
170 write.table(annotationMelange$COSY$liste_resultat, file=argLs[["annotationCOSY"]], quote=FALSE,
171 row.names=FALSE,sep="\t")
172 if (nrow(annotationMelange$COSY$listing_ppm_commun) != 0)
173 write.table(annotationMelange$COSY$listing_ppm_commun, file=argLs[["ppmCommunCOSY"]], quote=FALSE,
174 row.names=FALSE,sep="\t")
175 }
176
177 if (hmbc==1)
178 {
179 write.table(annotationMelange$HMBC$liste_resultat, file=argLs[["annotationHMBC"]], quote=FALSE,
180 row.names=FALSE,sep="\t")
181 if (nrow(annotationMelange$HMBC$listing_ppm_commun) != 0)
182 write.table(annotationMelange$HMBC$listing_ppm_commun, file=argLs[["ppmCommunHMBC"]], quote=FALSE,
183 row.names=FALSE,sep="\t")
184 }
185
186 if (hsqc==1)
187 {
188 write.table(annotationMelange$HSQC$liste_resultat, file=argLs[["annotationHSQC"]], quote=FALSE,
189 row.names=FALSE,sep="\t")
190 if (nrow(annotationMelange$HSQC$listing_ppm_commun) != 0)
191 write.table(annotationMelange$HSQC$listing_ppm_commun, file=argLs[["ppmCommunHSQC"]], quote=FALSE,
192 row.names=FALSE,sep="\t")
193 }
194
195 if (jres==1)
196 {
197 write.table(annotationMelange$JRES$liste_resultat, file=argLs[["annotationJRES"]], quote=FALSE,
198 row.names=FALSE,sep="\t")
199 if (nrow(annotationMelange$JRES$listing_ppm_commun) != 0)
200 write.table(annotationMelange$JRES$listing_ppm_commun, file=argLs[["ppmCommunJRES"]], quote=FALSE,
201 row.names=FALSE,sep="\t")
202 }
203
204 if (tocsy==1)
205 {
206 write.table(annotationMelange$TOCSY$liste_resultat, file=argLs[["annotationTOCSY"]], quote=FALSE,
207 row.names=FALSE,sep="\t")
208 if (nrow(annotationMelange$TOCSY$listing_ppm_commun) != 0)
209 write.table(annotationMelange$TOCSY$listing_ppm_commun, file=argLs[["ppmCommunTOCSY"]], quote=FALSE,
210 row.names=FALSE,sep="\t")
211 }
212
213 ## Combinaison de sequences
214 if (cosy + jres + hmbc + hsqc + tocsy > 1)
215 {
216 write.table(annotationMelange$combination, file=argLs[["annotationCombination"]], quote=FALSE,
217 row.names=FALSE,sep="\t")
218 }
219 st1=Sys.time()
220 print(st1-st0)
221
222 ## Ending
223 ##--------
224 cat("\nEnd of '2D NMR annotation' Galaxy module call: ", as.character(Sys.time()), sep = "")
225 sink()
226 options(stringsAsFactors = strAsFacL)
227 rm(list = ls())