Mercurial > repos > recetox > isolib

comparison isolib.R @ 3:6b0fef8a77c0 draft

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
planemo upload for repository https://github.com/RECETOX/galaxytools/tree/master/tools/isolib commit bc3445f7c41271b0062c7674108f57708d08dd28
author recetox
date Thu, 30 May 2024 14:52:02 +0000
parents b3251a7dae25
children 2b1118bce0b1
comparison
equal deleted inserted replaced
2:b3251a7dae25 3:6b0fef8a77c0
4 library(MetaboCoreUtils) 4 library(MetaboCoreUtils)
5 library(readr) 5 library(readr)
6 6
7 #' @param args A list of command line arguments. 7 #' @param args A list of command line arguments.
8 main <- function() { 8 main <- function() {
9 data(isotopes) 9 data(isotopes)
10 data(adducts) 10 data(adducts)
11 11
12 args <- commandArgs(trailingOnly = TRUE) 12 args <- commandArgs(trailingOnly = TRUE)
13 compound_table <- read_tsv( 13 compound_table <- read_tsv(
14 file = args[1], 14 file = args[1],
15 col_types = "ccd", 15 col_types = "ccd",
16 col_select = tidyselect::all_of(c("name", "formula")) | tidyselect::any_of("rt") 16 col_select = tidyselect::all_of(c("name", "formula")) | tidyselect::any_of("rt")
17 ) 17 )
18 adducts_to_use <- c(unlist(strsplit(args[2], ",", fixed = TRUE))) 18 adducts_to_use <- c(unlist(strsplit(args[2], ",", fixed = TRUE)))
19 19
20 chemforms <- compound_table$formula 20 chemforms <- compound_table$formula
21 chemforms <- check_chemform(isotopes, chemforms)[, 2] 21 chemforms <- check_chemform(isotopes, chemforms)[, 2]
22 22
23 spectra <- data.frame() 23 spectra <- data.frame()
24 24
25 for (current in adducts_to_use) { 25 for (current in adducts_to_use) {
26 adduct <- adducts[adducts$Name == current, ] 26 adduct <- adducts[adducts$Name == current, ]
27 multiplied_chemforms <- multiform(chemforms, adduct$Mult) 27 multiplied_chemforms <- multiform(chemforms, adduct$Mult)
28 28
29 if (adduct$Ion_mode == "negative") { 29 if (adduct$Ion_mode == "negative") {
30 merged_chemforms <- subform(multiplied_chemforms, adduct$Formula_ded) 30 merged_chemforms <- subform(multiplied_chemforms, adduct$Formula_ded)
31 } else { 31 } else {
32 merged_chemforms <- mergeform(multiplied_chemforms, adduct$Formula_add) 32 merged_chemforms <- mergeform(multiplied_chemforms, adduct$Formula_add)
33 }
34
35 charge_string <- paste0(if (adduct$Charge > 0) "+" else "-", if (abs(adduct$Charge) > 1) abs(adduct$Charge) else "")
36 adduct_string <- paste0("[", adduct$Name, "]", charge_string)
37 precursor_mz <- calculateMass(multiplied_chemforms) + adduct$Mass
38
39 if (args[4] == TRUE) {
40 names <- paste(compound_table$name, paste0("(", adduct$Name, ")"), sep = " ")
41 } else {
42 names <- compound_table$name
43 }
44
45 spectra_df <- data.frame(
46 name = names,
47 adduct = adduct_string,
48 formula = chemforms,
49 charge = adduct$Charge,
50 ionization_mode = adduct$Ion_mode,
51 precursor_mz = precursor_mz,
52 msLevel = as.integer(1)
53 )
54
55 if ("rt" %in% colnames(compound_table)) {
56 spectra_df$retention_time <- compound_table$rt
57 }
58
59 patterns <- enviPat::isopattern(
60 isotopes = isotopes,
61 chemforms = merged_chemforms,
62 charge = adduct$Charge,
63 threshold = as.numeric(args[3]),
64 )
65
66 mzs <- list()
67 intensities <- list()
68 for (i in seq_along(patterns)) {
69 mzs <- append(mzs, list(patterns[[i]][, 1]))
70 intensities <- append(intensities, list(patterns[[i]][, 2]))
71 }
72
73 spectra_df$mz <- mzs
74 spectra_df$intensity <- intensities
75 spectra <- rbind(spectra, spectra_df)
33 } 76 }
34 77
35 charge_string <- paste0(if (adduct$Charge > 0) "+" else "-", if (abs(adduct$Charge) > 1) abs(adduct$Charge) else "") 78 sps <- Spectra(spectra)
36 adduct_string <- paste0("[", adduct$Name, "]", charge_string) 79 export(sps, MsBackendMsp(), file = args[5])
37 precursor_mz <- calculateMass(multiplied_chemforms) + adduct$Mass
38
39 if (args[4] == TRUE) {
40 names <- paste(compound_table$name, paste0("(", adduct$Name, ")"), sep = " ")
41 } else {
42 names <- compound_table$name
43 }
44
45 spectra_df <- data.frame(
46 name = names,
47 adduct = adduct_string,
48 formula = chemforms,
49 charge = adduct$Charge,
50 ionization_mode = adduct$Ion_mode,
51 precursor_mz = precursor_mz,
52 msLevel = as.integer(1)
53 )
54
55 if ("rt" %in% colnames(compound_table)) {
56 spectra_df$retention_time <- compound_table$rt
57 }
58
59 patterns <- enviPat::isopattern(
60 isotopes = isotopes,
61 chemforms = merged_chemforms,
62 charge = adduct$Charge,
63 threshold = as.numeric(args[3]),
64 )
65
66 mzs <- list()
67 intensities <- list()
68 for (i in seq_along(patterns)) {
69 mzs <- append(mzs, list(patterns[[i]][, 1]))
70 intensities <- append(intensities, list(patterns[[i]][, 2]))
71 }
72
73 spectra_df$mz <- mzs
74 spectra_df$intensity <- intensities
75 spectra <- rbind(spectra, spectra_df)
76 }
77
78 sps <- Spectra(spectra)
79 export(sps, MsBackendMsp(), file = args[5])
80 } 80 }
81 81
82 # Get the command line arguments 82 # Get the command line arguments
83 args <- commandArgs(trailingOnly = TRUE) 83 args <- commandArgs(trailingOnly = TRUE)
84 # Call the main function 84 # Call the main function