Mercurial > repos > galaxyp > psm_to_sam
comparison PSM2SAM.R @ 0:90ecb65017a0 draft
planemo upload for repository https://github.com/galaxyproteomics/tools-galaxyp/tools/bumbershoot/psm2sam commit 9910ed076e4b8a3f083351b89fa861d0a4a93beb
author | galaxyp |
---|---|
date | Wed, 17 May 2017 20:23:27 -0400 |
parents | |
children | 757ecf27b4a9 |
comparison
equal
deleted
inserted
replaced
-1:000000000000 | 0:90ecb65017a0 |
---|---|
1 #!/usr/bin/env Rscript | |
2 | |
3 ## begin warning handler | |
4 withCallingHandlers({ | |
5 | |
6 library(methods) # Because Rscript does not always do this | |
7 | |
8 options('useFancyQuotes' = FALSE) | |
9 | |
10 suppressPackageStartupMessages(library("optparse")) | |
11 suppressPackageStartupMessages(library("RGalaxy")) | |
12 | |
13 | |
14 option_list <- list() | |
15 | |
16 option_list$exon_anno <- make_option('--exon_anno', type='character') | |
17 option_list$proteinseq <- make_option('--proteinseq', type='character') | |
18 option_list$procodingseq <- make_option('--procodingseq', type='character') | |
19 option_list$bam_file <- make_option('--bam', type='character') | |
20 option_list$idpDB_file <- make_option('--idpDB', type='character') | |
21 option_list$pepXmlTab_file <- make_option('--pepXmlTab', type='character') | |
22 option_list$peptideShakerPsmReport_file <- make_option('--peptideShakerPsmReport', type='character') | |
23 option_list$variantAnnotation_file <- make_option('--variantAnnotation', type='character') | |
24 option_list$searchEngineScore <- make_option('--searchEngineScore', type='character') | |
25 | |
26 | |
27 opt <- parse_args(OptionParser(option_list=option_list)) | |
28 | |
29 | |
30 psm2sam <- function( | |
31 exon_anno_file = GalaxyInputFile(required=TRUE), | |
32 proteinseq_file = GalaxyInputFile(required=TRUE), | |
33 procodingseq_file = GalaxyInputFile(required=TRUE), | |
34 bam_file = GalaxyInputFile(required=TRUE), | |
35 idpDB_file = GalaxyInputFile(required=FALSE), | |
36 pepXmlTab_file = GalaxyInputFile(required=FALSE), | |
37 peptideShakerPsmReport_file = GalaxyInputFile(required=FALSE), | |
38 variantAnnotation_file = GalaxyInputFile(required=FALSE), | |
39 searchEngineScore = GalaxyCharacterParam(required=FALSE) | |
40 ) | |
41 { | |
42 options(stringsAsFactors = FALSE) | |
43 | |
44 if (length(bam_file) == 0) | |
45 { | |
46 stop("BAM file must be specified to provide sequence headers") | |
47 } | |
48 | |
49 outputHeader = grep("^@(?!PG)", readLines(bam_file, n=500, warn=FALSE), value=TRUE, perl=TRUE) | |
50 if (length(outputHeader) == 0) | |
51 { | |
52 stop("failed to read header lines from bam_file") | |
53 } | |
54 | |
55 # load customProDB from GitHub (NOTE: downloading the zip is faster than cloning the repo with git2r or devtools::install_github) | |
56 download.file("https://github.com/chambm/customProDB/archive/master.zip", "customProDB.zip", quiet=TRUE) | |
57 unzip("customProDB.zip") | |
58 devtools::load_all("customProDB-master") | |
59 | |
60 # load proBAMr from GitHub | |
61 download.file("https://github.com/chambm/proBAMr/archive/master.zip", "proBAMr.zip", quiet=TRUE) | |
62 unzip("proBAMr.zip") | |
63 devtools::load_all("proBAMr-master") | |
64 | |
65 psmInputLength = length(idpDB_file)+length(pepXmlTab_file)+length(peptideShakerPsmReport_file) | |
66 if (psmInputLength == 0) | |
67 { | |
68 stop("one of the input PSM file parameters must be specified") | |
69 } | |
70 else if (psmInputLength > 1) | |
71 { | |
72 stop("only one of the input PSM file parameters can be specified") | |
73 } | |
74 | |
75 if (length(idpDB_file) > 0) | |
76 { | |
77 if (length(searchEngineScore) == 0) | |
78 stop("searchEngineScore parameter must be specified when reading IDPicker PSMs, e.g. 'MyriMatch:MVH'") | |
79 passedPSM = readIdpDB(idpDB_file, searchEngineScore) | |
80 } | |
81 else if (length(pepXmlTab_file) > 0) | |
82 { | |
83 if (length(searchEngineScore) == 0) | |
84 stop("searchEngineScore parameter must be specified when reading pepXmlTab PSMs, e.g. 'mvh'") | |
85 passedPSM = readPepXmlTab(pepXmlTab_file, searchEngineScore) | |
86 } | |
87 else if (length(peptideShakerPsmReport_file) > 0) | |
88 { | |
89 if (length(searchEngineScore) > 0) | |
90 warning("searchEngineScore parameter is ignored when reading PeptideShaker PSM report") | |
91 passedPSM = readPeptideShakerPsmReport(peptideShakerPsmReport_file) | |
92 } | |
93 | |
94 load(exon_anno_file) | |
95 load(proteinseq_file) | |
96 load(procodingseq_file) | |
97 | |
98 if (length(variantAnnotation_file) > 0) | |
99 { | |
100 load(variantAnnotation_file) # variantAnnotation list, with members snvprocoding/snvproseq and indelprocoding/indelproseq | |
101 | |
102 varprocoding = unique(rbind(variantAnnotation$snvprocoding, variantAnnotation$indelprocoding)) | |
103 varproseq = unique(rbind(variantAnnotation$snvproseq, variantAnnotation$indelproseq)) | |
104 } | |
105 else | |
106 { | |
107 varprocoding = NULL | |
108 varproseq = NULL | |
109 } | |
110 | |
111 # add proBAMr program key | |
112 outputHeader = c(outputHeader, paste0("@PG\tID:proBAMr\tVN:", packageVersion("proBAMr"))) | |
113 | |
114 # first write header lines to the output SAM | |
115 writeLines(outputHeader, "output.sam") | |
116 | |
117 # then write the PSM "reads" | |
118 PSMtab2SAM(passedPSM, exon, | |
119 proteinseq, procodingseq, | |
120 varproseq, varprocoding, | |
121 outfile = "output.sam", | |
122 show_progress = FALSE) | |
123 | |
124 invisible(NULL) | |
125 } | |
126 | |
127 params <- list() | |
128 for(param in names(opt)) | |
129 { | |
130 if (!param == "help") | |
131 params[param] <- opt[param] | |
132 } | |
133 | |
134 setClass("GalaxyRemoteError", contains="character") | |
135 wrappedFunction <- function(f) | |
136 { | |
137 tryCatch(do.call(f, params), | |
138 error=function(e) new("GalaxyRemoteError", conditionMessage(e))) | |
139 } | |
140 | |
141 | |
142 suppressPackageStartupMessages(library(RGalaxy)) | |
143 do.call(psm2sam, params) | |
144 | |
145 ## end warning handler | |
146 }, warning = function(w) { | |
147 cat(paste("Warning:", conditionMessage(w), "\n")) | |
148 invokeRestart("muffleWarning") | |
149 }) |