annotate graph_coverage_multi_png.R @ 1:7bcef9bf4f18 draft default tip

Uploaded wrapper
author p.lucas
date Thu, 13 Jun 2024 08:32:27 +0000
parents 5558479694c7
children
Ignore whitespace changes - Everywhere: Within whitespace: At end of lines:
rev   line source
0
5558479694c7 Uploaded script
p.lucas
parents:
diff changeset
1 ## draw coverage depth graphics for n references
5558479694c7 Uploaded script
p.lucas
parents:
diff changeset
2 ## if n <= 8: return png
5558479694c7 Uploaded script
p.lucas
parents:
diff changeset
3 ## if n > 8: return pdf
5558479694c7 Uploaded script
p.lucas
parents:
diff changeset
4
5558479694c7 Uploaded script
p.lucas
parents:
diff changeset
5 args = commandArgs(trailingOnly=TRUE)
5558479694c7 Uploaded script
p.lucas
parents:
diff changeset
6 library("reshape")
5558479694c7 Uploaded script
p.lucas
parents:
diff changeset
7 library("ggplot2")
5558479694c7 Uploaded script
p.lucas
parents:
diff changeset
8 library("gridExtra") # Load some extensions to have multiplot pdf
5558479694c7 Uploaded script
p.lucas
parents:
diff changeset
9 library("stringr") # To replace characters in strings
5558479694c7 Uploaded script
p.lucas
parents:
diff changeset
10
5558479694c7 Uploaded script
p.lucas
parents:
diff changeset
11 cov_file = args[1]
5558479694c7 Uploaded script
p.lucas
parents:
diff changeset
12 output = args[2]
5558479694c7 Uploaded script
p.lucas
parents:
diff changeset
13 coverage = read.table(cov_file,sep="\t", header=FALSE)
5558479694c7 Uploaded script
p.lucas
parents:
diff changeset
14 coverage = rename(coverage,c(V1="Chr", V2="locus",V3="depth"))
5558479694c7 Uploaded script
p.lucas
parents:
diff changeset
15
5558479694c7 Uploaded script
p.lucas
parents:
diff changeset
16 headers = unique(coverage$Chr)
5558479694c7 Uploaded script
p.lucas
parents:
diff changeset
17
5558479694c7 Uploaded script
p.lucas
parents:
diff changeset
18 ## round sup(length(each different ref)) / 8
5558479694c7 Uploaded script
p.lucas
parents:
diff changeset
19 ## nb_pages = ceiling(length(headers) / 8)
5558479694c7 Uploaded script
p.lucas
parents:
diff changeset
20 nb_pages = floor(length(headers) / 8)
5558479694c7 Uploaded script
p.lucas
parents:
diff changeset
21
5558479694c7 Uploaded script
p.lucas
parents:
diff changeset
22 print("nb_images:")
5558479694c7 Uploaded script
p.lucas
parents:
diff changeset
23 print(nb_pages)
5558479694c7 Uploaded script
p.lucas
parents:
diff changeset
24
5558479694c7 Uploaded script
p.lucas
parents:
diff changeset
25 output_p = ''
5558479694c7 Uploaded script
p.lucas
parents:
diff changeset
26
5558479694c7 Uploaded script
p.lucas
parents:
diff changeset
27 ## -------------------------------------------------------
5558479694c7 Uploaded script
p.lucas
parents:
diff changeset
28 ## plot each page (i more than 1 page)
5558479694c7 Uploaded script
p.lucas
parents:
diff changeset
29 if(nb_pages > 0){
5558479694c7 Uploaded script
p.lucas
parents:
diff changeset
30 for(n in 1:nb_pages-1){
5558479694c7 Uploaded script
p.lucas
parents:
diff changeset
31 print(n) ## start with 0
5558479694c7 Uploaded script
p.lucas
parents:
diff changeset
32
5558479694c7 Uploaded script
p.lucas
parents:
diff changeset
33 ## creates output file
5558479694c7 Uploaded script
p.lucas
parents:
diff changeset
34 if(n == 0){
5558479694c7 Uploaded script
p.lucas
parents:
diff changeset
35 output_p <- output
5558479694c7 Uploaded script
p.lucas
parents:
diff changeset
36 }else{
5558479694c7 Uploaded script
p.lucas
parents:
diff changeset
37 output_p <- str_replace(output, ".png", paste0("_", n ,".png") )
5558479694c7 Uploaded script
p.lucas
parents:
diff changeset
38 }
5558479694c7 Uploaded script
p.lucas
parents:
diff changeset
39 png(output_p)
5558479694c7 Uploaded script
p.lucas
parents:
diff changeset
40 print(paste0("output file:", output_p))
5558479694c7 Uploaded script
p.lucas
parents:
diff changeset
41 header_start_index = n*8 +1
5558479694c7 Uploaded script
p.lucas
parents:
diff changeset
42 header_end_index = header_start_index + 7
5558479694c7 Uploaded script
p.lucas
parents:
diff changeset
43 treated_headers = headers[header_start_index:header_end_index]
5558479694c7 Uploaded script
p.lucas
parents:
diff changeset
44 coverage_eight_headers <- coverage[coverage$Chr %in% treated_headers,]
5558479694c7 Uploaded script
p.lucas
parents:
diff changeset
45 print("header_start_index:")
5558479694c7 Uploaded script
p.lucas
parents:
diff changeset
46 print(header_start_index)
5558479694c7 Uploaded script
p.lucas
parents:
diff changeset
47 print("header_end_index:")
5558479694c7 Uploaded script
p.lucas
parents:
diff changeset
48 print(header_end_index)
5558479694c7 Uploaded script
p.lucas
parents:
diff changeset
49 print("treated_headers:")
5558479694c7 Uploaded script
p.lucas
parents:
diff changeset
50 print(treated_headers)
5558479694c7 Uploaded script
p.lucas
parents:
diff changeset
51 maxlocus = max( coverage_eight_headers$locus )
5558479694c7 Uploaded script
p.lucas
parents:
diff changeset
52 # # maxdepth = max( coverage_eight_headers$depth )
5558479694c7 Uploaded script
p.lucas
parents:
diff changeset
53 print("maxlocus:")
5558479694c7 Uploaded script
p.lucas
parents:
diff changeset
54 print(maxlocus)
5558479694c7 Uploaded script
p.lucas
parents:
diff changeset
55 print(paste0(output_p," file will be created"))
5558479694c7 Uploaded script
p.lucas
parents:
diff changeset
56 plot8cov <- ggplot(coverage_eight_headers, aes(x=locus, y=depth )) + geom_line(color="red", size=1) + facet_wrap(~coverage_eight_headers$Chr, scales="free", ncol=2, nrow=4) + scale_x_continuous(expand = c(0, 0), limits = c(0, NA))
5558479694c7 Uploaded script
p.lucas
parents:
diff changeset
57 print(plot8cov)
5558479694c7 Uploaded script
p.lucas
parents:
diff changeset
58 dev.off()
5558479694c7 Uploaded script
p.lucas
parents:
diff changeset
59 print(paste0(output_p," file created"))
5558479694c7 Uploaded script
p.lucas
parents:
diff changeset
60 }
5558479694c7 Uploaded script
p.lucas
parents:
diff changeset
61 }
5558479694c7 Uploaded script
p.lucas
parents:
diff changeset
62 ## creates output file
5558479694c7 Uploaded script
p.lucas
parents:
diff changeset
63 if(nb_pages == 0){
5558479694c7 Uploaded script
p.lucas
parents:
diff changeset
64 output_p <- output
5558479694c7 Uploaded script
p.lucas
parents:
diff changeset
65 }else{
5558479694c7 Uploaded script
p.lucas
parents:
diff changeset
66 output_p <- str_replace(output, ".png", paste0("_", nb_pages ,".png") )
5558479694c7 Uploaded script
p.lucas
parents:
diff changeset
67 }
5558479694c7 Uploaded script
p.lucas
parents:
diff changeset
68 png(output_p)
5558479694c7 Uploaded script
p.lucas
parents:
diff changeset
69
5558479694c7 Uploaded script
p.lucas
parents:
diff changeset
70 nb_elem_last_page = length(headers) %% 8
5558479694c7 Uploaded script
p.lucas
parents:
diff changeset
71 if(nb_elem_last_page > 0){
5558479694c7 Uploaded script
p.lucas
parents:
diff changeset
72 header_start_index = length(headers)-nb_elem_last_page+1
5558479694c7 Uploaded script
p.lucas
parents:
diff changeset
73 header_end_index = length(headers)
5558479694c7 Uploaded script
p.lucas
parents:
diff changeset
74 treated_headers = headers[header_start_index:header_end_index]
5558479694c7 Uploaded script
p.lucas
parents:
diff changeset
75 print("treated_headers:")
5558479694c7 Uploaded script
p.lucas
parents:
diff changeset
76 print(treated_headers)
5558479694c7 Uploaded script
p.lucas
parents:
diff changeset
77 coverage_last_headers= coverage[coverage$Chr %in% treated_headers,]
5558479694c7 Uploaded script
p.lucas
parents:
diff changeset
78 print("header_start_index:")
5558479694c7 Uploaded script
p.lucas
parents:
diff changeset
79 print(header_start_index)
5558479694c7 Uploaded script
p.lucas
parents:
diff changeset
80 print("header_end_index:")
5558479694c7 Uploaded script
p.lucas
parents:
diff changeset
81 print(header_end_index)
5558479694c7 Uploaded script
p.lucas
parents:
diff changeset
82 print("treated_headers:")
5558479694c7 Uploaded script
p.lucas
parents:
diff changeset
83 print(treated_headers)
5558479694c7 Uploaded script
p.lucas
parents:
diff changeset
84 maxlocus = max( coverage_last_headers$locus )
5558479694c7 Uploaded script
p.lucas
parents:
diff changeset
85 ## maxdepth = pmax( coverage_last_headers$depth )
5558479694c7 Uploaded script
p.lucas
parents:
diff changeset
86 print("maxlocus:")
5558479694c7 Uploaded script
p.lucas
parents:
diff changeset
87 print(maxlocus)
5558479694c7 Uploaded script
p.lucas
parents:
diff changeset
88 print(paste0(output_p," file will be created"))
5558479694c7 Uploaded script
p.lucas
parents:
diff changeset
89
5558479694c7 Uploaded script
p.lucas
parents:
diff changeset
90 plot8cov <- ggplot(coverage_last_headers, aes(x=locus, y=depth )) + geom_line(color="red", size=1) + facet_wrap(~coverage_last_headers$Chr, scales="free", ncol=2, nrow=4) + scale_x_continuous(expand = c(0, 0), limits = c(0, NA))
5558479694c7 Uploaded script
p.lucas
parents:
diff changeset
91 print(plot8cov)
5558479694c7 Uploaded script
p.lucas
parents:
diff changeset
92 dev.off()
5558479694c7 Uploaded script
p.lucas
parents:
diff changeset
93 print(paste0(output_p," file created"))
5558479694c7 Uploaded script
p.lucas
parents:
diff changeset
94 }
5558479694c7 Uploaded script
p.lucas
parents:
diff changeset
95 ##-------------------------------------------------------