Mercurial > repos > computational-metabolomics > mspurity_spectralmatching

comparison filterFragSpectra.R @ 6:fecfe8c80e25 draft

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
"planemo upload for repository https://github.com/computational-metabolomics/mspurity-galaxy commit 2579c8746819670348c378f86116f83703c493eb"
author computational-metabolomics
date Thu, 04 Mar 2021 12:32:11 +0000
parents 5ff9d40c7a42
children 20f4fdaf36bb
comparison
equal deleted inserted replaced
5:829af6fecd96 6:fecfe8c80e25
3 library(xcms) 3 library(xcms)
4 print(sessionInfo()) 4 print(sessionInfo())
5 5
6 6
7 option_list <- list( 7 option_list <- list(
8 make_option("--out_rdata", type="character"), 8 make_option("--out_rdata", type = "character"),
9 make_option("--out_peaklist_prec", type="character"), 9 make_option("--out_peaklist_prec", type = "character"),
10 make_option("--out_peaklist_frag", type="character"), 10 make_option("--out_peaklist_frag", type = "character"),
11 make_option("--pa", type="character"), 11 make_option("--pa", type = "character"),
12 12
13 make_option("--ilim", default=0.0), 13 make_option("--ilim", default = 0.0),
14 make_option("--plim", default=0.0), 14 make_option("--plim", default = 0.0),
15 15
16 make_option("--ra", default=0.0), 16 make_option("--ra", default = 0.0),
17 make_option("--snr", default=0.0), 17 make_option("--snr", default = 0.0),
18 18
19 make_option("--rmp", action="store_true"), 19 make_option("--rmp", action = "store_true"),
20 make_option("--snmeth", default="median", type="character") 20 make_option("--snmeth", default = "median", type = "character")
21 ) 21 )
22 22
23 opt <- parse_args(OptionParser(option_list=option_list)) 23 opt <- parse_args(OptionParser(option_list = option_list))
24 print(opt) 24 print(opt)
25 25
26 26
27 loadRData <- function(rdata_path, name){ 27 loadRData <- function(rdata_path, name) {
28 #loads an RData file, and returns the named xset object if it is there 28 #loads an RData file, and returns the named xset object if it is there
29 load(rdata_path) 29 load(rdata_path)
30 return(get(ls()[ls() %in% name])) 30 return(get(ls()[ls() %in% name]))
31 } 31 }
32 32
33 # Requires 33 # Requires
34 pa <- loadRData(opt$pa, 'pa') 34 pa <- loadRData(opt$pa, "pa")
35 35
36 if(is.null(opt$rmp)){ 36 if (is.null(opt$rmp)) {
37 opt$rmp = FALSE 37 opt$rmp <- FALSE
38 }else{ 38 }else{
39 opt$rmp = TRUE 39 opt$rmp <- TRUE
40 } 40 }
41 41
42 pa <- filterFragSpectra(pa, 42 pa <- filterFragSpectra(pa,
43 ilim=opt$ilim, 43 ilim = opt$ilim,
44 plim=opt$plim, 44 plim = opt$plim,
45 ra=opt$ra, 45 ra = opt$ra,
46 snr=opt$snr, 46 snr = opt$snr,
47 rmp=opt$rmp, 47 rmp = opt$rmp,
48 snmeth=opt$snmeth) 48 snmeth = opt$snmeth)
49 49
50 print(pa) 50 print(pa)
51 save(pa, file=opt$out_rdata) 51 save(pa, file = opt$out_rdata)
52 52
53 # get the msms data for grpid from the purityA object 53 # get the msms data for grpid from the purityA object
54 msmsgrp <- function(grpid, pa){ 54 msmsgrp <- function(grpid, pa) {
55 msms <- pa@grped_ms2[grpid] 55 msms <- pa@grped_ms2[grpid]
56 56
57 grpinfo <- pa@grped_df[pa@grped_df$grpid==grpid,] 57 grpinfo <- pa@grped_df[pa@grped_df$grpid == grpid, ]
58 58
59 grpinfo$subsetid <- 1:nrow(grpinfo) 59 grpinfo$subsetid <- seq_len(nrow(grpinfo))
60 result <- plyr::ddply(grpinfo, ~subsetid, setid, msms=msms) 60 result <- plyr::ddply(grpinfo, ~subsetid, setid, msms = msms)
61 return(result) 61 return(result)
62 } 62 }
63 63
64 # Set the relevant details 64 # Set the relevant details
65 setid <- function(grpinfo_i, msms){ 65 setid <- function(grpinfo_i, msms) {
66 msms_i <- msms[[1]][[grpinfo_i$subsetid]] 66 msms_i <- msms[[1]][[grpinfo_i$subsetid]]
67 n <- nrow(msms_i) 67 n <- nrow(msms_i)
68 msms_i <- data.frame(msms_i) 68 msms_i <- data.frame(msms_i)
69 colnames(msms_i)[1:2] <- c('mz', 'i') 69 colnames(msms_i)[1:2] <- c("mz", "i")
70 m <- cbind('grpid'=rep(grpinfo_i$grpid,n), 'pid'=rep(grpinfo_i$pid,n), 'fileid'=rep(grpinfo_i$fileid,n), msms_i) 70 m <- cbind("grpid" = rep(grpinfo_i$grpid, n), "pid" = rep(grpinfo_i$pid, n), "fileid" = rep(grpinfo_i$fileid, n), msms_i)
71 return(m) 71 return(m)
72 } 72 }
73 73
74 74
75 75
76 if (length(pa)>0){ 76 if (length(pa) > 0) {
77 77
78 if (length(pa@grped_ms2)==0){ 78 if (length(pa@grped_ms2) == 0) {
79 message('No spectra available') 79 message("No spectra available")
80 } else{ 80 } else {
81 81
82 # get group ids 82 # get group ids
83 grpids <- unique(as.character(pa@grped_df$grpid)) 83 grpids <- unique(as.character(pa@grped_df$grpid))
84 84
85 # loop through all the group ids 85 # loop through all the group ids
86 df_fragments = plyr::adply(grpids, 1, msmsgrp, pa=pa) 86 df_fragments <- plyr::adply(grpids, 1, msmsgrp, pa = pa)
87 df_fragments = merge(df_fragments, pa@puritydf[,c("pid", "acquisitionNum", "precursorScanNum")], by="pid") 87 df_fragments <- merge(df_fragments, pa@puritydf[, c("pid", "acquisitionNum", "precursorScanNum")], by = "pid")
88 df_fragments = df_fragments[order(df_fragments$grpid, df_fragments$pid, df_fragments$mz),] 88 df_fragments <- df_fragments[order(df_fragments$grpid, df_fragments$pid, df_fragments$mz), ]
89 #select and reorder columns 89 #select and reorder columns
90 df_fragments = df_fragments[,c("grpid", "pid", "precursorScanNum", "acquisitionNum", "fileid", "mz", "i", "snr", "ra", "purity_pass_flag", "intensity_pass_flag", "ra_pass_flag", "snr_pass_flag", "pass_flag")] 90 df_fragments <- df_fragments[, c("grpid", "pid", "precursorScanNum", "acquisitionNum", "fileid", "mz", "i", "snr", "ra", "purity_pass_flag", "intensity_pass_flag", "ra_pass_flag", "snr_pass_flag", "pass_flag")]
91 91
92 pa@grped_df$filename = sapply(pa@grped_df$fileid, function(x) names(pa@fileList)[as.integer(x)]) 92 pa@grped_df$filename <- sapply(pa@grped_df$fileid, function(x) names(pa@fileList)[as.integer(x)])
93 93
94 print(head(pa@grped_df)) 94 print(head(pa@grped_df))
95 write.table(pa@grped_df, opt$out_peaklist_prec, row.names=FALSE, sep='\t') 95 write.table(pa@grped_df, opt$out_peaklist_prec, row.names = FALSE, sep = "\t")
96 print(head(df_fragments)) 96 print(head(df_fragments))
97 write.table(df_fragments, opt$out_peaklist_frag, row.names=FALSE, sep='\t') 97 write.table(df_fragments, opt$out_peaklist_frag, row.names = FALSE, sep = "\t")
98 } 98 }
99 } 99 }
100