ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/osprai/osprai/trunk/mdl_module.py
(Generate patch)
# Line 2 | Line 2
2   mdl: Example model functions module for SPRI data.
3   Christopher Lausted, Institute for Systems Biology,
4   OSPRAI developers
5 < Last modified on 100421 (yymmdd)
5 > Last modified on 100518 (yymmdd)
6  
7 < Example:
7 > Examples:
8   #import mdl_module as mdl
9 < #params = {'Rate': {'value':1, 'min':-100.0, 'max':100.0, 'fixed':False} }
10 < #params = dict(Rate = dict(value=1, min=0, max=10, fixed=False))
11 < #times = range(100)
12 < #data1 = drift(time, data0, params)
13 < print data1
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__ = "100421"
25 > __version__ = "100518"
26  
27  
28   ## Import libraries
# Line 41 | Line 51
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
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 stock 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 >    It also can take an optional parameter useful for dilution series.
73 >    params['cofa']['value'] is concentration factor, (1/dilution factor)
74 >    """
75 >    ## Skip parameter validation steps for now.
76 >    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 >    
86 >    """
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 >    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 >        dt = time[i] - time[i-1]
102 >        ## Treat function as having three separate phases.
103 >        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 >            y[i] = y[i-1] + dy*dt
112 >        elif (t2 < time[i] <= t3):
113 >            ## Dissociation (conc=0)
114 >            yb = y[i-1]-base
115 >            dy = 0 - koff*yb
116 >            y[i] = y[i-1] + dy*dt
117 >            
118 >    return y
119 >    ## End of simple1to1() function
120 >
121 >
122 > 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
130 <    params['rmax']['value'] is maximum response.
131 <    params['conc']['value'] is time of concentration of analyte
132 <    params['kon']['value']  is on-rate of analyte
133 <    params['t2']['value']   is time of end binding / begin washing
134 <    params['koff']['value']  is off-rate of analyte
135 <    params['t3']['value']   is time end of washing / data fitting.
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 >    
138 >    It also can take an optional parameter useful for dilution series.
139 >    params['cofa']['value'] is concentration factor, (1/dilution factor)
140      """
141      ## Skip parameter validation steps for now.
142      t1 = params['t1']['value']
# Line 62 | Line 145
145      kon = params['kon']['value']
146      t2 = params['t2']['value']
147      koff = params['koff']['value']
148 <    t3 = params['t3']['value']
148 >    t3 = params['t3']['value']
149 >    kmtl = params['kmtl']['value']
150 >    if ('cofa' in params.keys()):
151 >        conc *= float(params['cofa']['value'])
152 >        
153 >    """
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
184 <            if (abs(y[i-1]) > 999999999): dy = 0  ## Is this useful?
185 <            y[i] = y[i-1] + dy * (time[i] - time[i-1])
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
189 <            yb = y[i-1]-base
190 <            dy = 0 - koff*yb
191 <            y[i] = y[i-1] + dy * (time[i] - time[i-1])
192 <            
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 + ################################# End of module #################################

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines