comparison averageFragSpectra.R @ 8:e98ae5c1f4e7 draft

planemo upload for repository https://github.com/computational-metabolomics/mspurity-galaxy commit 7e1748612a9f9dce11a9e54ff36752b600e7aea3
author computational-metabolomics
date Wed, 12 Jun 2024 16:01:28 +0000
parents 52d48bcd3608
children
comparison
equal deleted inserted replaced
7:94eb1fcfde88 8:e98ae5c1f4e7
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
9 if (length(x$av_intra) > 0) { 8 if (length(x$av_intra) > 0) {
10 av_intra_df <- plyr::ldply(x$av_intra) 9 av_intra_df <- plyr::ldply(x$av_intra)
11 10
12 if (nrow(av_intra_df) == 0) { 11 if (nrow(av_intra_df) == 0) {
13 av_intra_df <- NULL 12 av_intra_df <- NULL
14 }else{ 13 } else {
15 av_intra_df$method <- "intra" 14 av_intra_df$method <- "intra"
16 } 15 }
17 16 } else {
18 }else{
19 av_intra_df <- NULL 17 av_intra_df <- NULL
20 } 18 }
21 19
22 if ((is.null(x$av_inter)) || (nrow(x$av_inter) == 0)) { 20 if ((is.null(x$av_inter)) || (nrow(x$av_inter) == 0)) {
23 av_inter_df <- NULL 21 av_inter_df <- NULL
24 }else{ 22 } else {
25 av_inter_df <- x$av_inter 23 av_inter_df <- x$av_inter
26 av_inter_df$method <- "inter" 24 av_inter_df$method <- "inter"
27 } 25 }
28 26
29 if ((is.null(x$av_all)) || (nrow(x$av_all) == 0)) { 27 if ((is.null(x$av_all)) || (nrow(x$av_all) == 0)) {
30 av_all_df <- NULL 28 av_all_df <- NULL
31 }else{ 29 } else {
32 av_all_df <- x$av_all 30 av_all_df <- x$av_all
33 av_all_df$method <- "all" 31 av_all_df$method <- "all"
34 } 32 }
35 33
36 combined <- plyr::rbind.fill(av_intra_df, av_inter_df, av_all_df) 34 combined <- plyr::rbind.fill(av_intra_df, av_inter_df, av_all_df)
58 opt <- parse_args(OptionParser(option_list = option_list)) 56 opt <- parse_args(OptionParser(option_list = option_list))
59 print(opt) 57 print(opt)
60 58
61 59
62 load_r_data <- function(rdata_path, name) { 60 load_r_data <- function(rdata_path, name) {
63 #loads an RData file, and returns the named xset object if it is there 61 # loads an RData file, and returns the named xset object if it is there
64 load(rdata_path) 62 load(rdata_path)
65 return(get(ls()[ls() %in% name])) 63 return(get(ls()[ls() %in% name]))
66 } 64 }
67 65
68 # Requires 66 # Requires
69 pa <- load_r_data(opt$pa, "pa") 67 pa <- load_r_data(opt$pa, "pa")
70 68
71 pa@cores <- opt$cores 69 pa@cores <- opt$cores
72 70
73 if (is.null(opt$rmp)) { 71 if (is.null(opt$rmp)) {
74 rmp <- FALSE 72 rmp <- FALSE
75 }else{ 73 } else {
76 rmp <- TRUE 74 rmp <- TRUE
77 } 75 }
78 76
79 if (is.null(opt$sumi)) { 77 if (is.null(opt$sumi)) {
80 sumi <- FALSE 78 sumi <- FALSE
81 }else{ 79 } else {
82 sumi <- TRUE 80 sumi <- TRUE
83 } 81 }
84 82
85 if (opt$av_level == "intra") { 83 if (opt$av_level == "intra") {
86 pa <- msPurity::averageIntraFragSpectra(pa, 84 pa <- msPurity::averageIntraFragSpectra(pa,
87 minfrac = opt$minfrac, 85 minfrac = opt$minfrac,
88 minnum = opt$minnum, 86 minnum = opt$minnum,
89 ppm = opt$ppm, 87 ppm = opt$ppm,
90 snr = opt$snr, 88 snr = opt$snr,
91 ra = opt$ra, 89 ra = opt$ra,
92 av = opt$av, 90 av = opt$av,
93 sumi = sumi, 91 sumi = sumi,
94 rmp = rmp, 92 rmp = rmp,
95 cores = opt$cores) 93 cores = opt$cores
96 94 )
97 } else if (opt$av_level == "inter") { 95 } else if (opt$av_level == "inter") {
98
99 pa <- msPurity::averageInterFragSpectra(pa, 96 pa <- msPurity::averageInterFragSpectra(pa,
100 minfrac = opt$minfrac, 97 minfrac = opt$minfrac,
101 minnum = opt$minnum, 98 minnum = opt$minnum,
102 ppm = opt$ppm, 99 ppm = opt$ppm,
103 snr = opt$snr, 100 snr = opt$snr,
104 ra = opt$ra, 101 ra = opt$ra,
105 av = opt$av, 102 av = opt$av,
106 sumi = sumi, 103 sumi = sumi,
107 rmp = rmp, 104 rmp = rmp,
108 cores = opt$cores) 105 cores = opt$cores
106 )
109 } else if (opt$av_level == "all") { 107 } else if (opt$av_level == "all") {
110
111 pa <- msPurity::averageAllFragSpectra(pa, 108 pa <- msPurity::averageAllFragSpectra(pa,
112 minfrac = opt$minfrac, 109 minfrac = opt$minfrac,
113 minnum = opt$minnum, 110 minnum = opt$minnum,
114 ppm = opt$ppm, 111 ppm = opt$ppm,
115 snr = opt$snr, 112 snr = opt$snr,
116 ra = opt$ra, 113 ra = opt$ra,
117 av = opt$av, 114 av = opt$av,
118 sumi = sumi, 115 sumi = sumi,
119 rmp = rmp, 116 rmp = rmp,
120 cores = opt$cores) 117 cores = opt$cores
118 )
121 } 119 }
122 120
123 print(pa) 121 print(pa)
124 save(pa, file = opt$out_rdata) 122 save(pa, file = opt$out_rdata)
125 123
126 if (length(pa) > 0) { 124 if (length(pa) > 0) {
127
128 av_spectra <- plyr::ldply(pa@av_spectra, get_av_spectra) 125 av_spectra <- plyr::ldply(pa@av_spectra, get_av_spectra)
129 126
130 if (nrow(av_spectra) == 0) { 127 if (nrow(av_spectra) == 0) {
131 message("No average spectra available") 128 message("No average spectra available")
132 } else { 129 } else {
136 if ((length(pa@av_intra_params) > 0) || (length(pa@av_inter_params) > 0)) { 133 if ((length(pa@av_intra_params) > 0) || (length(pa@av_inter_params) > 0)) {
137 # Add some extra info (only required if av_intra or av_inter performed) 134 # Add some extra info (only required if av_intra or av_inter performed)
138 colnames(av_spectra)[2] <- "fileid" 135 colnames(av_spectra)[2] <- "fileid"
139 av_spectra$avid <- seq_len(nrow(av_spectra)) 136 av_spectra$avid <- seq_len(nrow(av_spectra))
140 137
141 filenames <- sapply(av_spectra$fileid, 138 filenames <- sapply(
142 function(x) names(pa@fileList)[as.integer(x)]) 139 av_spectra$fileid,
140 function(x) names(pa@fileList)[as.integer(x)]
141 )
143 # filenames_galaxy <- sapply( 142 # filenames_galaxy <- sapply(
144 # av_spectra$fileid, function(x) basename(pa@fileList[as.integer(x)])) 143 # av_spectra$fileid, function(x) basename(pa@fileList[as.integer(x)]))
145 144
146 av_spectra <- as.data.frame( 145 av_spectra <- as.data.frame(
147 append(av_spectra, list(filename = filenames), after = 2)) 146 append(av_spectra, list(filename = filenames), after = 2)
147 )
148 } 148 }
149 149
150 150
151 print(head(av_spectra)) 151 print(head(av_spectra))
152 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")
153
154 } 153 }
155 } 154 }