Mercurial > repos > iuc > purge_dups
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 ) |
