comparison hist_plot.py @ 3:76d4cbefff85 draft

"planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/purge_dups commit 5d56aa02b0f905507e1d98a2d74f0629b7591cd3"
author iuc
date Mon, 14 Jun 2021 18:01:05 +0000
parents
children
comparison
equal deleted inserted replaced
2:17b378303f2d 3:76d4cbefff85
1 #!/usr/bin/env python3
2 # read depth histogram plot
3
4 # imported from https://github.com/dfguan/purge_dups/blob/master/scripts/hist_plot.py
5
6 import argparse
7
8 import matplotlib as mpl
9 import matplotlib.pyplot as plt
10
11 mpl.use("Agg")
12
13
14 def col_hist(stat_fn, delim):
15 hists = []
16 # we consider the coverage histogram start with 0
17 with open(stat_fn) as f:
18 for ln in f:
19 lnlist = ln.strip().split(delim)
20 hists.append(int(lnlist[1]))
21 return hists
22
23
24 def get_cutoffs(con):
25 if con:
26 lnlst = []
27 with open(con) as f:
28 lnlst = f.readline().strip().split("\t")
29 if len(lnlst):
30 return [int(lnlst[0]), int(lnlst[3]), int(lnlst[5])]
31 else:
32 return []
33 else:
34 return []
35
36
37 def mk_plot(hists, cutoffs, ttle, xm, xM, ym, yM, out_fl):
38
39 if ttle is None:
40 ttle = "read depth histogram"
41 if xm is None:
42 xm = 0
43 if xM is None:
44 xM = len(hists) - 2 # ignore the last read depth count
45 if ym is None:
46 ym = 0
47 if yM is None:
48 yM = 1.2 * max(hists)
49 x = [t for t in range(xm, xM)]
50 plt.axis([xm, xM, ym, yM])
51 width = 8
52 height = 6
53 plt.figure(num=None, figsize=(width, height))
54 plt.plot(x, hists[xm:xM], label="l", color="blue") #
55 plt.xticks([z for z in range(xm, xM, 10)], fontsize=3)
56 # cutoffs
57 colors = ["r", "g", "c"]
58 if len(cutoffs):
59 for i in range(len(cutoffs)):
60 plt.text(cutoffs[i], 0, str(cutoffs[i]), fontsize=5, color=colors[i])
61 plt.axvline(x=cutoffs[i], linewidth=1, color=colors[i])
62
63 plt.title(ttle)
64 plt.gca().xaxis.grid(True, color="black", alpha=0.2)
65
66 # plt.grid(True, color="black", alpha=0.2)
67 # plt.gca().get_legend().remove()
68
69 plt.tight_layout()
70 plt.savefig(out_fl, dpi=300)
71
72
73 if __name__ == "__main__":
74 parser = argparse.ArgumentParser(description="read depth histogram plot")
75
76 parser.add_argument(
77 "-c",
78 "--cutoffs",
79 type=str,
80 action="store",
81 dest="con",
82 help="read depth cutoffs",
83 )
84 parser.add_argument(
85 "-y", "--ymin", type=int, action="store", dest="ymin", help="set ymin"
86 )
87 parser.add_argument(
88 "-x", "--xmin", type=int, action="store", dest="xmin", help="set xmin"
89 )
90 parser.add_argument(
91 "-Y", "--ymax", type=int, action="store", dest="ymax", help="set ymax"
92 )
93 parser.add_argument(
94 "-X", "--xmax", type=int, action="store", dest="xmax", help="set xmax"
95 )
96 parser.add_argument(
97 "-t",
98 "--title",
99 type=str,
100 action="store",
101 dest="title",
102 help="figure title [NULL]",
103 default="",
104 )
105 parser.add_argument(
106 "-d",
107 "--delim",
108 type=str,
109 action="store",
110 dest="delim",
111 help="delimiter",
112 default="\t",
113 )
114 parser.add_argument("-v", "--version", action="version", version="hist_plot 0.0.0")
115 parser.add_argument("stat_fn", type=str, action="store", help="stat file")
116 parser.add_argument("out_fn", type=str, action="store", help="output file")
117 opts = parser.parse_args()
118 hists = col_hist(opts.stat_fn, opts.delim)
119 cutoffs = get_cutoffs(opts.con)
120 mk_plot(
121 hists,
122 cutoffs,
123 opts.title,
124 opts.xmin,
125 opts.xmax,
126 opts.ymin,
127 opts.ymax,
128 opts.out_fn,
129 )