1 import math
3 """Returns DNA/DNA tm using nearest neighbor thermodynamics.
4
5 dnac is DNA concentration [nM]
6 saltc is salt concentration [mM].
7 rna=0 is for DNA/DNA (default), for RNA, rna should be 1.
8
9 Sebastian Bassi <sbassi@genesdigitales.com>"""
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28 dh=0
29 ds=0
30
31 def tercorr(stri):
32 deltah=0
33 deltas=0
34 if rna==0:
35
36
37 if stri.startswith('G') or stri.startswith('C'):
38 deltah=deltah-0.1
39 deltas=deltas+2.8
40 elif stri.startswith('A') or stri.startswith('T'):
41 deltah=deltah-2.3
42 deltas=deltas-4.1
43 if stri.endswith('G') or stri.endswith('C'):
44 deltah=deltah-0.1
45 deltas=deltas+2.8
46 elif stri.endswith('A') or stri.endswith('T'):
47 deltah=deltah-2.3
48 deltas=deltas-4.1
49 dhL=dh+deltah
50 dsL=ds+deltas
51 return dsL,dhL
52 elif rna==1:
53
54 if stri.startswith('G') or stri.startswith('C'):
55 deltah=deltah-3.61
56 deltas=deltas-1.5
57 elif stri.startswith('A') or stri.startswith('T') or \
58 stri.startswith('U'):
59 deltah=deltah-3.72
60 deltas=deltas+10.5
61 if stri.endswith('G') or stri.endswith('C'):
62 deltah=deltah-3.61
63 deltas=deltas-1.5
64 elif stri.endswith('A') or stri.endswith('T') or \
65 stri.endswith('U'):
66 deltah=deltah-3.72
67 deltas=deltas+10.5
68 dhL=dh+deltah
69 dsL=ds+deltas
70
71 return dsL,dhL
72
73 def overcount(st,p):
74 """Returns how many p are on st, works even for overlapping"""
75 ocu=0
76 x=0
77 while 1:
78 try:
79 i=st.index(p,x)
80 except ValueError:
81 break
82 ocu=ocu+1
83 x=i+1
84 return ocu
85
86 R=1.987
87 sup=s.upper()
88 vsTC,vh=tercorr(sup)
89 vs=vsTC
90
91 k=(dnac/4.0)*1e-9
92
93
94 if rna==0:
95
96
97 vh=vh+(overcount(sup,"AA"))*7.9+(overcount(sup,"TT"))*7.9+(overcount(sup,"AT"))*7.2+(overcount(sup,"TA"))*7.2+(overcount(sup,"CA"))*8.5+(overcount(sup,"TG"))*8.5+(overcount(sup,"GT"))*8.4+(overcount(sup,"AC"))*8.4
98 vh=vh+(overcount(sup,"CT"))*7.8+(overcount(sup,"AG"))*7.8+(overcount(sup,"GA"))*8.2+(overcount(sup,"TC"))*8.2
99 vh=vh+(overcount(sup,"CG"))*10.6+(overcount(sup,"GC"))*9.8+(overcount(sup,"GG"))*8+(overcount(sup,"CC"))*8
100 vs=vs+(overcount(sup,"AA"))*22.2+(overcount(sup,"TT"))*22.2+(overcount(sup,"AT"))*20.4+(overcount(sup,"TA"))*21.3
101 vs=vs+(overcount(sup,"CA"))*22.7+(overcount(sup,"TG"))*22.7+(overcount(sup,"GT"))*22.4+(overcount(sup,"AC"))*22.4
102 vs=vs+(overcount(sup,"CT"))*21.0+(overcount(sup,"AG"))*21.0+(overcount(sup,"GA"))*22.2+(overcount(sup,"TC"))*22.2
103 vs=vs+(overcount(sup,"CG"))*27.2+(overcount(sup,"GC"))*24.4+(overcount(sup,"GG"))*19.9+(overcount(sup,"CC"))*19.9
104 ds=vs
105 dh=vh
106
107 else:
108
109
110 vh=vh+(overcount(sup,"AA"))*6.82+(overcount(sup,"TT"))*6.6+(overcount(sup,"AT"))*9.38+(overcount(sup,"TA"))*7.69+(overcount(sup,"CA"))*10.44+(overcount(sup,"TG"))*10.5+(overcount(sup,"GT"))*11.4+(overcount(sup,"AC"))*10.2
111 vh=vh+(overcount(sup,"CT"))*10.48+(overcount(sup,"AG"))*7.6+(overcount(sup,"GA"))*12.44+(overcount(sup,"TC"))*13.3
112 vh=vh+(overcount(sup,"CG"))*10.64+(overcount(sup,"GC"))*14.88+(overcount(sup,"GG"))*13.39+(overcount(sup,"CC"))*12.2
113 vs=vs+(overcount(sup,"AA"))*19.0+(overcount(sup,"TT"))*18.4+(overcount(sup,"AT"))*26.7+(overcount(sup,"TA"))*20.5
114 vs=vs+(overcount(sup,"CA"))*26.9+(overcount(sup,"TG"))*27.8+(overcount(sup,"GT"))*29.5+(overcount(sup,"AC"))*26.2
115 vs=vs+(overcount(sup,"CT"))*27.1+(overcount(sup,"AG"))*19.2+(overcount(sup,"GA"))*32.5+(overcount(sup,"TC"))*35.5
116 vs=vs+(overcount(sup,"CG"))*26.7+(overcount(sup,"GC"))*36.9+(overcount(sup,"GG"))*32.7+(overcount(sup,"CC"))*29.7
117 ds=vs
118 dh=vh
119
120 ds=ds-0.368*(len(s)-1)*math.log(saltc/1e3)
121 tm=((1000* (-dh))/(-ds+(R * (math.log(k)))))-273.15
122
123
124 return tm
125