Mercurial > repos > tomnl > msp_split

comparison flag-remove-peaks.R @ 2:89f33758ad22 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 f79fa34772bbab836d89cf8bad52d49285409a98
author tomnl
date Thu, 14 Jun 2018 09:09:58 -0400
parents aa55ff1d76d4
children
comparison
equal deleted inserted replaced
1:0db4c33f5f5b 2:89f33758ad22
1 library(XCMSwrapper) 1 library(msPurity)
2 library(optparse) 2 library(optparse)
3 3 print(sessionInfo())
4 option_list <- list( 4 option_list <- list(
5 make_option(c("-o", "--out_dir"), type="character", default=getwd(), 5 make_option(c("-o", "--out_dir"), type="character", default=getwd(),
6 help="Output folder for resulting files [default = %default]" 6 help="Output folder for resulting files [default = %default]"
7 ), 7 ),
8 make_option(c("-x", "--xset_path"), type="character", default=file.path(getwd(),"xset.rds"), 8 make_option(c("-x", "--xset_path"), type="character", default=file.path(getwd(),"xset.rds"),
9 help="The path to the xcmsSet object [default = %default]" 9 help="The path to the xcmsSet object [default = %default]"
10 ), 10 ),
11 make_option("--polarity", default='NA', 11 make_option("--polarity", default=NA,
12 help="polarity (just used for naming purpose for files being saved) [positive, negative, NA] [default %default]" 12 help="polarity (just used for naming purpose for files being saved) [positive, negative, NA] [default %default]"
13 ), 13 ),
14 make_option("--rsd_i_blank", default=NA, 14 make_option("--rsd_i_blank", default=100,
15 help="RSD threshold for the blank [default = %default]" 15 help="RSD threshold for the blank [default = %default]"
16 ), 16 ),
17 make_option("--minfrac_blank", default=0.5, 17 make_option("--minfrac_blank", default=0.5,
18 help="minimum fraction of files for features needed for the blank [default = %default]" 18 help="minimum fraction of files for features needed for the blank [default = %default]"
19 ), 19 ),
20 make_option("--rsd_rt_blank", default=NA, 20 make_option("--rsd_rt_blank", default=100,
21 help="RSD threshold for the RT of the blank [default = %default]" 21 help="RSD threshold for the RT of the blank [default = %default]"
22 ), 22 ),
23 23
24 make_option("--ithres_blank", default=NA, 24 make_option("--ithres_blank", default=0,
25 help="Intensity threshold for the blank [default = %default]" 25 help="Intensity threshold for the blank [default = %default]"
26 ), 26 ),
27 make_option("--s2b", default=10, 27 make_option("--s2b", default=10,
28 help="fold change (sample/blank) needed for sample peak to be allowed. e.g. 28 help="fold change (sample/blank) needed for sample peak to be allowed. e.g.
29 if s2b set to 10 and the recorded sample 'intensity' value was 100 and blank was 10. 29 if s2b set to 10 and the recorded sample 'intensity' value was 100 and blank was 10.
36 make_option("--egauss_thr", default=NA, 36 make_option("--egauss_thr", default=NA,
37 help="Threshold for filtering out non gaussian shaped peaks. Note this only works 37 help="Threshold for filtering out non gaussian shaped peaks. Note this only works
38 if the 'verbose columns' and 'fit gauss' was used with xcms 38 if the 'verbose columns' and 'fit gauss' was used with xcms
39 [default = %default]" 39 [default = %default]"
40 ), 40 ),
41 make_option("--rsd_i_sample", default=NA, 41 make_option("--rsd_i_sample", default=100,
42 help="RSD threshold for the samples [default = %default]" 42 help="RSD threshold for the samples [default = %default]"
43 ), 43 ),
44 make_option("--minfrac_sample", default=0.8, 44 make_option("--minfrac_sample", default=0.8,
45 help="minimum fraction of files for features needed for the samples [default = %default]" 45 help="minimum fraction of files for features needed for the samples [default = %default]"
46 ), 46 ),
47 make_option("--rsd_rt_sample", default=NA, 47 make_option("--rsd_rt_sample", default=100,
48 help="RSD threshold for the RT of the samples [default = %default]" 48 help="RSD threshold for the RT of the samples [default = %default]"
49 ), 49 ),
50 make_option("--ithres_sample", default=5000, 50 make_option("--ithres_sample", default=5000,
51 help="Intensity threshold for the sample [default = %default]" 51 help="Intensity threshold for the sample [default = %default]"
52 ), 52 ),
86 #), 86 #),
87 87
88 # store options 88 # store options
89 opt<- parse_args(OptionParser(option_list=option_list)) 89 opt<- parse_args(OptionParser(option_list=option_list))
90 90
91 opt <- replace(opt, opt == "NA", NA)
92
93
94
95
91 if (is.null(opt$temp_save)){ 96 if (is.null(opt$temp_save)){
92 temp_save<-FALSE 97 temp_save<-FALSE
93 }else{ 98 }else{
94 temp_save<-TRUE 99 temp_save<-TRUE
95 } 100 }
126 print(blank_class) 131 print(blank_class)
127 } 132 }
128 133
129 134
130 if (is.null(opt$multilist)){ 135 if (is.null(opt$multilist)){
131 ffrm_out <- XCMSwrapper::flag_remove(xset, 136 ffrm_out <- flag_remove(xset,
132 pol=opt$polarity, 137 pol=opt$polarity,
133 rsd_i_blank=opt$rsd_i_blank, 138 rsd_i_blank=opt$rsd_i_blank,
134 minfrac_blank=opt$minfrac_blank, 139 minfrac_blank=opt$minfrac_blank,
135 rsd_rt_blank=opt$rsd_rt_blank, 140 rsd_rt_blank=opt$rsd_rt_blank,
136 ithres_blank=opt$ithres_blank, 141 ithres_blank=opt$ithres_blank,
146 bw=opt$bw, 151 bw=opt$bw,
147 out_dir=opt$out_dir, 152 out_dir=opt$out_dir,
148 temp_save=temp_save, 153 temp_save=temp_save,
149 remove_spectra=remove_spectra, 154 remove_spectra=remove_spectra,
150 grp_rm_ids=unlist(strsplit(as.character(opt$grp_rm_ids), split=", "))[[1]]) 155 grp_rm_ids=unlist(strsplit(as.character(opt$grp_rm_ids), split=", "))[[1]])
151 156 print('flag remove finished')
152 xset <- ffrm_out[[1]] 157 xset <- ffrm_out[[1]]
153 grp_peaklist <- ffrm_out[[2]] 158 grp_peaklist <- ffrm_out[[2]]
154 removed_peaks <- ffrm_out[[3]] 159 removed_peaks <- ffrm_out[[3]]
155 160
156 save.image(file=file.path(opt$out_dir, 'xset_filtered.RData')) 161 save.image(file=file.path(opt$out_dir, 'xset_filtered.RData'))
157 162
158 # grpid needed for mspurity ID needed for deconrank... (will clean up at some up) 163 # grpid needed for mspurity ID needed for deconrank... (will clean up at some up)
164 peak_pth <- file.path(opt$out_dir, 'peaklist_filtered.tsv')
165 print(peak_pth)
159 write.table(data.frame('grpid'=rownames(grp_peaklist), 'ID'=rownames(grp_peaklist), grp_peaklist), 166 write.table(data.frame('grpid'=rownames(grp_peaklist), 'ID'=rownames(grp_peaklist), grp_peaklist),
160 file.path(opt$out_dir, 'peaklist_filtered.txt'), row.names=FALSE, sep='\t') 167 peak_pth, row.names=FALSE, sep='\t')
161 168
162 removed_peaks <- data.frame(removed_peaks) 169 removed_peaks <- data.frame(removed_peaks)
163 write.table(data.frame('ID'=rownames(removed_peaks),removed_peaks), 170 write.table(data.frame('ID'=rownames(removed_peaks),removed_peaks),
164 file.path(opt$out_dir, 'removed_peaks.txt'), row.names=FALSE, sep='\t') 171 file.path(opt$out_dir, 'removed_peaks.tsv'), row.names=FALSE, sep='\t')
165 172
166 }else{ 173 }else{
167 174
168 175
169 # TODO 176 # TODO