Mercurial > repos > computational-metabolomics > mspurity_createmsp

comparison averageFragSpectra.R @ 6:d25273689e04 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:31:13 +0000
parents 35898942bfbb
children b91b9492a4bf
comparison
equal deleted inserted replaced
5:50409aeeeac5 6:d25273689e04
2 library(msPurity) 2 library(msPurity)
3 library(xcms) 3 library(xcms)
4 print(sessionInfo()) 4 print(sessionInfo())
5 5
6 6
7 get_av_spectra <- function(x){ 7 get_av_spectra <- function(x) {
8 8
9 if (length(x$av_intra)>0){ 9 if (length(x$av_intra) > 0) {
10 av_intra_df <- plyr::ldply(x$av_intra) 10 av_intra_df <- plyr::ldply(x$av_intra)
11 11
12 if (nrow(av_intra_df)==0){ 12 if (nrow(av_intra_df) == 0) {
13 av_intra_df <- NULL 13 av_intra_df <- NULL
14 }else{ 14 }else{
15 av_intra_df$method <- 'intra' 15 av_intra_df$method <- "intra"
16 } 16 }
17 17
18 }else{ 18 }else{
19 av_intra_df <- NULL 19 av_intra_df <- NULL
20 } 20 }
21 21
22 if ((is.null(x$av_inter)) || (nrow(x$av_inter)==0)){ 22 if ((is.null(x$av_inter)) || (nrow(x$av_inter) == 0)) {
23 av_inter_df <- NULL 23 av_inter_df <- NULL
24 }else{ 24 }else{
25 av_inter_df <- x$av_inter 25 av_inter_df <- x$av_inter
26 av_inter_df$method <- 'inter' 26 av_inter_df$method <- "inter"
27 } 27 }
28 28
29 if ((is.null(x$av_all)) || (nrow(x$av_all)==0)){ 29 if ((is.null(x$av_all)) || (nrow(x$av_all) == 0)) {
30 av_all_df <- NULL 30 av_all_df <- NULL
31 }else{ 31 }else{
32 av_all_df <- x$av_all 32 av_all_df <- x$av_all
33 av_all_df$method <- 'all' 33 av_all_df$method <- "all"
34 } 34 }
35 35
36 combined <- plyr::rbind.fill(av_intra_df, av_inter_df, av_all_df) 36 combined <- plyr::rbind.fill(av_intra_df, av_inter_df, av_all_df)
37 37
38 return(combined) 38 return(combined)
39 } 39 }
40 40
41 41
42 option_list <- list( 42 option_list <- list(
43 make_option("--out_rdata", type="character"), 43 make_option("--out_rdata", type = "character"),
44 make_option("--out_peaklist", type="character"), 44 make_option("--out_peaklist", type = "character"),
45 make_option("--pa", type="character"), 45 make_option("--pa", type = "character"),
46 46 make_option("--av_level", type = "character"),
47 make_option("--av_level", type="character"), 47 make_option("--minfrac", default = 0.5),
48 48 make_option("--minnum", default = 1),
49 make_option("--minfrac", default=0.5), 49 make_option("--ppm", default = 5.0),
50 make_option("--minnum", default=1), 50 make_option("--snr", default = 0),
51 make_option("--ppm", default=5.0), 51 make_option("--ra", default = 0),
52 52 make_option("--av", default = "median", type = "character"),
53 make_option("--snr", default=0), 53 make_option("--sumi", action = "store_true"),
54 54 make_option("--rmp", action = "store_true"),
55 make_option("--ra", default=0), 55 make_option("--cores", default = 1)
56
57 make_option("--av", default="median", type="character"),
58 make_option("--sumi", action="store_true"),
59
60 make_option("--rmp", action="store_true"),
61 make_option("--cores", default=1)
62 ) 56 )
63 57
64 opt <- parse_args(OptionParser(option_list=option_list)) 58 opt <- parse_args(OptionParser(option_list = option_list))
65 print(opt) 59 print(opt)
66 60
67 61
68 loadRData <- function(rdata_path, name){ 62 load_r_data <- function(rdata_path, name) {
69 #loads an RData file, and returns the named xset object if it is there 63 #loads an RData file, and returns the named xset object if it is there
70 load(rdata_path) 64 load(rdata_path)
71 return(get(ls()[ls() %in% name])) 65 return(get(ls()[ls() %in% name]))
72 } 66 }
73 67
74 # Requires 68 # Requires
75 pa <- loadRData(opt$pa, 'pa') 69 pa <- load_r_data(opt$pa, "pa")
76 70
77 pa@cores <- opt$cores 71 pa@cores <- opt$cores
78 72
79 if(is.null(opt$rmp)){ 73 if (is.null(opt$rmp)) {
80 rmp = FALSE 74 rmp <- FALSE
81 }else{ 75 }else{
82 rmp = TRUE 76 rmp <- TRUE
83 } 77 }
84 78
85 if(is.null(opt$sumi)){ 79 if (is.null(opt$sumi)) {
86 80 sumi <- FALSE
87 sumi = FALSE
88 }else{ 81 }else{
89 sumi = TRUE 82 sumi <- TRUE
90
91 } 83 }
92 84
85 if (opt$av_level == "intra") {
86 pa <- msPurity::averageIntraFragSpectra(pa,
87 minfrac = opt$minfrac,
88 minnum = opt$minnum,
89 ppm = opt$ppm,
90 snr = opt$snr,
91 ra = opt$ra,
92 av = opt$av,
93 sumi = sumi,
94 rmp = rmp,
95 cores = opt$cores)
93 96
94 if(opt$av_level=="intra"){ 97 } else if (opt$av_level == "inter") {
95
96 pa <- msPurity::averageIntraFragSpectra(pa,
97 minfrac=opt$minfrac,
98 minnum=opt$minnum,
99 ppm=opt$ppm,
100 snr=opt$snr,
101 ra=opt$ra,
102 av=opt$av,
103 sumi=sumi,
104 rmp=rmp,
105 cores=opt$cores)
106
107 } else if(opt$av_level=="inter"){
108 98
109 pa <- msPurity::averageInterFragSpectra(pa, 99 pa <- msPurity::averageInterFragSpectra(pa,
110 minfrac=opt$minfrac, 100 minfrac = opt$minfrac,
111 minnum=opt$minnum, 101 minnum = opt$minnum,
112 ppm=opt$ppm, 102 ppm = opt$ppm,
113 snr=opt$snr, 103 snr = opt$snr,
114 ra=opt$ra, 104 ra = opt$ra,
115 av=opt$av, 105 av = opt$av,
116 sumi=sumi, 106 sumi = sumi,
117 rmp=rmp, 107 rmp = rmp,
118 cores=opt$cores) 108 cores = opt$cores)
119 } else if(opt$av_level=="all"){ 109 } else if (opt$av_level == "all") {
120 110
121 pa <- msPurity::averageAllFragSpectra(pa, 111 pa <- msPurity::averageAllFragSpectra(pa,
122 minfrac=opt$minfrac, 112 minfrac = opt$minfrac,
123 minnum=opt$minnum, 113 minnum = opt$minnum,
124 ppm=opt$ppm, 114 ppm = opt$ppm,
125 snr=opt$snr, 115 snr = opt$snr,
126 ra=opt$ra, 116 ra = opt$ra,
127 av=opt$av, 117 av = opt$av,
128 sumi=sumi, 118 sumi = sumi,
129 rmp=rmp, 119 rmp = rmp,
130 cores=opt$cores) 120 cores = opt$cores)
131
132 } 121 }
133 122
134 print(pa) 123 print(pa)
135 save(pa, file=opt$out_rdata) 124 save(pa, file = opt$out_rdata)
136 125
137 126 if (length(pa) > 0) {
138 if (length(pa)>0){
139 127
140 av_spectra <- plyr::ldply(pa@av_spectra, get_av_spectra) 128 av_spectra <- plyr::ldply(pa@av_spectra, get_av_spectra)
141 129
142 if (nrow(av_spectra)==0){ 130 if (nrow(av_spectra) == 0) {
143 message('No average spectra available') 131 message("No average spectra available")
144 } else{ 132 } else {
145 colnames(av_spectra)[1] <- 'grpid' 133 colnames(av_spectra)[1] <- "grpid"
146 av_spectra$grpid <- names(pa@av_spectra)[av_spectra$grpid] 134 av_spectra$grpid <- names(pa@av_spectra)[av_spectra$grpid]
147 135
148 if((length(pa@av_intra_params)>0) || (length(pa@av_inter_params)>0) ){ 136 if ((length(pa@av_intra_params) > 0) || (length(pa@av_inter_params) > 0)) {
149 # Add some extra info (only required if av_intra or av_inter performed) 137 # Add some extra info (only required if av_intra or av_inter performed)
150 colnames(av_spectra)[2] <- 'fileid' 138 colnames(av_spectra)[2] <- "fileid"
151 av_spectra$avid <- 1:nrow(av_spectra) 139 av_spectra$avid <- seq_len(nrow(av_spectra))
152 140
153 filenames <- sapply(av_spectra$fileid, function(x) names(pa@fileList)[as.integer(x)]) 141 filenames <- sapply(av_spectra$fileid,
154 # filenames_galaxy <- sapply(av_spectra$fileid, function(x) basename(pa@fileList[as.integer(x)])) 142 function(x) names(pa@fileList)[as.integer(x)])
155 143 # filenames_galaxy <- sapply(
156 av_spectra = as.data.frame(append(av_spectra, list(filename = filenames), after=2)) 144 # av_spectra$fileid, function(x) basename(pa@fileList[as.integer(x)]))
145
146 av_spectra <- as.data.frame(
147 append(av_spectra, list(filename = filenames), after = 2))
157 } 148 }
158 149
159 150
160 print(head(av_spectra)) 151 print(head(av_spectra))
161 write.table(av_spectra, opt$out_peaklist, row.names=FALSE, sep='\t') 152 write.table(av_spectra, opt$out_peaklist, row.names = FALSE, sep = "\t")
162 153
163 } 154 }
164 } 155 }
165