forked from combogenomics/regtools
-
Notifications
You must be signed in to change notification settings - Fork 0
/
plotMethodDist
executable file
·97 lines (78 loc) · 2.18 KB
/
plotMethodDist
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
#!/usr/bin/python
import sys
if len(sys.argv) < 3:
print 'USAGE: plotMethodDist REGULATOR INFILE [[INFILE] ...]'
sys.exit(1)
reg = sys.argv[1]
infiles = sys.argv[2:]
logmethods = ['rsat', 'nhmmer', 'mast']
import math
meth = {}
pdist = set()
for f in infiles:
for l in open(f):
s = l.rstrip().split('\t')
try:
st = int(s[6])
en = int(s[7])
except:continue
sc = float(s[10])
me = s[9]
meth[me] = meth.get(me, set())
if me in logmethods:
meth[me].add(((st+en)/2.0, math.log10(sc)))
else:
meth[me].add(((st+en)/2.0, sc))
pdist.add((st+en)/2.0)
import matplotlib
matplotlib.use('Agg')
import matplotlib.pyplot as plt
fig = plt.figure(figsize=(8,8))
#matplotlib.rcParams['axes.color_cycle'] = ['r', 'b', 'r', 'g', 'k']
#plt.rc('axes', color_cycle=['b', 'r', 'g', 'k'])
colors = ['b', 'r', 'g', 'k']
import numpy as np
nmeth = np.array(meth.keys())
# Array reshape (tricky!)
# 1- Finger crossed for a perfect square
square = np.sqrt( len(nmeth) )
if square.is_integer():
# Yep
rows = cols = int(square)
else:
# Oh snap!
cols = int(square)
rows = len(nmeth) / float(cols)
if not rows.is_integer():
# Some fake signals will be added
rows = int(rows) + 1
rows = int(rows)
i = 0
for m in sorted(meth.keys()):
ax = fig.add_subplot(rows, cols, i)
j = i
while True:
if j > len(colors)-1:
j = i-len(colors)
else:
c = colors[j]
break
if m in logmethods:
ms = max([-x[1] for x in meth[m]])
mins = min([-x[1] for x in meth[m]])
try:
ax.plot([x[0] for x in meth[m]], [(-x[1]-mins)/(ms-mins) for x in meth[m]], '+%s'%c)
except:pass
else:
ms = max([x[1] for x in meth[m]])
mins = min([x[1] for x in meth[m]])
try:
ax.plot([x[0] for x in meth[m]], [(x[1]-mins)/(ms-mins) for x in meth[m]], '+%s'%c)
except:pass
ax.set_title('%s'%m)
ax.set_ylim(0, 1)
ax.set_xlim(min(pdist),0)
i += 1
plt.tight_layout()
plt.savefig('%s.png'%reg)
plt.savefig('%s.pdf'%reg)