comparison CAMERA_findAdducts.R @ 0:3b2f7cb74c33 draft default tip

planemo upload commit 24d44ee26b7c23380c2b42fae2f7f6e58472100d
author workflow4metabolomics
date Sun, 24 Nov 2024 21:31:14 +0000
parents
children
comparison
equal deleted inserted replaced
-1:000000000000 0:3b2f7cb74c33
1 #!/usr/bin/env Rscript
2
3 # ----- PACKAGE -----
4 cat("\tSESSION INFO\n")
5
6 # Import the different functions
7 source_local <- function(fname) {
8 argv <- commandArgs(trailingOnly = FALSE)
9 base_dir <- dirname(substring(argv[grep("--file=", argv)], 8))
10 source(paste(base_dir, fname, sep = "/"))
11 }
12 source_local("lib.r")
13
14 pkgs <- c("CAMERA", "xcms", "multtest", "batch")
15 loadAndDisplayPackages(pkgs)
16 cat("\n\n")
17 # ----- ARGUMENTS -----
18 cat("\tARGUMENTS INFO\n")
19
20 args <- parseCommandArgs(evaluate = FALSE) # interpretation of arguments given in command line as an R list of objects
21 write.table(as.matrix(args), col.names = FALSE, quote = FALSE, sep = "\t")
22
23 cat("\n\n")
24
25 print("Arguments retrieved from the command line:")
26 print(args)
27
28 # Function to convert "NULL" strings to actual NULL values
29 convertNullString <- function(x) {
30 if (x == "NULL") {
31 return(NULL)
32 }
33 return(x)
34 }
35
36 # Function to convert string to numeric lists
37 convert_psg_list <- function(x) {
38 # Check if x is NULL or has zero length before further processing
39 if (is.null(x) || length(x) == 0) {
40 return(NULL)
41 }
42
43 # Force conversion to character
44 x <- as.character(x)
45
46 if (grepl("^[0-9]+$", x)) {
47 # If the string represents a single numeric value
48 return(as.numeric(x))
49 } else {
50 # Convert string representation of a list to a numeric vector
51 # Use a regular expression to split by common separators
52 return(as.numeric(unlist(strsplit(x, "[,;\\s]+"))))
53 }
54 }
55
56 for (arg in names(args)) {
57 args[[arg]] <- convertNullString(args[[arg]])
58 }
59
60 # Convert only the 'psg_list' element in args
61 args$psg_list <- convert_psg_list(args$psg_list)
62
63 print("Argument types:")
64 print(sapply(args, class))
65
66 # Check if the image file exists
67 if (!file.exists(args$image)) {
68 stop("The RData file does not exist: ", args$image)
69 }
70
71 # ----- PROCESSING INFILE -----
72
73 # Load the RData file
74 load(args$image)
75 args$image <- NULL
76
77 # Save arguments for report generation
78 if (!exists("listOFlistArguments")) listOFlistArguments <- list()
79 listOFlistArguments[[format(Sys.time(), "%y%m%d-%H:%M:%S_findAdducts")]] <- args
80
81 # We unzip automatically the chromatograms from the zip files.
82 if (!exists("zipfile")) zipfile <- NULL
83 if (!exists("singlefile")) singlefile <- NULL
84 rawFilePath <- getRawfilePathFromArguments(singlefile, zipfile, args)
85 zipfile <- rawFilePath$zipfile
86 singlefile <- rawFilePath$singlefile
87 args <- rawFilePath$args
88
89 print(paste("singlefile :", singlefile))
90 if (!is.null(singlefile)) {
91 directory <- retrieveRawfileInTheWorkingDir(singlefile, zipfile)
92 }
93
94 # Check if the 'rules' argument in 'args' is NULL
95 if (is.null(args$rules)) {
96 # If 'args$rules' is NULL, set 'rulset' to NULL
97 args$rulset <- NULL
98 } else {
99 # Try to read the rules file with different delimiters
100 delimiters <- c(";", "\t", ",") # List of possible delimiters
101 success <- FALSE # Flag to check if reading was successful
102
103 for (sep in delimiters) {
104 # Attempt to read the rules file with the current separator
105 args$rulset <- read.table(args$rules, header = TRUE, sep = sep)
106
107 # Check if the number of columns is at least 4
108 if (ncol(args$rulset) >= 4) {
109 success <- TRUE # Mark success if the format is correct
110 break # Exit the loop if the file was read successfully
111 }
112 }
113
114 # If reading the rules file failed for all delimiters
115 if (!success) {
116 # Display an error message if the file is not well formatted
117 error_message <- "The rules file appears to be improperly formatted. Accepted column separators are ;, tab, and ,."
118 print(error_message)
119 stop(error_message) # Stop execution with an error
120 }
121 }
122
123 # Verify that the object xa is loaded
124 if (!exists("xa")) {
125 stop("The object xa was not found in the RData file.")
126 }
127
128 print("Loaded xa object:")
129 print(xa)
130
131 # Apply the findAdducts function on the xsAnnotate object
132 print("Calling findAdducts function:")
133 xa <- findAdducts(xa, ppm = args$ppm, mzabs = args$mzabs, multiplier = args$multiplier, polarity = args$polarity, rules = args$rulset, max_peaks = args$max_peaks, psg_list = args$psg_list, intval = args$intval)
134
135 print("Result of findAdducts function:")
136 print(xa)
137
138 # Extract the list of annotated peaks
139 peakList <- getPeaklist(xa, intval = args$intval)
140
141 if (length(phenoData@data$sample_name) == 1) {
142 peakList$name <- make.unique(paste0("M", round(peakList[, "mz"], 0), "T", round(peakList[, "rt"], 0)), "_")
143 variableMetadata <- peakList[, c("name", setdiff(names(peakList), "name"))]
144 variableMetadata <- formatIonIdentifiers(variableMetadata, numDigitsRT = args$numDigitsRT, numDigitsMZ = args$numDigitsMZ)
145 } else {
146 names_default <- groupnames(xa@xcmsSet, mzdec = 0, rtdec = 0) # Names without decimals
147 names_custom <- groupnames(xa@xcmsSet, mzdec = args$numDigitsMZ, rtdec = args$numDigitsRT) # Names with "x" decimals
148
149 variableMetadata <- data.frame(
150 name = names_default,
151 name_custom = names_custom,
152 stringsAsFactors = FALSE
153 )
154 variableMetadata <- cbind(variableMetadata, peakList[, !(make.names(colnames(peakList)) %in% c(make.names(sampnames(xa@xcmsSet))))])
155 }
156
157 if (!exists("RTinMinute")) RTinMinute <- FALSE
158
159 if (args$convertRTMinute && RTinMinute == FALSE) {
160 RTinMinute <- TRUE
161 variableMetadata <- RTSecondToMinute(variableMetadata = variableMetadata, convertRTMinute = args$convertRTMinute)
162 }
163
164 # Save the extracted peak list as a TSV file named 'variableMetadata.tsv'
165 output_file_tsv <- "variableMetadata.tsv"
166 write.table(variableMetadata, file = output_file_tsv, sep = "\t", row.names = FALSE, quote = FALSE)
167
168 # Save the updated xsAnnotate object
169 output_file_RData <- "camera_findAdducts.RData"
170 objects2save <- c("xa", "variableMetadata", "listOFlistArguments", "zipfile", "singlefile", "RTinMinute", "phenoData")
171 save(list = objects2save[objects2save %in% ls()], file = output_file_RData)
172
173 cat("Output files generated:", output_file_tsv, "and", output_file_RData, "\n")