瀏覽代碼

First working version of Hg model (no pregnancy yet)

Andrej 2 年之前
父節點
當前提交
8b7ee24093
共有 3 個文件被更改,包括 185 次插入162 次删除
  1. 23 14
      models/humanHG.json
  2. 51 16
      pythonScripts/cModel.py
  3. 111 132
      pythonScripts/compartmentModel.ipynb

+ 23 - 14
models/humanHG.json

@@ -1,10 +1,12 @@
 {
+   "timeUnit":"min",
    "compartments":[
-      "redBloodCells","arterial","venous","kidney","urine",
+      "redBloodCells","plasma","venous","kidney","urine",
       "richlyPerfused","fat","slowlyPerfused","hair",
       "brainBlood","brain","placenta",
       "liver","gut","intestine","feces","inorganicMercury"],
    "commentParameters":"every parameter must have a value, name is optional and used in SE computation",
+   "volumeComment":"in a compartment w/o volume, value is mass",
    "volumes":{
       "plasma":"plasmaVolume",
       "venous":"plasmaVolume",
@@ -15,7 +17,12 @@
       "brainBlood":"brainBloodVolume",
       "placenta":"placentaVolume",
       "gut":"gutVolume",
-      "intestine":"intestineVolume"},
+      "intestine":"intestineVolume",
+      "redBloodCells":"redBloodCellsVolume",
+      "hair":"hairVolume",
+      "brain":"brainVolume",
+      "liver":"liverVolume"
+      },
    "flows":{
       "(plasma,venous):kidney":"kidneyFlow",
       "(plasma,venous):richlyPerfused":"richlyPerfusedFlow",
@@ -30,8 +37,7 @@
    "commentSources":"constants in concentration units per s",
    "commentHeavyside":"value is mass of exogene/volume of container/duration",
    "sources":{
-      "arterial":{"name":"exponential","k":-0.005,"constant":0},
-      "venous":{"name":"Heavyside","limit":3,"value":1.22e4}},
+      "intestine":{"name":"constant","value":0.0486,"unit":"ug/min"}},
    "commentBindings":"split to flow and diffusion dominated",
    "bindings": {
       "flowComments":"group common sources or targets with brackets",
@@ -42,8 +48,10 @@
          "venous->plasma"],
       "diffusion":{
          "plasma->redBloodCells":"kRBC",
+         "redBloodCells->plasma":"kRBC",
          "kidney->urine":"kU",
          "slowlyPerfused->hair":"kH",
+         "hair->slowlyPerfused":"kH",
          "brainBlood->brain":"kBR",
          "liver->intestine":"kB",
          "liver->inorganicMercury":"kI",
@@ -62,8 +70,8 @@
        "liver":"liverPC",
        "gut":"gutPC",
        "brainBlood->brain":"brainPC",
-       "slowlyPerfused->hair":"hairPC",
-       "plasma->redBloodCells":"rbcPC"},
+       "hair->slowlyPerfused":"hairPC",
+       "redBloodCells->plasma":"rbcPC"},
      "parameters":{
         "kidneyPC":{"value":4.0},
         "richlyPerfusedPC":{"value":1},
@@ -77,6 +85,7 @@
         "gutPC":{"value":1},
         "rbcPC":{"value":12},
         "flowNotes":"assuming 70kg individual, fetus at 0.1kg",
+        "plasmaFlow":{"value":8.067,"unit":"l/min"},
         "kidneyFlow":{"value":1.412,"unit":"l/min"},
         "richlyPerfusedFlow":{"value":1.484,"unit":"l/min"},
         "fatFlow":{"value":0.419,"unit":"l/min"},
@@ -99,14 +108,14 @@
         "liverVolume":{"value":1.82,"unit":"kg"},
         "gutVolume":{"value":1.19,"unit":"kg"},
         "intestineVolume":{"value":0.98,"unit":"kg"},
-        "kRBC":{"value":0.025,"unit":"l/min"},
+        "kRBC":{"value":0.60,"unit":"l/min"},
         "kU":{"value":0,"unit":"l/min"},
-        "kH":{"value":1.167e-7,"unit":"l/min"},
-        "kBR":{"value":1.67e-4,"unit":"l/min"},
-        "kI":{"value":1.67e-7,"unit":"l/min"},
-        "kB":{"value":1.67e-6,"unit":"l/min"},
-        "kD":{"value":1.67e-6,"unit":"l/min"},
-        "kF":{"value":3.33e-6,"unit":"l/min"},
+        "kH":{"value":2.823e-6,"unit":"l/min"},
+        "kBR":{"value":4.033e-3,"unit":"l/min"},
+        "kI":{"value":4.033e-6,"unit":"l/min"},
+        "kB":{"value":4.033e-5,"unit":"l/min"},
+        "kD":{"value":4.033e-5,"unit":"l/min"},
+        "kF":{"value":8.066e-5,"unit":"l/min"},
         "kF1":{"value":0,"unit":"l/min"},
-        "kR":{"value":8.33e-5,"unit":"l/min"}
+        "kR":{"value":2.016e-3,"unit":"l/min"}}
 }

+ 51 - 16
pythonScripts/cModel.py

@@ -17,28 +17,46 @@ class model:
       self.compartments[compartmentName]['targets']={}
       self.compartments[compartmentName]['sensTargets']={}
 
+   def getTimeUnit(self):
+
+      try: 
+         return self.mod['timeUnit']
+      except KeyError:
+         return 's'
+
    def bind(self,sourceCompartment,targetCompartment,k,pc,scaleToVolume=0):
+      
+
 #establish a flow from one compartment to the other
       cSrc=self.compartments[sourceCompartment]
       cTarg=self.compartments[targetCompartment]
 
       vSPar=self.getVolumePar(sourceCompartment,scaleToVolume)
       vTPar=self.getVolumePar(targetCompartment,scaleToVolume)
+      
+      tu=self.getTimeUnit()
+      fk=get(tu,k)
+      fpc=get(tu,pc)
+      fvs=get(tu,vSPar)
+      fvt=get(tu,vTPar)
 
       #the source equation (where we subtract the current)
-      addValue(cSrc['targets'],sourceCompartment,-get(k)/get(pc)/get(vSPar))
+      addValue(cSrc['targets'],sourceCompartment,-fk/fpc/fvs)
       #the target equation (where we add the current)
-      addValue(cTarg['targets'],sourceCompartment,get(k)/get(pc)/get(vTPar))
+      addValue(cTarg['targets'],sourceCompartment,fk/fpc/fvt)
         
                     
    def getDerivative(self,variable, sign, qPar, pcPar, vPar):
         
-#for flow based transfer, k=Q/vS/pc, 
-#and jacobian derivative is -Q/V/pc/pc
-#for diffusion dominated transfer, k=k1/pc
-      q=get(qPar)
-      pc=get(pcPar)
-      v=get(vPar)
+      #for flow based transfer, k=Q/vS/pc, 
+      #and jacobian derivative is -Q/V/pc/pc
+      #for diffusion dominated transfer, k=k1/pc
+      
+      tu=self.getTimeUnit()
+      q=get(tu,qPar)
+      pc=get(tu,pcPar)
+      v=get(tu,vPar)
+
       if variable=="partitionCoefficient":
          return sign*(-1)*q/v/pc/pc
       if variable=="diffusionCoefficient":
@@ -120,6 +138,9 @@ class model:
    def inspect(self):
       comps=self.compartments
       pars=self.mod['parameters']
+      
+      tu=self.getTimeUnit()
+      print('Time unit: {}'.format(tu))
       print('Compartments')
       for c in comps:
          print('{}/{}:'.format(c,self.lut[c]))
@@ -135,19 +156,19 @@ class model:
          fName=self.flows[f]
          fParName=self.mod['flows'][fName]
          fPar=pars[fParName]
-         print('\t{}[{}]:{} [{}]'.format(f,fName,fParName,get(fPar)))
+         print('\t{}[{}]:{} [{}]'.format(f,fName,fParName,get(tu,fPar)))
 
       print('Volumes')
       for v in self.mod['volumes']:
          vParName=self.mod['volumes'][v]
          vPar=pars[vParName]
-         print('\t{}:{} [{}]'.format(v,vParName,get(vPar)))
+         print('\t{}:{} [{}]'.format(v,vParName,get(tu,vPar)))
 
       print('Partition coefficients')
       for pc in self.mod['partitionCoefficients']:
          pcParName=self.mod['partitionCoefficients'][pc]
          pcPar=pars[pcParName]
-         print('\t{}:{} [{}]'.format(pc,pcParName,get(pcPar)))
+         print('\t{}:{} [{}]'.format(pc,pcParName,get(tu,pcPar)))
 
 
       print('SE parameters')
@@ -169,6 +190,7 @@ class model:
          #src=mod['sources'][s]
          self.add_source(s,self.mod['sources'][s])
       self.flows={}
+      print('timeUnit: {}'.format(self.getTimeUnit()))
       pars=self.mod['parameters']
       for f in self.mod['flows']:
          #skip comments
@@ -357,15 +379,28 @@ def addValue(qdict,compName,v):
       qdict[compName]=v
 
 
-def get(par):
+def get(timeUnit,par):
    v=par["value"]
 #convert to seconds
    try:
-      if par['unit'].split('/')[1]=='min':
-         return v/60
+      parUnits=par['unit'].split('/')
    except (KeyError,IndexError):
-      pass
-
+      #no unit given
+      return v
+   
+   try:
+      if parUnits[1]==timeUnit:
+         return v
+   except IndexError:
+      #no / in unit name
+      return v
+   if parUnits[1]=='min' and timeUnit=='s':
+      return v/60
+   
+   if parUnits[1]=='s' and timeUnit=='min':
+      return 60*v
+
+   #no idea what to do
    return v
 
 def calculateDerivative(par):

File diff suppressed because it is too large
+ 111 - 132
pythonScripts/compartmentModel.ipynb


Some files were not shown because too many files changed in this diff