comparison small_rna_maps.r @ 12:d33263e6e812 draft

planemo upload for repository https://github.com/ARTbio/tools-artbio/tree/master/tools/small_rna_maps commit b0676fd329c2ca50002f9f2fede531d8e550569f
author artbio
date Sat, 07 Apr 2018 06:14:50 -0400
parents a561a71bd7d7
children cd75c72e1d75
comparison
equal deleted inserted replaced
11:a561a71bd7d7 12:d33263e6e812
1 ## Setup R error handling to go to stderr 1 ## Setup R error handling to go to stderr
2 options( show.error.messages=F, 2 options( show.error.messages=F,
3 error = function () { cat( geterrmessage(), file=stderr() ); q( "no", 1, F ) } ) 3 error = function () { cat( geterrmessage(), file=stderr() ); q( "no", 1, F ) } )
4 options(warn = -1) 4 options(warn = -1)
5 library(RColorBrewer) 5 library(RColorBrewer)
6 library(lattice) 6 library(lattice)
7 library(latticeExtra) 7 library(latticeExtra)
8 library(grid) 8 library(grid)
9 library(gridExtra) 9 library(gridExtra)
10 library(optparse) 10 library(optparse)
11 11
12 12
13 option_list <- list( 13 option_list <- list(
14 make_option(c("-f", "--first_dataframe"), type="character", help="path to first dataframe"), 14 make_option(c("-i", "--ymin"), type="double", help="set min ylimit. e.g. '-100.0'"),
15 make_option(c("-e", "--extra_dataframe"), type="character", help="path to additional dataframe"), 15 make_option(c("-a", "--ymax"), type="double", help="set max ylimit. e.g. '100.0'"),
16 make_option(c("-n", "--normalization"), type="character", help="space-separated normalization/size factors"), 16 make_option(c("-f", "--first_dataframe"), type="character", help="path to first dataframe"),
17 make_option("--first_plot_method", type = "character", help="How additional data should be plotted"), 17 make_option(c("-e", "--extra_dataframe"), type="character", help="path to additional dataframe"),
18 make_option("--extra_plot_method", type = "character", help="How additional data should be plotted"), 18 make_option(c("-n", "--normalization"), type="character", help="space-separated normalization/size factors"),
19 make_option("--global", type = "character", help="data should be plotted as global size distribution"), 19 make_option("--first_plot_method", type = "character", help="How additional data should be plotted"),
20 make_option("--output_pdf", type = "character", help="path to the pdf file with plots") 20 make_option("--extra_plot_method", type = "character", help="How additional data should be plotted"),
21 ) 21 make_option("--global", type = "character", help="data should be plotted as global size distribution"),
22 22 make_option("--output_pdf", type = "character", help="path to the pdf file with plots")
23 )
24
23 parser <- OptionParser(usage = "%prog [options] file", option_list = option_list) 25 parser <- OptionParser(usage = "%prog [options] file", option_list = option_list)
24 args = parse_args(parser) 26 args = parse_args(parser)
25 27
26 # data frames implementation 28 # data frames implementation
27 ## first table 29 ## first table
28 Table = read.delim(args$first_dataframe, header=T, row.names=NULL) 30 Table = read.delim(args$first_dataframe, header=T, row.names=NULL)
29 if (args$first_plot_method == "Counts" | args$first_plot_method == "Size") { 31 if (args$first_plot_method == "Counts" | args$first_plot_method == "Size") {
30 Table <- within(Table, Counts[Polarity=="R"] <- (Counts[Polarity=="R"]*-1)) 32 Table <- within(Table, Counts[Polarity=="R"] <- (Counts[Polarity=="R"]*-1))
31 } 33 }
32 n_samples=length(unique(Table$Dataset)) 34 n_samples=length(unique(Table$Dataset))
33 samples = unique(Table$Dataset) 35 samples = unique(Table$Dataset)
34 if (args$normalization != "") { 36 if (args$normalization != "") {
35 norm_factors = as.numeric(unlist(strsplit(args$normalization, " "))) 37 norm_factors = as.numeric(unlist(strsplit(args$normalization, " ")))
36 } else { 38 } else {
37 norm_factors = rep(1, n_samples) 39 norm_factors = rep(1, n_samples)
38 } 40 }
39 if (args$first_plot_method == "Counts" | args$first_plot_method == "Size" | args$first_plot_method == "Coverage") { 41 if (args$first_plot_method == "Counts" | args$first_plot_method == "Size" | args$first_plot_method == "Coverage") {
40 i = 1 42 i = 1
41 for (sample in samples) { 43 for (sample in samples) {
42 Table[, length(Table)][Table$Dataset==sample] <- Table[, length(Table)][Table$Dataset==sample]*norm_factors[i] 44 Table[, length(Table)][Table$Dataset==sample] <- Table[, length(Table)][Table$Dataset==sample]*norm_factors[i]
43 i = i + 1 45 i = i + 1
44 } 46 }
45 } 47 }
46 genes=unique(Table$Chromosome) 48 genes=unique(Table$Chromosome)
47 per_gene_readmap=lapply(genes, function(x) subset(Table, Chromosome==x)) 49 per_gene_readmap=lapply(genes, function(x) subset(Table, Chromosome==x))
48 per_gene_limit=lapply(genes, function(x) c(1, unique(subset(Table, Chromosome==x)$Chrom_length)) ) 50 per_gene_limit=lapply(genes, function(x) c(1, unique(subset(Table, Chromosome==x)$Chrom_length)) )
49 n_genes=length(per_gene_readmap) 51 n_genes=length(per_gene_readmap)
50 # second table 52 # second table
51 if (args$extra_plot_method != '') { 53 if (args$extra_plot_method != '') {
52 ExtraTable=read.delim(args$extra_dataframe, header=T, row.names=NULL) 54 ExtraTable=read.delim(args$extra_dataframe, header=T, row.names=NULL)
53 if (args$extra_plot_method == "Counts" | args$extra_plot_method=='Size') { 55 if (args$extra_plot_method == "Counts" | args$extra_plot_method=='Size') {
54 ExtraTable <- within(ExtraTable, Counts[Polarity=="R"] <- (Counts[Polarity=="R"]*-1)) 56 ExtraTable <- within(ExtraTable, Counts[Polarity=="R"] <- (Counts[Polarity=="R"]*-1))
57 }
58 if (args$extra_plot_method == "Counts" | args$extra_plot_method == "Size" | args$extra_plot_method == "Coverage") {
59 i = 1
60 for (sample in samples) {
61 ExtraTable[, length(ExtraTable)][ExtraTable$Dataset==sample] <- ExtraTable[, length(ExtraTable)][ExtraTable$Dataset==sample]*norm_factors[i]
62 i = i + 1
55 } 63 }
56 if (args$extra_plot_method == "Counts" | args$extra_plot_method == "Size" | args$extra_plot_method == "Coverage") { 64 }
57 i = 1 65 per_gene_size=lapply(genes, function(x) subset(ExtraTable, Chromosome==x))
58 for (sample in samples) {
59 ExtraTable[, length(ExtraTable)][ExtraTable$Dataset==sample] <- ExtraTable[, length(ExtraTable)][ExtraTable$Dataset==sample]*norm_factors[i]
60 i = i + 1
61 }
62 }
63 per_gene_size=lapply(genes, function(x) subset(ExtraTable, Chromosome==x))
64 } 66 }
65 67
66 ## functions 68 ## functions
67 globalbc = function(df, global="", ...) { 69 globalbc = function(df, global="", ...) {
68 if (global == "yes") { 70 if (global == "yes") {
69 bc <- barchart(Counts~as.factor(Size)|factor(Dataset, levels=unique(Dataset)), 71 bc <- barchart(Counts~as.factor(Size)|factor(Dataset, levels=unique(Dataset)),
70 data = df, origin = 0, 72 data = df, origin = 0,
71 horizontal=FALSE, 73 horizontal=FALSE,
72 col=c("darkblue"), 74 col=c("darkblue"),
73 scales=list(y=list(tick.number=4, rot=90, relation="free", cex=0.5, alternating=T), x=list(rot=0, cex=0.6, tck=0.5, alternating=c(3,3))), 75 scales=list(y=list(tick.number=4, rot=90, relation="free", cex=0.5, alternating=T), x=list(rot=0, cex=0.6, tck=0.5, alternating=c(3,3))),
74 xlab=list(label=bottom_first_method[[args$first_plot_method]], cex=.85), 76 xlab=list(label=bottom_first_method[[args$first_plot_method]], cex=.85),
75 ylab=list(label=legend_first_method[[args$first_plot_method]], cex=.85), 77 ylab=list(label=legend_first_method[[args$first_plot_method]], cex=.85),
76 main=title_first_method[[args$first_plot_method]], 78 main=title_first_method[[args$first_plot_method]],
77 layout = c(2, 6), newpage=T, 79 layout = c(2, 6), newpage=T,
78 as.table=TRUE, 80 as.table=TRUE,
79 aspect=0.5, 81 aspect=0.5,
80 strip = strip.custom(par.strip.text = list(cex = 1), which.given=1, bg="lightblue"), 82 strip = strip.custom(par.strip.text = list(cex = 1), which.given=1, bg="lightblue"),
81 ... 83 ...
82 ) 84 )
83 } else { 85 } else {
84 bc <- barchart(Counts~as.factor(Size)|factor(Dataset, levels=unique(Dataset)), 86 bc <- barchart(Counts~as.factor(Size)|factor(Dataset, levels=unique(Dataset)),
85 data = df, origin = 0, 87 data = df, origin = 0,
86 horizontal=FALSE, 88 horizontal=FALSE,
87 group=Polarity, 89 group=Polarity,
88 stack=TRUE, 90 stack=TRUE,
89 col=c('red', 'blue'), 91 col=c('red', 'blue'),
90 scales=list(y=list(tick.number=4, rot=90, relation="free", cex=0.5, alternating=T), x=list(rot=0, cex=0.6, tck=0.5, alternating=c(3,3))), 92 scales=list(y=list(tick.number=4, rot=90, relation="free", cex=0.5, alternating=T), x=list(rot=0, cex=0.6, tck=0.5, alternating=c(3,3))),
91 xlab=list(label=bottom_first_method[[args$first_plot_method]], cex=.85), 93 xlab=list(label=bottom_first_method[[args$first_plot_method]], cex=.85),
92 ylab=list(label=legend_first_method[[args$first_plot_method]], cex=.85), 94 ylab=list(label=legend_first_method[[args$first_plot_method]], cex=.85),
93 main=title_first_method[[args$first_plot_method]], 95 main=title_first_method[[args$first_plot_method]],
94 layout = c(2, 6), newpage=T, 96 layout = c(2, 6), newpage=T,
95 as.table=TRUE, 97 as.table=TRUE,
96 aspect=0.5, 98 aspect=0.5,
97 strip = strip.custom(par.strip.text = list(cex = 1), which.given=1, bg="lightblue"), 99 strip = strip.custom(par.strip.text = list(cex = 1), which.given=1, bg="lightblue"),
98 ... 100 ...
99 ) 101 )
100 } 102 }
101 return(bc) 103 return(bc)
102 } 104 }
103 plot_unit = function(df, method=args$first_plot_method, ...) { 105 plot_unit = function(df, method=args$first_plot_method, ...) {
104 if (method == 'Counts') { 106 if (exists('ymin', where=args)){
105 p = xyplot(Counts~Coordinate|factor(Dataset, levels=unique(Dataset))+factor(Chromosome, levels=unique(Chromosome)), 107 min=args$ymin
106 data=df, 108 }else{
107 type='h', 109 min=''
108 lwd=1.5, 110 }
109 scales= list(relation="free", x=list(rot=0, cex=0.7, axs="i", tck=0.5), y=list(tick.number=4, rot=90, cex=0.7)), 111 if ((exists('ymax', where=args))){
110 xlab=NULL, main=NULL, ylab=NULL, 112 max=args$ymax
111 as.table=T, 113 }else{
112 origin = 0, 114 max=''
113 horizontal=FALSE, 115 }
114 group=Polarity, 116 ylimits=c(min,max)
115 col=c("red","blue"), 117 if (method == 'Counts') {
116 par.strip.text = list(cex=0.7), 118 p = xyplot(Counts~Coordinate|factor(Dataset, levels=unique(Dataset))+factor(Chromosome, levels=unique(Chromosome)),
117 ...) 119 data=df,
118 } else if (method != "Size") { 120 type='h',
119 p = xyplot(eval(as.name(method))~Coordinate|factor(Dataset, levels=unique(Dataset))+factor(Chromosome, levels=unique(Chromosome)), 121 lwd=1.5,
120 data=df, 122 scales= list(relation="free", x=list(rot=0, cex=0.7, axs="i", tck=0.5), y=list(tick.number=4, rot=90, cex=0.7)),
121 type='p', 123 xlab=NULL, main=NULL, ylab=NULL, ylim=ylimits,
122 pch=19, 124 as.table=T,
123 cex=0.35, 125 origin = 0,
124 scales= list(relation="free", x=list(rot=0, cex=0.7, axs="i", tck=0.5), y=list(tick.number=4, rot=90, cex=0.7)), 126 horizontal=FALSE,
125 xlab=NULL, main=NULL, ylab=NULL, 127 group=Polarity,
126 as.table=T, 128 col=c("red","blue"),
127 origin = 0, 129 par.strip.text = list(cex=0.7),
128 horizontal=FALSE, 130 ...)
129 group=Polarity, 131 p=combineLimits(p)
130 col=c("red","blue"), 132 } else if (method != "Size") {
131 par.strip.text = list(cex=0.7), 133 p = xyplot(eval(as.name(method))~Coordinate|factor(Dataset, levels=unique(Dataset))+factor(Chromosome, levels=unique(Chromosome)),
132 ...) 134 data=df,
133 } else { 135 type='p',
134 p = barchart(Counts~as.factor(Size)|factor(Dataset, levels=unique(Dataset))+Chromosome, data = df, origin = 0, 136 pch=19,
135 horizontal=FALSE, 137 cex=0.35,
136 group=Polarity, 138 scales= list(relation="free", x=list(rot=0, cex=0.7, tck=0.5), y=list(tick.number=4, rot=90, cex=0.7)),
137 stack=TRUE, 139 xlab=NULL, main=NULL, ylab=NULL, ylim=ylimits,
138 col=c('red', 'blue'), 140 as.table=T,
139 scales=list(y=list(tick.number=4, rot=90, relation="free", cex=0.7), x=list(rot=0, cex=0.7, axs="i", tck=0.5)), 141 origin = 0,
140 xlab = NULL, 142 horizontal=FALSE,
141 ylab = NULL, 143 group=Polarity,
142 main = NULL, 144 col=c("red","blue"),
143 as.table=TRUE, 145 par.strip.text = list(cex=0.7),
144 par.strip.text = list(cex=0.6), 146 ...)
145 ...) 147 } else {
146 } 148 p = barchart(Counts~as.factor(Size)|factor(Dataset, levels=unique(Dataset))+Chromosome, data = df, origin = 0,
147 combineLimits(p) 149 horizontal=FALSE,
148 } 150 group=Polarity,
151 stack=TRUE,
152 col=c('red', 'blue'),
153 scales=list(y=list(rot=90, relation="free", cex=0.7), x=list(rot=0, cex=0.7, axs="i", tck=c(1,0))),
154 xlab = NULL,
155 ylab = NULL,
156 main = NULL,
157 as.table=TRUE,
158 par.strip.text = list(cex=0.6),
159 ...)
160 p=combineLimits(p)
161 }
162 return(p)
163 }
164
149 165
150 ## function parameters 166 ## function parameters
151 167
152 #par.settings.firstplot = list(layout.heights=list(top.padding=11, bottom.padding = -14)) 168 #par.settings.firstplot = list(layout.heights=list(top.padding=11, bottom.padding = -14))
153 #par.settings.secondplot=list(layout.heights=list(top.padding=11, bottom.padding = -15), strip.background=list(col=c("lavender","deepskyblue"))) 169 #par.settings.secondplot=list(layout.heights=list(top.padding=11, bottom.padding = -15), strip.background=list(col=c("lavender","deepskyblue")))
154 par.settings.firstplot = list(layout.heights=list(top.padding=-2, bottom.padding=-2)) 170 par.settings.firstplot = list(layout.heights=list(top.padding=-2, bottom.padding=-2),strip.background=list(col=c("lightblue","lightgreen")))
155 par.settings.secondplot=list(layout.heights=list(top.padding=-1, bottom.padding=-1), strip.background=list(col=c("lavender","deepskyblue"))) 171 par.settings.secondplot=list(layout.heights=list(top.padding=-1, bottom.padding=-1),strip.background=list(col=c("lightblue","lightgreen")))
156 par.settings.single_plot=list(strip.background = list(col = c("lightblue", "lightgreen")))
157 title_first_method = list(Counts="Read Counts", Coverage="Coverage depths", Median="Median sizes", Mean="Mean sizes", Size="Size Distributions") 172 title_first_method = list(Counts="Read Counts", Coverage="Coverage depths", Median="Median sizes", Mean="Mean sizes", Size="Size Distributions")
158 title_extra_method = list(Counts="Read Counts", Coverage="Coverage depths", Median="Median sizes", Mean="Mean sizes", Size="Size Distributions") 173 title_extra_method = list(Counts="Read Counts", Coverage="Coverage depths", Median="Median sizes", Mean="Mean sizes", Size="Size Distributions")
159 legend_first_method =list(Counts="Read count", Coverage="Coverage depth", Median="Median size", Mean="Mean size", Size="Read count") 174 legend_first_method =list(Counts="Read count", Coverage="Coverage depth", Median="Median size", Mean="Mean size", Size="Read count")
160 legend_extra_method =list(Counts="Read count", Coverage="Coveragedepth", Median="Median size", Mean="Mean size", Size="Read count") 175 legend_extra_method =list(Counts="Read count", Coverage="Coverage depth", Median="Median size", Mean="Mean size", Size="Read count")
161 bottom_first_method =list(Counts="Coordinates (nbre of bases)",Coverage="Coordinates (nbre of bases)", Median="Coordinates (nbre of bases)", Mean="Coordinates (nbre of bases)", Size="Sizes of reads") 176 bottom_first_method =list(Counts="Coordinates (nbre of bases)",Coverage="Coordinates (nbre of bases)", Median="Coordinates (nbre of bases)", Mean="Coordinates (nbre of bases)", Size="Sizes of reads")
162 bottom_extra_method =list(Counts="Coordinates (nbre of bases)",Coverage="Coordinates (nbre of bases)", Median="Coordinates (nbre of bases)", Mean="Coordinates (nbre of bases)", Size="Sizes of reads") 177 bottom_extra_method =list(Counts="Coordinates (nbre of bases)",Coverage="Coordinates (nbre of bases)", Median="Coordinates (nbre of bases)", Mean="Coordinates (nbre of bases)", Size="Sizes of reads")
163 178
164 ## Plotting Functions 179 ## Plotting Functions
165 180
166 double_plot <- function(...) { 181 double_plot <- function(...) {
167 page_height = 15 182 page_height = 15
168 rows_per_page = 10 183 rows_per_page = 10
169 graph_heights=c(40,30,40,30,40,30,40,30,40,30,10) 184 graph_heights=c(40,30,40,30,40,30,40,30,40,30,10)
170 if (n_samples > 4) {page_width = 8.2677*n_samples/4} else {page_width = 2.3*n_samples +2.5} 185 page_width=8.2677 * n_samples / 2
171 pdf(file=args$output_pdf, paper="special", height=page_height, width=page_width) 186 pdf(file=args$output_pdf, paper="special", height=page_height, width=page_width)
172 for (i in seq(1,n_genes,rows_per_page/2)) { 187 for (i in seq(1,n_genes,rows_per_page/2)) {
173 start=i 188 start=i
174 end=i+rows_per_page/2-1 189 end=i+rows_per_page/2-1
175 if (end>n_genes) {end=n_genes} 190 if (end>n_genes) {end=n_genes}
176 if (end-start+1 < 5) {graph_heights=c(rep(c(40,30),end-start+1),10,rep(c(40,30),5-(end-start+1)))} 191 if (end-start+1 < 5) {graph_heights=c(rep(c(40,30),end-start+1),10,rep(c(40,30),5-(end-start+1)))}
177 first_plot.list = lapply(per_gene_readmap[start:end], function(x) plot_unit(x, strip=FALSE, par.settings=par.settings.firstplot)) 192 first_plot.list = lapply(per_gene_readmap[start:end], function(x) update(useOuterStrips(plot_unit(x, par.settings=par.settings.secondplot), strip.left=strip.custom(par.strip.text = list(cex=0.5)))))
178 second_plot.list = lapply(per_gene_size[start:end], function(x) plot_unit(x, method=args$extra_plot_method, par.settings=par.settings.secondplot)) 193 second_plot.list = lapply(per_gene_size[start:end], function(x) update(useOuterStrips(plot_unit(x, method=args$extra_plot_method, par.settings=par.settings.firstplot), strip.left=strip.custom(par.strip.text = list(cex=0.5)), strip=FALSE)))
179 plot.list=rbind(second_plot.list, first_plot.list) 194 plot.list=rbind(first_plot.list, second_plot.list)
180 args_list=c(plot.list, list( nrow=rows_per_page+1, ncol=1, heights=unit(graph_heights, rep("mm", 11)), 195 args_list=c(plot.list, list( nrow=rows_per_page+1, ncol=1, heights=unit(graph_heights, rep("mm", 11)),
181 top=textGrob(paste(title_first_method[[args$first_plot_method]], "and", title_extra_method[[args$extra_plot_method]]), gp=gpar(cex=1), vjust=0, just="top"), 196 top=textGrob(paste(title_first_method[[args$first_plot_method]], "and", title_extra_method[[args$extra_plot_method]]), gp=gpar(cex=1), vjust=0, just="top"),
182 left=textGrob(paste(legend_first_method[[args$first_plot_method]], "/", legend_extra_method[[args$extra_plot_method]]), gp=gpar(cex=1), just=0.675*(end-start-(2.2*(4/2.7))),vjust=2, rot=90), 197 left=textGrob(paste(legend_first_method[[args$first_plot_method]], "/", legend_extra_method[[args$extra_plot_method]]), gp=gpar(cex=1), vjust=0, hjust=0, x=1, y=(-0.38/4)*(end-start-(3.28/0.38)), rot=90),
183 sub=textGrob(paste(bottom_first_method[[args$first_plot_method]], "/", bottom_extra_method[[args$extra_plot_method]]), gp=gpar(cex=1), just="bottom", vjust=2) 198 sub=textGrob(paste(bottom_first_method[[args$first_plot_method]], "/", bottom_extra_method[[args$extra_plot_method]]), gp=gpar(cex=1), just="bottom", vjust=2)
184 ) 199 )
185 ) 200 )
186 do.call(grid.arrange, args_list) 201 do.call(grid.arrange, args_list)
187 } 202 }
188 devname=dev.off() 203 devname=dev.off()
189 } 204 }
190 205
191 206
192 single_plot <- function(...) { 207 single_plot <- function(...) {
193 width = 8.2677 * n_samples / 2 208 width = 8.2677 * n_samples / 2
194 rows_per_page=8 209 rows_per_page=8
195 graph_heights=c(rep(40,8),10) 210 graph_heights=c(rep(40,8),10)
196 pdf(file=args$output_pdf, paper="special", height=15, width=width) 211 pdf(file=args$output_pdf, paper="special", height=15, width=width)
197 for (i in seq(1,n_genes,rows_per_page)) { 212 for (i in seq(1,n_genes,rows_per_page)) {
198 start=i 213 start=i
199 end=i+rows_per_page-1 214 end=i+rows_per_page-1
200 if (end>n_genes) {end=n_genes} 215 if (end>n_genes) {end=n_genes}
201 if (end-start+1 < 8) {graph_heights=c(rep(c(40),end-start+1),10,rep(c(40),8-(end-start+1)))} 216 if (end-start+1 < 8) {graph_heights=c(rep(c(40),end-start+1),10,rep(c(40),8-(end-start+1)))}
202 first_plot.list = lapply(per_gene_readmap[start:end], function(x) plot_unit(x, par.settings=par.settings.firstplot)) 217 first_plot.list = lapply(per_gene_readmap[start:end], function(x) update(useOuterStrips(plot_unit(x, par.settings=par.settings.firstplot),strip.left=strip.custom(par.strip.text = list(cex=0.5)))))
203 plot.list=rbind(first_plot.list) 218 plot.list=rbind(first_plot.list)
204 args_list=c(plot.list, list( nrow=rows_per_page+1, ncol=1, heights=unit(graph_heights, rep("mm", 9)), 219 args_list=c(plot.list, list( nrow=rows_per_page+1, ncol=1, heights=unit(graph_heights, rep("mm", 9)),
205 top=textGrob(title_first_method[[args$first_plot_method]], gp=gpar(cex=1), vjust=0, just="top"), 220 top=textGrob(title_first_method[[args$first_plot_method]], gp=gpar(cex=1), vjust=0, just="top"),
206 left=textGrob(legend_first_method[[args$first_plot_method]], gp=gpar(cex=1), just=(6.4/7)*(end-start-(6.2*(7/6.4))),vjust=2, rot=90), 221 left=textGrob(legend_first_method[[args$first_plot_method]], gp=gpar(cex=1), vjust=0, hjust=0, x=1, y=(-0.41/7)*(end-start-(6.23/0.41)), rot=90),
207 sub=textGrob(bottom_first_method[[args$first_plot_method]], gp=gpar(cex=1), just="bottom", vjust=2) 222 sub=textGrob(bottom_first_method[[args$first_plot_method]], gp=gpar(cex=1), just="bottom", vjust=2)
208 ) 223 )
209 ) 224 )
210 do.call(grid.arrange, args_list) 225 do.call(grid.arrange, args_list)
211 } 226 }
212 devname=dev.off() 227 devname=dev.off()
213 } 228 }
214 229
215 # main 230 # main
216 231
217 if (args$extra_plot_method != '') { double_plot() } 232 if (args$extra_plot_method != '') { double_plot() }
218 if (args$extra_plot_method == '' & !exists('global', where=args)) { 233 if (args$extra_plot_method == '' & !exists('global', where=args)) {
219 single_plot() 234 single_plot()
220 } 235 }
221 if (exists('global', where=args)) { 236 if (exists('global', where=args)) {
222 pdf(file=args$output, paper="special", height=11.69) 237 pdf(file=args$output, paper="special", height=11.69)
223 Table <- within(Table, Counts[Polarity=="R"] <- abs(Counts[Polarity=="R"])) # retropedalage 238 Table <- within(Table, Counts[Polarity=="R"] <- abs(Counts[Polarity=="R"])) # retropedalage
224 library(reshape2) 239 library(reshape2)
225 ml = melt(Table, id.vars = c("Dataset", "Chromosome", "Polarity", "Size")) 240 ml = melt(Table, id.vars = c("Dataset", "Chromosome", "Polarity", "Size"))
226 if (args$global == "nomerge") { 241 if (args$global == "nomerge") {
227 castml = dcast(ml, Dataset+Polarity+Size ~ variable, function(x) sum(x)) 242 castml = dcast(ml, Dataset+Polarity+Size ~ variable, function(x) sum(x))
228 castml <- within(castml, Counts[Polarity=="R"] <- (Counts[Polarity=="R"]*-1)) 243 castml <- within(castml, Counts[Polarity=="R"] <- (Counts[Polarity=="R"]*-1))
229 bc = globalbc(castml, global="no") 244 bc = globalbc(castml, global="no")
230 } else { 245 } else {
231 castml = dcast(ml, Dataset+Size ~ variable, function(x) sum(x)) 246 castml = dcast(ml, Dataset+Size ~ variable, function(x) sum(x))
232 bc = globalbc(castml, global="yes") 247 bc = globalbc(castml, global="yes")
233 } 248 }
234 plot(bc) 249 plot(bc)
235 devname=dev.off() 250 devname=dev.off()
236 } 251 }
237 252
238
239
240
241
242
243
244
245
246