|
0
|
1 ## draw coverage depth graphics for n references
|
|
|
2 ## if n <= 8: return png
|
|
|
3 ## if n > 8: return pdf
|
|
|
4
|
|
|
5 args = commandArgs(trailingOnly=TRUE)
|
|
|
6 library("reshape")
|
|
|
7 library("ggplot2")
|
|
|
8 library("gridExtra") # Load some extensions to have multiplot pdf
|
|
|
9 library("stringr") # To replace characters in strings
|
|
|
10
|
|
|
11 cov_file = args[1]
|
|
|
12 output = args[2]
|
|
|
13 coverage = read.table(cov_file,sep="\t", header=FALSE)
|
|
|
14 coverage = rename(coverage,c(V1="Chr", V2="locus",V3="depth"))
|
|
|
15
|
|
|
16 headers = unique(coverage$Chr)
|
|
|
17
|
|
|
18 ## round sup(length(each different ref)) / 8
|
|
|
19 ## nb_pages = ceiling(length(headers) / 8)
|
|
|
20 nb_pages = floor(length(headers) / 8)
|
|
|
21
|
|
|
22 print("nb_images:")
|
|
|
23 print(nb_pages)
|
|
|
24
|
|
|
25 output_p = ''
|
|
|
26
|
|
|
27 ## -------------------------------------------------------
|
|
|
28 ## plot each page (i more than 1 page)
|
|
|
29 if(nb_pages > 0){
|
|
|
30 for(n in 1:nb_pages-1){
|
|
|
31 print(n) ## start with 0
|
|
|
32
|
|
|
33 ## creates output file
|
|
|
34 if(n == 0){
|
|
|
35 output_p <- output
|
|
|
36 }else{
|
|
|
37 output_p <- str_replace(output, ".png", paste0("_", n ,".png") )
|
|
|
38 }
|
|
|
39 png(output_p)
|
|
|
40 print(paste0("output file:", output_p))
|
|
|
41 header_start_index = n*8 +1
|
|
|
42 header_end_index = header_start_index + 7
|
|
|
43 treated_headers = headers[header_start_index:header_end_index]
|
|
|
44 coverage_eight_headers <- coverage[coverage$Chr %in% treated_headers,]
|
|
|
45 print("header_start_index:")
|
|
|
46 print(header_start_index)
|
|
|
47 print("header_end_index:")
|
|
|
48 print(header_end_index)
|
|
|
49 print("treated_headers:")
|
|
|
50 print(treated_headers)
|
|
|
51 maxlocus = max( coverage_eight_headers$locus )
|
|
|
52 # # maxdepth = max( coverage_eight_headers$depth )
|
|
|
53 print("maxlocus:")
|
|
|
54 print(maxlocus)
|
|
|
55 print(paste0(output_p," file will be created"))
|
|
|
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))
|
|
|
57 print(plot8cov)
|
|
|
58 dev.off()
|
|
|
59 print(paste0(output_p," file created"))
|
|
|
60 }
|
|
|
61 }
|
|
|
62 ## creates output file
|
|
|
63 if(nb_pages == 0){
|
|
|
64 output_p <- output
|
|
|
65 }else{
|
|
|
66 output_p <- str_replace(output, ".png", paste0("_", nb_pages ,".png") )
|
|
|
67 }
|
|
|
68 png(output_p)
|
|
|
69
|
|
|
70 nb_elem_last_page = length(headers) %% 8
|
|
|
71 if(nb_elem_last_page > 0){
|
|
|
72 header_start_index = length(headers)-nb_elem_last_page+1
|
|
|
73 header_end_index = length(headers)
|
|
|
74 treated_headers = headers[header_start_index:header_end_index]
|
|
|
75 print("treated_headers:")
|
|
|
76 print(treated_headers)
|
|
|
77 coverage_last_headers= coverage[coverage$Chr %in% treated_headers,]
|
|
|
78 print("header_start_index:")
|
|
|
79 print(header_start_index)
|
|
|
80 print("header_end_index:")
|
|
|
81 print(header_end_index)
|
|
|
82 print("treated_headers:")
|
|
|
83 print(treated_headers)
|
|
|
84 maxlocus = max( coverage_last_headers$locus )
|
|
|
85 ## maxdepth = pmax( coverage_last_headers$depth )
|
|
|
86 print("maxlocus:")
|
|
|
87 print(maxlocus)
|
|
|
88 print(paste0(output_p," file will be created"))
|
|
|
89
|
|
|
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))
|
|
|
91 print(plot8cov)
|
|
|
92 dev.off()
|
|
|
93 print(paste0(output_p," file created"))
|
|
|
94 }
|
|
|
95 ##-------------------------------------------------------
|