Mercurial > repos > tyty > structurefold

comparison reactivity_cal/react_norm_function.py @ 5:7a8ddf1819b1 draft

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
Uploaded
author tyty
date Mon, 15 Sep 2014 14:52:52 -0400
parents
children
comparison
equal deleted inserted replaced
4:a292aaf51735 5:7a8ddf1819b1
1 #!/usr/bin/env python
2 # -*- coding: utf-8 -*-
3 import sys
4 from Bio import SeqIO
5 import math
6 from parse_dis_react import *
7
8 def cap(a,value):
9 if a>=value:
10 return value
11 else:
12 return a
13
14 def react_norm(react_file, result_file, capped_value):
15 print("Normalizing.....")
16 react1 = parse_dist(react_file)
17 react = react1[1]
18 h = file(result_file, 'w')
19
20 capped = int(capped_value)
21
22 all_react = []
23
24
25 for t in react:
26 if react[t]!='null':
27 for i in range(len(react[t])):
28 if react[t][i]!='NA':
29 all_react.append(float(react[t][i]))
30 # except:
31 # print(react[t][i])
32 # print(t)
33 # print(i)
34
35 all_react.sort(reverse = True)
36 #print((all_react))
37 #print(all_react[int(len(all_react)*0.02)])
38 #print(all_react[int(len(all_react)*0.03)])
39 #print(all_react[int(len(all_react)*0.025)])
40 #print(all_react[int(len(all_react)*0.04)])
41 #print(all_react[int(len(all_react)*0.05)])
42 '''
43 mean = sum(all_react)/len(all_react)
44 print(mean)
45 temp = 0
46
47 for i in range(len(all_react)):
48 temp = temp+all_react[i]*all_react[i]
49 temp = temp/len(all_react)
50 sd = math.sqrt(temp-mean*mean)
51 '''
52 eight = all_react[int(len(all_react)*0.02):int(len(all_react)*0.1)]
53 meight = sum(eight)/len(eight)
54
55 for t in react:
56 h.write(t)
57 h.write('\n')
58 if react[t]!='null':
59 if (t.find('AT1G29930')==-1) and (t.find('At1g29930')==-1):
60 for i in range((len(react[t])-1)):
61 if react[t][i]!='NA':
62 h.write(str(cap((float(react[t][i])/meight),capped)))
63 else:
64 h.write('NA')
65 h.write('\t')
66 if react[t][i+1]!='NA':
67 h.write(str(cap((float(react[t][i+1])/meight),capped)))
68 else:
69 h.write('NA')
70 h.write('\n')
71 else:
72 for i in range((len(react[t])-1)):
73 if react[t][i]!='NA':
74 h.write(str(float(react[t][i])*2.6))
75 else:
76 h.write('NA')
77 h.write('\t')
78 if react[t][i+1]!='NA':
79 h.write(str(float(react[t][i])*2.6))
80 else:
81 h.write('NA')
82 h.write('\n')
83
84
85
86 h.close()
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114