ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/osprai/osprai/trunk/mdl_module.py
Revision: 26
Committed: Fri May 14 22:36:12 2010 UTC (9 years, 3 months ago) by clausted
File size: 6413 byte(s)
Log Message:
Added mass-transport limitation model simple1to1_mtl.

Line File contents
1 """
2 mdl: Example model functions module for SPRI data.
3 Christopher Lausted, Institute for Systems Biology,
4 OSPRAI developers
5 Last modified on 100511 (yymmdd)
6
7 Examples:
8 #import mdl_module as mdl
9 #import numpy as np
10 #times = np.arange(100)
11 #data = np.zeros(100)
12 #
13 #params1 = dict(rate=dict(value=1.0, min=-100.0, max=100.0, fixed=True))
14 #data1 = drift(time, data, param1)
15 #
16 #param2 = {'t1': {'value':30.0, 'min':30.0, 'max':30.0, 'fixed':True} }
17 #param2['rmax'] = {'value': 100.0}
18 #param2['conc'] = {'value': 1e-6}
19 #param2['kon'] = {'value': 2e4}
20 #param2['t2'] = {'value': 150.0}
21 #param2['koff'] = {'value': 1e-3}
22 #param2['t3'] = {'value': 270.0}
23 #data2 = simple1to1(time, data, param2)
24 """
25 __version__ = "100511"
26
27
28 ## Import libraries
29 import numpy as np
30 from copy import deepcopy
31
32
33 def drift(time, data, params):
34 """
35 This function simply models a constant signal drift in units/second.
36 It requires numpy arrays of times and starting data values,
37 It only requires one parameter in the params list.
38 params['rate']['value']
39 """
40 y = np.zeros(len(time), dtype=float)
41 try:
42 rate = params['rate']['value']
43 except KeyError:
44 print "Error: Rate parameter not supplied to drift() model."
45 return y
46
47 ## Numerical integration.
48 y[0] = data[0]
49 for i in range(1,len(time)):
50 dy = rate * (time[i]-time[i-1])
51 y[i] = y[i-1] + dy
52
53 return y
54 ## End of drift() function.
55
56
57 def simple1to1(time, data, params):
58 """
59 This function simply models a 1:1 interaction.
60 [A] + [L] --kon--> [AL]
61 [A] + [L] <-koff-- [AL]
62
63 It requires numpy arrays of times and starting data values.
64 params['t1']['value'] is time of injection for binding, (s)
65 params['rmax']['value'] is maximum response, (RIU)
66 params['conc']['value'] is concentration of analyte [A], (M)
67 params['kon']['value'] is on-rate of analyte, (1/Ms)
68 params['t2']['value'] is time of end binding & begin washing, (s)
69 params['koff']['value'] is off-rate of analyte, (1/s)
70 params['t3']['value'] is time end of washing & data fitting, (s)
71 """
72 ## Skip parameter validation steps for now.
73 t1 = params['t1']['value']
74 rmax = params['rmax']['value']
75 conc = params['conc']['value']
76 kon = params['kon']['value']
77 t2 = params['t2']['value']
78 koff = params['koff']['value']
79 t3 = params['t3']['value']
80
81 """
82 Derivation:
83 d[AL]/dt = kon*[A]*[L] - koff*[AL]
84 y = [AL]
85 (rmax-y) = [L]
86 conc = [A]
87 rmax = [AL] + [L]
88 dy/dt = conc*kon*(rmax-y) - koff*y
89 """
90
91 y = np.zeros(len(time), dtype=float)
92 base = y[0] = data[0] ## Baseline SPR signal.
93
94 ## Must iterate through data, numerical integration.
95 for i in range(1,len(time)):
96 dt = time[i] - time[i-1]
97 ## Treat function as having three separate phases.
98 if (time[i] <= t1):
99 ## Pre-binding phase.
100 base = y[i] = data[i]
101 elif (t1 < time[i] <= t2):
102 ## Binding phase
103 yb = y[i-1] - base
104 dy = conc*kon*(rmax-yb) - koff*yb
105 if (abs(y[i-1]) > 999999999): dy = 0 ## Is this useful?
106 y[i] = y[i-1] + dy*dt
107 elif (t2 < time[i] <= t3):
108 ## Dissociation (conc=0)
109 yb = y[i-1]-base
110 dy = 0 - koff*yb
111 y[i] = y[i-1] + dy*dt
112
113 return y
114 ## End of simple1to1() function
115
116
117 def simple1to1_mtl(time, data, params):
118 """
119 This function simply models a 1:1 interaction with mass transport limitation.
120 [Abulk] --km-> [Asurf] + [L] --kon--> [AL]
121 [Abulk] <-km-- [Asurf] + [L] <-koff-- [AL]
122
123 It requires numpy arrays of times and starting data values,
124 params['t1']['value'] is time of injection for binding, (s)
125 params['rmax']['value'] is maximum response, (RIU)
126 params['conc']['value'] is concentration of analyte [Abulk], (M)
127 params['kon']['value'] is on-rate of analyte, (1/Ms)
128 params['t2']['value'] is time of end binding & begin washing, (s)
129 params['koff']['value'] is off-rate of analyte, (1/s)
130 params['t3']['value'] is time end of washing & data fitting, (s)
131 params['kmtl']['value'] is rate of diffusion, (RIU/Ms)
132 """
133 ## Skip parameter validation steps for now.
134 t1 = params['t1']['value']
135 rmax = params['rmax']['value']
136 conc = params['conc']['value']
137 kon = params['kon']['value']
138 t2 = params['t2']['value']
139 koff = params['koff']['value']
140 t3 = params['t3']['value']
141 kmtl = params['kmtl']['value']
142
143 """
144 Derivation:
145 d[AL]/dt = (kon*[A]*[L] - koff*[AL]) / (1 + kon*[L]/kmtl)
146 y = [AL]
147 (rmax-y) = [L]
148 conc = [Abulk]
149 rmax = [AL] + [L]
150 dy/dt = (kon*conc*(rmax-y) - koff*y) / (1 + kon*(rmax-y)/kmtl)
151 """
152
153 stat = {} ## Error status dictionary.
154 y = np.zeros(len(time), dtype=float)
155 base = y[0] = data[0] ## Baseline SPR signal.
156
157 ## Error checks.
158 if (kmtl<=0):
159 ## Avoid div/0, assume user wants no mass transport limitation.
160 stat['kmtl must be > 0'] = True
161 kmtl = 1e40 ## An arbitrary very large number.
162
163 ## Must iterate through data, numerical integration.
164 for i in range(1,len(time)):
165 dt = time[i] - time[i-1]
166 ## Treat function as having three separate phases.
167 if (time[i] <= t1):
168 ## Pre-binding phase.
169 base = y[i] = data[i]
170 elif (t1 < time[i] <= t2):
171 ## Binding phase
172 yb = y[i-1] - base
173 dy = (conc*kon*(rmax-yb) - koff*yb) / (1 + kon*(rmax-yb)/kmtl)
174 ## Check for integration errors producing wild results.
175 if (abs(y[i-1]) > 9e9): dy = 0; stat['y>9e9'] = True
176 y[i] = y[i-1] + dy*dt
177 elif (t2 < time[i] <= t3):
178 ## Dissociation (conc=0)
179 yb = y[i-1] - base
180 dy = (0 - koff*yb) / (1 + kon*(rmax-yb)/kmtl)
181 y[i] = y[i-1] + dy*dt
182
183 if (len(stat.keys()) > 0): print "Errors in simple1to1:", stat.keys()
184
185 return y
186 ## End of simple1to1_mtl() function
187
188
189
190 ################################# End of module #################################