ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/osprai/osprai/trunk/modelbasic.py
Revision: 25
Committed: Wed Apr 28 20:22:24 2010 UTC (9 years, 6 months ago) by rjaynes
File size: 4038 byte(s)
Log Message:
Add py and obj files to allow modeling of more SPR experiments with converter and curvefitting modules.  This is the work of Yuhang Wan and Rui Hou.

1. In "converter.py": 
      Add the saving and reading function for the sprclass data object.
      Also add function "keyfile_read_fake" to provide default information for SPRit and ICM formats in case of the bug when do background_subtract.
      Fix the bugs in "background_subtract".
      Tested by DAM and ICM formats.
2. In model modules:
      "modelclass.py" is the parent class for all the other model classes that performs the theoretical simulating, loading and saving of the parameter or simulated data. Rui and I also add some other model modules like competing model, twostate model, parallel model, and the time variable concentrated models, where the simulated result is compared with Clamp's simulation to make sure the equations are correct. 
       The basicmodel and basicmodel_varyC class are tested. 
3. In "curvefitting.py":
      Add typical pipeline for operation. The examples are packed with the file. 
      Add function to show the Elapsed time for each fitting.
Line File contents
1 """
2 Provide the model for curvefitting.
3 This is a BASIC 1:1 REACTION model, where the reaction equations are:
4 L + A ==ka=> LA
5 L + A <=kd== LA
6
7 Yuhang Wan, Rui Hou
8 Last modified on 100427 (yymmdd) by YW
9
10 Typical Pipeline:(work together with "modelclass.py" the father class):
11 >import modelbasic as mb
12 >m = mb.basicmodel()
13 >m.check_attribute()
14 >m.wizard() or m.load_model('model.obj')
15 >m.simulate(t=range(1200))
16 >m.save_model('newmodel')
17 >m.save_parainfo('parainfo.txt')
18 >m.save_sim_data('sim.dat')
19
20 """
21 __version__ = "100427"
22
23 import numpy as np
24 import pylab as plt
25 import copy
26 import os
27 import pickle
28 from modelclass import *
29
30 class basicmodel(modelclass):
31 """The data attributes of the class:
32 parainfo: the parameter information
33 sim_data: the simulated data
34 """
35 def __init__ (self, parainfo=[], sim_data=[] ):
36 modelclass.__init__(self,parainfo,sim_data)
37
38 def wizard(self, ):
39 """ This function helps you to create a parameter information list.
40 In this basic 1:1 reaction model, there are 7 parameters:
41 rmax: Maximum analyte binding capacity(RU),
42 ka: Association rate constant(M-1S-1),
43 kd: Dissociation rate constant(S-1),
44 ca: Analyte concentration(M),
45 ton: Starting time for sample injection(s),
46 toff: Ending time for sample injection(s),
47 tout: Total time(s).
48 """
49 print ('')
50 pname = ['rmax','ka','kd','ca','ton','toff','tout']
51 parainfo = []
52 for i in range(7):
53 print i
54 tmp={}
55 tmp['name'] = pname[i]
56 vstr = raw_input("input the value of parameter '%s', if a series, seperate with ',': " %pname[i])
57 vtmp = vstr.strip().split(',')
58 tmp['number'] = len(vtmp)
59 if len(vtmp) == 1:
60 vtmp = float(vtmp[0])
61 else:
62 vtmp = map(float, vtmp)
63 tmp['value'] = vtmp
64 tmp['fixed'] = bool(input("Is '%s' fixed? (1/0): " %pname[i]))
65
66 parainfo.append(tmp)
67 ##parainfo.append({'name':'', 'value':0., 'fixed':0, 'number':0})
68 self.parainfo = parainfo
69 return
70
71 def function(self, t, paralist):
72 """ In this basic 1:1 model, this function calculates the theoretical
73 curve through the parameter list you give.
74 """
75 ## Assign the parameters to calculate the curve
76 for p in paralist:
77 if p['name'] == 'rmax': rmax = p['value']
78 elif p['name'] == 'ka': ka = p['value']
79 elif p['name'] == 'kd': kd = p['value']
80 elif p['name'] == 'ca': ca = p['value']
81 elif p['name'] == 'ton': ton = p['value']
82 elif p['name'] == 'toff': toff = p['value']
83 elif p['name'] == 'tout': tout = p['value']
84 else: print p['name'], p['value']
85 if type(ca) == list:
86 print "Error: This function can only generate data for a single concentration."
87 return
88
89 t1 = ton ##100 ##234 # Initial injection
90 t2 = toff ##178 ##570 # Wash, begin of the dissociation
91 t3 = tout ##400 ##840 # End wash
92
93 ## Must iterate through data, numerical integration.
94 g = np.zeros(len(t), dtype=float)
95
96 for i in range(1,len(t)):
97 if (t[i] > t3): break #Speed things up.
98 if (t1 < t[i] < t2):
99 ## Association
100 dG = ka*ca*(rmax-g[i-1]) - kd*g[i-1]
101 elif (t2 < t[i] < t3):
102 ## Dissociation
103 dG = 0 - kd*g[i-1]
104 else:
105 dG = 0
106 if (abs(g[i]) > 999999999): dG = 0
107
108 g[i] = g[i-1] + dG * (t[i] - t[i-1])
109
110 return g
111
112
113 ##### End of basicmodel class definition.