|
0
|
1 ##!/usr/bin/Rscript
|
|
|
2
|
|
|
3 ## Setup R error handling to go to stderr
|
|
|
4 options( show.error.messages = FALSE, error = function () { cat( geterrmessage(), file=stderr() ); q( "no", 1, F ) } )
|
|
|
5
|
|
|
6 # we need that to not crash galaxy with an UTF8 error on German LC settings.
|
|
|
7 Sys.setlocale("LC_MESSAGES", "en_US.UTF-8")
|
|
|
8
|
|
|
9 suppressMessages(library('getopt'))
|
|
|
10
|
|
|
11 #get options, using the spec as defined by the enclosed list.
|
|
|
12 #we read the options from the default: commandArgs(TRUE).
|
|
|
13 spec = matrix(c(
|
|
|
14 'treated', 't', 1, "character",
|
|
|
15 'control', 'c', 2, "character",
|
|
|
16 'method', 'm', 1, "character",
|
|
|
17 'fuComplexity', 'b', 2, "integer",
|
|
|
18 'k2nTreated', 'k', 2, "character",
|
|
|
19 'k2nControl', 'a', 2, "character",
|
|
|
20 'cutoff', 'g', 1, "integer",
|
|
|
21 'reference', 'r', 2, "character",
|
|
|
22 'compdata', 'h', 2, "character",
|
|
|
23 'dtcr', 'd', 0, "logical",
|
|
|
24 'dtcrWindow', 'e', 2, "integer",
|
|
|
25 'dtcrToZero', 'f', 2, "character",
|
|
|
26 'slograt', 's', 0, "logical",
|
|
|
27 'slogratWindow', 'p', 2, "integer",
|
|
|
28 'depthCorrection', 'q', 2, "character",
|
|
|
29 'pseudocount', 'l', 2, "integer",
|
|
|
30 'swinsor', 'w', 0, "logical",
|
|
|
31 'swinsorWindow', 'x', 2, "integer",
|
|
|
32 'winsorLevel', 'y', 2, "double",
|
|
|
33 'fixQuantile', 'z', 2, "character",
|
|
|
34 'ntOffset', 'n', 2, "integer",
|
|
|
35 'outputDir', 'o', 1, "character",
|
|
|
36 'bedgraph', 'v', 2, "character",
|
|
|
37 'bed', 'j', 2, "character",
|
|
|
38 'genome', 'i', 2, "character",
|
|
|
39 'trackName', 'u', 2, "character"
|
|
|
40 ), byrow=TRUE, ncol=4);
|
|
|
41
|
|
|
42 opt = getopt(spec);
|
|
|
43
|
|
|
44 suppressMessages(require(RNAprobR))
|
|
|
45
|
|
|
46
|
|
|
47 #Create output dir
|
|
|
48 dir.create(opt$outputDir, showWarnings = FALSE, recursive = TRUE)
|
|
|
49
|
|
|
50
|
|
|
51 if (opt$method=="counts"){
|
|
|
52 treated_euc <- readsamples(opt$treated, euc=opt$method)
|
|
|
53 } else if (opt$method=="Fu") {
|
|
|
54 treated_euc <- readsamples(opt$treated, euc=opt$method, m = opt$fuComplexity)
|
|
|
55 } else {
|
|
|
56 treated_euc <- readsamples(opt$treated, euc=opt$method, k2n_files=opt$k2nTreated)
|
|
|
57 }
|
|
|
58
|
|
|
59
|
|
|
60 if ( !is.null(opt$reference)){
|
|
|
61 comp_treated <- comp(treated_euc, cutoff=opt$cutoff,
|
|
|
62 fasta_file=opt$reference)
|
|
|
63 } else {
|
|
|
64 comp_treated <- comp(treated_euc, cutoff=opt$cutoff)
|
|
|
65 }
|
|
|
66
|
|
|
67
|
|
|
68 #If present, read control file
|
|
|
69 if ( !is.null(opt$control) && opt$control != 'None'){
|
|
|
70 if (opt$method=="counts"){
|
|
|
71 control_euc <- readsamples(opt$control, euc=opt$method)
|
|
|
72 } else if (opt$method=="Fu") {
|
|
|
73 control_euc <- readsamples(opt$control, euc=opt$method, m = opt$fuComplexity)
|
|
|
74 } else {
|
|
|
75 control_euc <- readsamples(opt$control, euc=opt$method, k2n_files=opt$k2nControl)
|
|
|
76 }
|
|
|
77 if ( !is.null(opt$reference)){
|
|
|
78 comp_control <- comp(control_euc, cutoff=opt$cutoff,
|
|
|
79 fasta_file=opt$reference)
|
|
|
80 } else {
|
|
|
81 comp_control <- comp(control_euc, cutoff=opt$cutoff)
|
|
|
82 }
|
|
|
83 }
|
|
|
84
|
|
|
85 #compdata
|
|
|
86 if ( !is.null(opt$compdata) ) {
|
|
|
87
|
|
|
88 normalized <- compdata(comp_treated, nt_offset=opt$ntOffset)
|
|
|
89 colnames(mcols(normalized))[1:4] <- paste(colnames(mcols(normalized))[1:4], ".treated", sep="")
|
|
|
90
|
|
|
91 if ( !is.null(opt$control) && opt$control != 'None') {
|
|
|
92 normalized <- compdata(comp_control, nt_offset=opt$ntOffset, add_to=normalized)
|
|
|
93 colnames(mcols(normalized))[6:9] <- c("TC.control", "TCR.control", "Cover.control", "PC.control")
|
|
|
94 }
|
|
|
95 }
|
|
|
96
|
|
|
97 #dtcr
|
|
|
98 if ( !is.null(opt$dtcr) ) {
|
|
|
99 if ( exists("normalized")) {
|
|
|
100 normalized <- dtcr(comp_control, comp_treated, window_size=opt$dtcrWindow,
|
|
|
101 nt_offset=opt$ntOffset, bring_to_zero=opt$dtcrToZero,
|
|
|
102 add_to=normalized)
|
|
|
103 }
|
|
|
104 else {
|
|
|
105 normalized <- dtcr(comp_control, comp_treated, window_size=opt$dtcrWindow,
|
|
|
106 nt_offset=opt$ntOffset, bring_to_zero=opt$dtcrToZero)
|
|
|
107 }
|
|
|
108 }
|
|
|
109
|
|
|
110 #slograt
|
|
|
111 if ( !is.null(opt$slograt) ) {
|
|
|
112 if ( exists("normalized")) {
|
|
|
113 normalized <- slograt(comp_control, comp_treated, window_size=opt$slogratWindow,
|
|
|
114 nt_offset=opt$ntOffset, depth_correction = opt$depthCorrection,
|
|
|
115 pseudocount=opt$pseudocount, add_to=normalized)
|
|
|
116 }
|
|
|
117 else {
|
|
|
118 normalized <- slograt(comp_control, comp_treated, window_size=opt$slogratWindow,
|
|
|
119 nt_offset=opt$ntOffset, depth_correction = opt$depthCorrection,
|
|
|
120 pseudocount=opt$pseudocount)
|
|
|
121 }
|
|
|
122 }
|
|
|
123
|
|
|
124 #swinsor
|
|
|
125 if ( !is.null(opt$swinsor) ) {
|
|
|
126 if ( exists("normalized")) {
|
|
|
127 normalized <- swinsor(comp_treated, winsor_level = opt$winsorLevel,
|
|
|
128 window_size=opt$swinsorWindow,
|
|
|
129 only_top=opt$fixQuantile, nt_offset=opt$ntOffset,
|
|
|
130 add_to=normalized)
|
|
|
131 }
|
|
|
132 else {
|
|
|
133 normalized <- swinsor(comp_treated, winsor_level = opt$winsorLevel,
|
|
|
134 window_size=opt$swinsorWindow,
|
|
|
135 only_top=opt$fixQuantile, nt_offset=opt$ntOffset)
|
|
|
136 }
|
|
|
137 }
|
|
|
138
|
|
|
139 #bedgraph output
|
|
|
140 if ( !is.null(opt$bedgraph)) {
|
|
|
141 if ( !is.null(opt$dtcr)){
|
|
|
142 norm2bedgraph(normalized, bed_file = opt$bed, norm_method = "dtcr",
|
|
|
143 genome_build = opt$genome,
|
|
|
144 bedgraph_out_file = paste(opt$outputDir, "/dtcr",
|
|
|
145 sep=""),
|
|
|
146 track_name = opt$trackName,
|
|
|
147 track_description = opt$trackDesc)
|
|
|
148 }
|
|
|
149 if ( !is.null(opt$slograt)){
|
|
|
150 norm2bedgraph(normalized, bed_file = opt$bed, norm_method = "slograt",
|
|
|
151 genome_build = opt$genome,
|
|
|
152 bedgraph_out_file = paste(opt$outputDir, "/slograt",
|
|
|
153 sep=""),
|
|
|
154 track_name = opt$trackName,
|
|
|
155 track_description = opt$trackDesc)
|
|
|
156 }
|
|
|
157 if ( !is.null(opt$swinsor)){
|
|
|
158 norm2bedgraph(normalized, bed_file = opt$bed, norm_method = "swinsor",
|
|
|
159 genome_build = opt$genome,
|
|
|
160 bedgraph_out_file = paste(opt$outputDir, "/swinsor",
|
|
|
161 sep=""),
|
|
|
162 track_name = opt$trackName,
|
|
|
163 track_description = opt$trackDesc)
|
|
|
164 }
|
|
|
165 }
|
|
|
166
|
|
|
167 output <- GR2norm_df(normalized)
|
|
|
168 write.table( output, paste(opt$outputDir, "/norm_df.txt", sep = ""), sep = "\t",
|
|
|
169 quote = F, row.names = F)
|
