ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/osprai/osprai/trunk/mdl_module.py
Revision: 29
Committed: Thu May 20 01:14:55 2010 UTC (9 years, 3 months ago) by clausted
File size: 6940 byte(s)
Log Message:
Add a parameter ('cofa') to the 1to1 models that represents a concentration factor.  This allows for easier analysis of a dilution series.  For example, if two ROIs contain data for the same sample but at a dilutions of 1:2 and 1:4, the set 'conc' value to the molar stock concentration and set the two 'cofa' values to 0.5 and 0.25.
Line User Rev File contents
1 clausted 17 """
2     mdl: Example model functions module for SPRI data.
3     Christopher Lausted, Institute for Systems Biology,
4     OSPRAI developers
5 clausted 29 Last modified on 100518 (yymmdd)
6 clausted 17
7 clausted 22 Examples:
8 clausted 17 #import mdl_module as mdl
9 clausted 22 #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 clausted 17 """
25 clausted 29 __version__ = "100518"
26 clausted 17
27    
28     ## Import libraries
29     import numpy as np
30 clausted 18 from copy import deepcopy
31 clausted 17
32    
33 clausted 18 def drift(time, data, params):
34 clausted 17 """
35     This function simply models a constant signal drift in units/second.
36 clausted 18 It requires numpy arrays of times and starting data values,
37 clausted 17 It only requires one parameter in the params list.
38 clausted 18 params['rate']['value']
39 clausted 17 """
40 clausted 18 y = np.zeros(len(time), dtype=float)
41 clausted 17 try:
42 clausted 18 rate = params['rate']['value']
43 clausted 17 except KeyError:
44     print "Error: Rate parameter not supplied to drift() model."
45 clausted 18 return y
46 clausted 17
47 clausted 18 ## 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 clausted 19 ## End of drift() function.
55 clausted 18
56    
57     def simple1to1(time, data, params):
58     """
59 clausted 26 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 clausted 29 params['conc']['value'] is stock concentration of analyte [A], (M)
67 clausted 26 params['kon']['value'] is on-rate of analyte, (1/Ms)
68     params['t2']['value'] is time of end binding & begin washing, (s)
69 clausted 29 params['koff']['value'] is off-rate of analyte, (1/s)
70 clausted 26 params['t3']['value'] is time end of washing & data fitting, (s)
71 clausted 29
72     It also can take an optional parameter useful for dilution series.
73     params['cofa']['value'] is concentration factor, (1/dilution factor)
74 clausted 18 """
75     ## Skip parameter validation steps for now.
76 clausted 29 t1 = float(params['t1']['value'])
77     rmax = float(params['rmax']['value'])
78     conc = float(params['conc']['value'])
79     kon = float(params['kon']['value'])
80     t2 = float(params['t2']['value'])
81     koff = float(params['koff']['value'])
82     t3 = float(params['t3']['value'])
83     if ('cofa' in params.keys()):
84     conc *= float(params['cofa']['value'])
85 clausted 17
86 clausted 26 """
87     Derivation:
88     d[AL]/dt = kon*[A]*[L] - koff*[AL]
89     y = [AL]
90     (rmax-y) = [L]
91     conc = [A]
92     rmax = [AL] + [L]
93     dy/dt = conc*kon*(rmax-y) - koff*y
94     """
95    
96 clausted 18 y = np.zeros(len(time), dtype=float)
97     base = y[0] = data[0] ## Baseline SPR signal.
98    
99     ## Must iterate through data, numerical integration.
100     for i in range(1,len(time)):
101 clausted 26 dt = time[i] - time[i-1]
102     ## Treat function as having three separate phases.
103 clausted 18 if (time[i] <= t1):
104     ## Pre-binding phase.
105     base = y[i] = data[i]
106     elif (t1 < time[i] <= t2):
107     ## Binding phase
108     yb = y[i-1] - base
109     dy = conc*kon*(rmax-yb) - koff*yb
110     if (abs(y[i-1]) > 999999999): dy = 0 ## Is this useful?
111 clausted 26 y[i] = y[i-1] + dy*dt
112 clausted 18 elif (t2 < time[i] <= t3):
113 clausted 26 ## Dissociation (conc=0)
114 clausted 18 yb = y[i-1]-base
115     dy = 0 - koff*yb
116 clausted 26 y[i] = y[i-1] + dy*dt
117 clausted 18
118     return y
119 clausted 19 ## End of simple1to1() function
120 clausted 18
121 clausted 19
122 clausted 26 def simple1to1_mtl(time, data, params):
123     """
124     This function simply models a 1:1 interaction with mass transport limitation.
125     [Abulk] --km-> [Asurf] + [L] --kon--> [AL]
126     [Abulk] <-km-- [Asurf] + [L] <-koff-- [AL]
127    
128     It requires numpy arrays of times and starting data values,
129     params['t1']['value'] is time of injection for binding, (s)
130     params['rmax']['value'] is maximum response, (RIU)
131     params['conc']['value'] is concentration of analyte [Abulk], (M)
132     params['kon']['value'] is on-rate of analyte, (1/Ms)
133     params['t2']['value'] is time of end binding & begin washing, (s)
134     params['koff']['value'] is off-rate of analyte, (1/s)
135     params['t3']['value'] is time end of washing & data fitting, (s)
136     params['kmtl']['value'] is rate of diffusion, (RIU/Ms)
137 clausted 29
138     It also can take an optional parameter useful for dilution series.
139     params['cofa']['value'] is concentration factor, (1/dilution factor)
140 clausted 26 """
141     ## Skip parameter validation steps for now.
142     t1 = params['t1']['value']
143     rmax = params['rmax']['value']
144     conc = params['conc']['value']
145     kon = params['kon']['value']
146     t2 = params['t2']['value']
147     koff = params['koff']['value']
148     t3 = params['t3']['value']
149     kmtl = params['kmtl']['value']
150 clausted 29 if ('cofa' in params.keys()):
151     conc *= float(params['cofa']['value'])
152    
153 clausted 26 """
154     Derivation:
155     d[AL]/dt = (kon*[A]*[L] - koff*[AL]) / (1 + kon*[L]/kmtl)
156     y = [AL]
157     (rmax-y) = [L]
158     conc = [Abulk]
159     rmax = [AL] + [L]
160     dy/dt = (kon*conc*(rmax-y) - koff*y) / (1 + kon*(rmax-y)/kmtl)
161     """
162    
163     stat = {} ## Error status dictionary.
164     y = np.zeros(len(time), dtype=float)
165     base = y[0] = data[0] ## Baseline SPR signal.
166    
167     ## Error checks.
168     if (kmtl<=0):
169     ## Avoid div/0, assume user wants no mass transport limitation.
170     stat['kmtl must be > 0'] = True
171     kmtl = 1e40 ## An arbitrary very large number.
172    
173     ## Must iterate through data, numerical integration.
174     for i in range(1,len(time)):
175     dt = time[i] - time[i-1]
176     ## Treat function as having three separate phases.
177     if (time[i] <= t1):
178     ## Pre-binding phase.
179     base = y[i] = data[i]
180     elif (t1 < time[i] <= t2):
181     ## Binding phase
182     yb = y[i-1] - base
183     dy = (conc*kon*(rmax-yb) - koff*yb) / (1 + kon*(rmax-yb)/kmtl)
184     ## Check for integration errors producing wild results.
185     if (abs(y[i-1]) > 9e9): dy = 0; stat['y>9e9'] = True
186     y[i] = y[i-1] + dy*dt
187     elif (t2 < time[i] <= t3):
188     ## Dissociation (conc=0)
189     yb = y[i-1] - base
190     dy = (0 - koff*yb) / (1 + kon*(rmax-yb)/kmtl)
191     y[i] = y[i-1] + dy*dt
192    
193     if (len(stat.keys()) > 0): print "Errors in simple1to1:", stat.keys()
194    
195     return y
196     ## End of simple1to1_mtl() function
197    
198    
199    
200 clausted 19 ################################# End of module #################################