Mercurial > repos > devteam > weightedaverage
annotate WeightedAverage.py @ 2:efa2b391e887 draft default tip
planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/weightedaverage commit f770c3c58f1e7e1fa5ed22d7f7aca856d36729e8
author | devteam |
---|---|
date | Wed, 05 Oct 2016 13:39:38 -0400 |
parents | 90611e86a998 |
children |
rev | line source |
---|---|
0 | 1 #!/usr/bin/env python |
2 """ | |
3 usage: %prog bed_file_1 bed_file_2 out_file | |
4 -1, --cols1=N,N,N,N: Columns for chr, start, end, strand in first file | |
5 -2, --cols2=N,N,N,N,N: Columns for chr, start, end, strand, name/value in second file | |
2
efa2b391e887
planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/weightedaverage commit f770c3c58f1e7e1fa5ed22d7f7aca856d36729e8
devteam
parents:
1
diff
changeset
|
6 -z, --allow_zeros: Include zeros in calculations |
0 | 7 """ |
8 | |
9 import collections | |
10 import sys | |
11 from galaxy.tools.util.galaxyops import * | |
12 from bx.cookbook import doc_optparse | |
13 | |
14 | |
15 #export PYTHONPATH=~/galaxy/lib/ | |
16 #running command python WeightedAverage.py interval_interpolate.bed value_interpolate.bed interpolate_result.bed | |
17 | |
18 def stop_err(msg): | |
19 sys.stderr.write(msg) | |
20 sys.exit() | |
21 | |
22 | |
23 def FindRate(chromosome, start_stop, dictType): | |
24 OverlapList = [] | |
25 for tempO in dictType[chromosome]: | |
26 DatabaseInterval = [tempO[0], tempO[1]] | |
27 Overlap = GetOverlap( start_stop, DatabaseInterval ) | |
28 if Overlap > 0: | |
29 OverlapList.append([Overlap, tempO[2]]) | |
30 | |
31 if len(OverlapList) > 0: | |
32 SumRecomb = 0 | |
33 SumOverlap = 0 | |
34 for member in OverlapList: | |
35 SumRecomb += member[0]*member[1] | |
36 SumOverlap += member[0] | |
37 averageRate = SumRecomb/SumOverlap | |
38 return averageRate | |
39 else: | |
40 return 'NA' | |
41 | |
42 | |
43 def GetOverlap(a, b): | |
44 return min(a[1], b[1])-max(a[0], b[0]) | |
45 | |
2
efa2b391e887
planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/weightedaverage commit f770c3c58f1e7e1fa5ed22d7f7aca856d36729e8
devteam
parents:
1
diff
changeset
|
46 def get_float_no_zero( field ): |
efa2b391e887
planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/weightedaverage commit f770c3c58f1e7e1fa5ed22d7f7aca856d36729e8
devteam
parents:
1
diff
changeset
|
47 rval = float( field ) |
efa2b391e887
planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/weightedaverage commit f770c3c58f1e7e1fa5ed22d7f7aca856d36729e8
devteam
parents:
1
diff
changeset
|
48 assert rval |
efa2b391e887
planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/weightedaverage commit f770c3c58f1e7e1fa5ed22d7f7aca856d36729e8
devteam
parents:
1
diff
changeset
|
49 return rval |
efa2b391e887
planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/weightedaverage commit f770c3c58f1e7e1fa5ed22d7f7aca856d36729e8
devteam
parents:
1
diff
changeset
|
50 |
0 | 51 |
52 options, args = doc_optparse.parse( __doc__ ) | |
53 | |
54 try: | |
55 chr_col_1, start_col_1, end_col_1, strand_col1 = parse_cols_arg( options.cols1 ) | |
56 chr_col_2, start_col_2, end_col_2, strand_col2, name_col_2 = parse_cols_arg( options.cols2 ) | |
57 input1, input2, input3 = args | |
58 except Exception, eee: | |
59 print eee | |
60 stop_err( "Data issue: click the pencil icon in the history item to correct the metadata attributes." ) | |
61 | |
2
efa2b391e887
planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/weightedaverage commit f770c3c58f1e7e1fa5ed22d7f7aca856d36729e8
devteam
parents:
1
diff
changeset
|
62 if options.allow_zeros: |
efa2b391e887
planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/weightedaverage commit f770c3c58f1e7e1fa5ed22d7f7aca856d36729e8
devteam
parents:
1
diff
changeset
|
63 get_value = float |
efa2b391e887
planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/weightedaverage commit f770c3c58f1e7e1fa5ed22d7f7aca856d36729e8
devteam
parents:
1
diff
changeset
|
64 else: |
efa2b391e887
planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/weightedaverage commit f770c3c58f1e7e1fa5ed22d7f7aca856d36729e8
devteam
parents:
1
diff
changeset
|
65 get_value = get_float_no_zero |
0 | 66 RecombChrDict = collections.defaultdict(list) |
67 | |
68 skipped = 0 | |
2
efa2b391e887
planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/weightedaverage commit f770c3c58f1e7e1fa5ed22d7f7aca856d36729e8
devteam
parents:
1
diff
changeset
|
69 for line in open( input2 ): |
0 | 70 temp = line.strip().split() |
71 try: | |
2
efa2b391e887
planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/weightedaverage commit f770c3c58f1e7e1fa5ed22d7f7aca856d36729e8
devteam
parents:
1
diff
changeset
|
72 value = get_value( temp[ name_col_2 ] ) |
efa2b391e887
planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/weightedaverage commit f770c3c58f1e7e1fa5ed22d7f7aca856d36729e8
devteam
parents:
1
diff
changeset
|
73 except Exception: |
0 | 74 skipped += 1 |
75 continue | |
2
efa2b391e887
planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/weightedaverage commit f770c3c58f1e7e1fa5ed22d7f7aca856d36729e8
devteam
parents:
1
diff
changeset
|
76 tempIndex = [ int( temp[ start_col_2 ] ), int( temp[ end_col_2 ] ), value ] |
efa2b391e887
planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/weightedaverage commit f770c3c58f1e7e1fa5ed22d7f7aca856d36729e8
devteam
parents:
1
diff
changeset
|
77 RecombChrDict[ temp[ chr_col_2 ] ].append( tempIndex ) |
0 | 78 |
79 print "Skipped %d features with invalid values" % (skipped) | |
80 | |
2
efa2b391e887
planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/weightedaverage commit f770c3c58f1e7e1fa5ed22d7f7aca856d36729e8
devteam
parents:
1
diff
changeset
|
81 fdd = open( input3, 'w' ) |
efa2b391e887
planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/weightedaverage commit f770c3c58f1e7e1fa5ed22d7f7aca856d36729e8
devteam
parents:
1
diff
changeset
|
82 for line in open( input1 ): |
efa2b391e887
planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/weightedaverage commit f770c3c58f1e7e1fa5ed22d7f7aca856d36729e8
devteam
parents:
1
diff
changeset
|
83 line = line.strip() |
efa2b391e887
planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/weightedaverage commit f770c3c58f1e7e1fa5ed22d7f7aca856d36729e8
devteam
parents:
1
diff
changeset
|
84 temp = line.split('\t') |
efa2b391e887
planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/weightedaverage commit f770c3c58f1e7e1fa5ed22d7f7aca856d36729e8
devteam
parents:
1
diff
changeset
|
85 chromosome = temp[ chr_col_1 ] |
efa2b391e887
planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/weightedaverage commit f770c3c58f1e7e1fa5ed22d7f7aca856d36729e8
devteam
parents:
1
diff
changeset
|
86 start = int( temp[ start_col_1 ] ) |
efa2b391e887
planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/weightedaverage commit f770c3c58f1e7e1fa5ed22d7f7aca856d36729e8
devteam
parents:
1
diff
changeset
|
87 stop = int( temp[ end_col_1 ] ) |
0 | 88 start_stop = [start, stop] |
89 RecombRate = FindRate( chromosome, start_stop, RecombChrDict ) | |
90 try: | |
91 RecombRate = "%.4f" % (float(RecombRate)) | |
92 except: | |
93 RecombRate = RecombRate | |
2
efa2b391e887
planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/weightedaverage commit f770c3c58f1e7e1fa5ed22d7f7aca856d36729e8
devteam
parents:
1
diff
changeset
|
94 fdd.write( "%s\t%s\n" % ( line, RecombRate ) ) |
0 | 95 fdd.close() |