Module pk2Comp
Expand source code
#!/usr/bin/env python
# coding: utf-8
# In[51]:
# Import commands
from scipy.stats import gamma
import numpy as np
import matplotlib.pyplot as plt
#%matplotlib inline
from scipy.integrate import odeint
import math as math
import os
import pandas as pd
class pk2Comp:
"""The pk2Comp object is a two compartment PK model that outputs graphs of concentration of tracer over time."""
def __init__ (self, numParam = 4, Flow = 1/60, Vp = 0.05, Visf = 0.15, PS = 1/60):
"""Initializes the model with default parameter values for flow, Vp, Visf, and PS.
Parameters
----------
numParam: int
numParam is the number of parameters you want to optimize for the model. Defaults to 4.
Flow : double
Flow is the flow of plasma through the blood vessel in mL/(mL*min). Defaults to 1/60.
Vp : double
Vp is the volume of plasma in mL. Defaults to 0.05.
Visf : double
Visf is the volume of interstitial fluid in mL. Defaults to 0.15.
PS : double
PS is the permeability-surface area constant in mL/(g*min). Defaults to 1/60.
"""
# Declare Variables for initial conditions
self.numParam = numParam
self.flow = Flow
self.Vp = Vp
self.Visf = Visf
self.PS = PS
self.sol = []
self.Cp0 = 0 # Initial concentration of tracer in plasma
self.Cisf0 = 0 # Initial concentration of tracer in interstitial space
self.tmax = 10 #Time in seconds
self.dt = 0.1 #Time step
self.a = 2. # Alpha for gamma distribution
self.rv = gamma(self.a, loc = 2, scale = 0.55) #input function
self.sol = []
self.Mass_plasma = [] #mass of tracer in plasma
self.Mass_isf = []
self.Q = []
# Define the time array
self.time = np.arange(0, self.tmax + self.dt, self.dt)
# Derivative function
def derivs(self, curr_vals, time):
"""Finds derivatives of ODEs.
Parameters
----------
curr_vals : double[]
curr_vals it he current values of the variables we wish to "update" from the curr_vals list.
time : double[]
time is our time array from 0 to tmax with timestep dt.
Returns
-------
dCp_dt : double[]
contains the derivative of concentration in plasma with respect to time.
dCisf_dt : double[]
contains the derivative of concentration in interstitial fluid with respect to time.
"""
# Unpack the current values of the variables we wish to "update" from the curr_vals list
Cp, Cisf = curr_vals
# Define value of input function Cin
Cin = self.rv.pdf(time)
# Right-hand side of odes, which are used to computer the derivative
dCp_dt = (self.flow/self.Vp)*(Cin - Cp) + (self.PS/self.Vp)*(Cisf - Cp)
dCisf_dt = (self.PS/self.Visf)*(Cp - Cisf)
return dCp_dt, dCisf_dt
def main(self):
"""Main function to solve ODEs"""
# Store the initial values in a list
init = [self.Cp0, self.Cisf0]
# Solve the odes with odeint
self.sol = odeint(self.derivs, init, self.time)
self.Mass_plasma = self.Vp * self.sol[:,0] #mass of tracer in plasma
self.Mass_isf = self.Visf * self.sol[:,1] #mass of tracer in isf
#Tp = Vp/(flow + PS) # mean transit time
#E = 1 - np.exp(-PS/flow) #extraction fraction
self.Q = self.Mass_plasma + self.Mass_isf
#print('The mean transit time is ' + str(Tp))
#print('The extraction fraction is ' + str(E))
def getPlot(self):
"""Plots the solution of the solved ODEs.
Attributes
----------
sol : double[]
contains the solutions of our ODE functions.
"""
# Plot the results using the values stored in the solution variable, "sol"
# Plot Cp using the "0" element from the solution
plt.figure(1)
plt.plot(self.time, self.rv.pdf(self.time), color = 'blue', label = 'Input Function')
plt.plot(self.time, self.sol[:,0],color="green", label = 'Cp')
# Plot Cisf using the "1" element from the solution
plt.plot(self.time, self.sol[:,1],color="purple", label = 'Cisf')
plt.xlabel('Time [s]')
plt.ylabel('Concentration [mM]')
plt.legend(loc = 'best')
plt.grid()
# Plot mass of tracer using the "2" element from the solution
plt.figure(2)
plt.plot(self.time, self.Mass_plasma,color="red", label = 'Plasma')
# Plot mass of tracer in tissue using the "3" element from the solution
plt.plot(self.time, self.Mass_isf,color="black", label = 'Interstitial Space')
plt.plot(self.time, self.Q, color="blue", label = 'Total mass')
plt.xlabel('Time [s]')
plt.ylabel('Mass [mg]')
plt.legend(loc = 'best')
plt.grid()
print('Cp at 10 sec is ' + str(self.sol[100,0]))
print('Cisf at 10 sec is ' + str(self.sol[100,1]))
# In[52]:
test = pk2Comp(4,.2,.2,.7,1/60)
test.main()
test.getPlot()
# In[ ]:
Classes
class pk2Comp (numParam=4, Flow=0.016666666666666666, Vp=0.05, Visf=0.15, PS=0.016666666666666666)
-
The pk2Comp object is a two compartment PK model that outputs graphs of concentration of tracer over time.
Initializes the model with default parameter values for flow, Vp, Visf, and PS. Parameters
numParam: int numParam is the number of parameters you want to optimize for the model. Defaults to 4.
Flow : double Flow is the flow of plasma through the blood vessel in mL/(mL*min). Defaults to 1/60.
Vp : double Vp is the volume of plasma in mL. Defaults to 0.05.
Visf : double Visf is the volume of interstitial fluid in mL. Defaults to 0.15.
PS : double PS is the permeability-surface area constant in mL/(g*min). Defaults to 1/60.
Expand source code
class pk2Comp: """The pk2Comp object is a two compartment PK model that outputs graphs of concentration of tracer over time.""" def __init__ (self, numParam = 4, Flow = 1/60, Vp = 0.05, Visf = 0.15, PS = 1/60): """Initializes the model with default parameter values for flow, Vp, Visf, and PS. Parameters ---------- numParam: int numParam is the number of parameters you want to optimize for the model. Defaults to 4. Flow : double Flow is the flow of plasma through the blood vessel in mL/(mL*min). Defaults to 1/60. Vp : double Vp is the volume of plasma in mL. Defaults to 0.05. Visf : double Visf is the volume of interstitial fluid in mL. Defaults to 0.15. PS : double PS is the permeability-surface area constant in mL/(g*min). Defaults to 1/60. """ # Declare Variables for initial conditions self.numParam = numParam self.flow = Flow self.Vp = Vp self.Visf = Visf self.PS = PS self.sol = [] self.Cp0 = 0 # Initial concentration of tracer in plasma self.Cisf0 = 0 # Initial concentration of tracer in interstitial space self.tmax = 10 #Time in seconds self.dt = 0.1 #Time step self.a = 2. # Alpha for gamma distribution self.rv = gamma(self.a, loc = 2, scale = 0.55) #input function self.sol = [] self.Mass_plasma = [] #mass of tracer in plasma self.Mass_isf = [] self.Q = [] # Define the time array self.time = np.arange(0, self.tmax + self.dt, self.dt) # Derivative function def derivs(self, curr_vals, time): """Finds derivatives of ODEs. Parameters ---------- curr_vals : double[] curr_vals it he current values of the variables we wish to "update" from the curr_vals list. time : double[] time is our time array from 0 to tmax with timestep dt. Returns ------- dCp_dt : double[] contains the derivative of concentration in plasma with respect to time. dCisf_dt : double[] contains the derivative of concentration in interstitial fluid with respect to time. """ # Unpack the current values of the variables we wish to "update" from the curr_vals list Cp, Cisf = curr_vals # Define value of input function Cin Cin = self.rv.pdf(time) # Right-hand side of odes, which are used to computer the derivative dCp_dt = (self.flow/self.Vp)*(Cin - Cp) + (self.PS/self.Vp)*(Cisf - Cp) dCisf_dt = (self.PS/self.Visf)*(Cp - Cisf) return dCp_dt, dCisf_dt def main(self): """Main function to solve ODEs""" # Store the initial values in a list init = [self.Cp0, self.Cisf0] # Solve the odes with odeint self.sol = odeint(self.derivs, init, self.time) self.Mass_plasma = self.Vp * self.sol[:,0] #mass of tracer in plasma self.Mass_isf = self.Visf * self.sol[:,1] #mass of tracer in isf #Tp = Vp/(flow + PS) # mean transit time #E = 1 - np.exp(-PS/flow) #extraction fraction self.Q = self.Mass_plasma + self.Mass_isf #print('The mean transit time is ' + str(Tp)) #print('The extraction fraction is ' + str(E)) def getPlot(self): """Plots the solution of the solved ODEs. Attributes ---------- sol : double[] contains the solutions of our ODE functions. """ # Plot the results using the values stored in the solution variable, "sol" # Plot Cp using the "0" element from the solution plt.figure(1) plt.plot(self.time, self.rv.pdf(self.time), color = 'blue', label = 'Input Function') plt.plot(self.time, self.sol[:,0],color="green", label = 'Cp') # Plot Cisf using the "1" element from the solution plt.plot(self.time, self.sol[:,1],color="purple", label = 'Cisf') plt.xlabel('Time [s]') plt.ylabel('Concentration [mM]') plt.legend(loc = 'best') plt.grid() # Plot mass of tracer using the "2" element from the solution plt.figure(2) plt.plot(self.time, self.Mass_plasma,color="red", label = 'Plasma') # Plot mass of tracer in tissue using the "3" element from the solution plt.plot(self.time, self.Mass_isf,color="black", label = 'Interstitial Space') plt.plot(self.time, self.Q, color="blue", label = 'Total mass') plt.xlabel('Time [s]') plt.ylabel('Mass [mg]') plt.legend(loc = 'best') plt.grid() print('Cp at 10 sec is ' + str(self.sol[100,0])) print('Cisf at 10 sec is ' + str(self.sol[100,1]))
Methods
def derivs(self, curr_vals, time)
-
Finds derivatives of ODEs.
Parameters
curr_vals : double[] curr_vals it he current values of the variables we wish to "update" from the curr_vals list.
time : double[] time is our time array from 0 to tmax with timestep dt.
Returns
dCp_dt
:double
[]- contains the derivative of concentration in plasma with respect to time.
dCisf_dt
:double
[]- contains the derivative of concentration in interstitial fluid with respect to time.
Expand source code
def derivs(self, curr_vals, time): """Finds derivatives of ODEs. Parameters ---------- curr_vals : double[] curr_vals it he current values of the variables we wish to "update" from the curr_vals list. time : double[] time is our time array from 0 to tmax with timestep dt. Returns ------- dCp_dt : double[] contains the derivative of concentration in plasma with respect to time. dCisf_dt : double[] contains the derivative of concentration in interstitial fluid with respect to time. """ # Unpack the current values of the variables we wish to "update" from the curr_vals list Cp, Cisf = curr_vals # Define value of input function Cin Cin = self.rv.pdf(time) # Right-hand side of odes, which are used to computer the derivative dCp_dt = (self.flow/self.Vp)*(Cin - Cp) + (self.PS/self.Vp)*(Cisf - Cp) dCisf_dt = (self.PS/self.Visf)*(Cp - Cisf) return dCp_dt, dCisf_dt
def getPlot(self)
-
Plots the solution of the solved ODEs.
Attributes
sol : double[] contains the solutions of our ODE functions.
Expand source code
def getPlot(self): """Plots the solution of the solved ODEs. Attributes ---------- sol : double[] contains the solutions of our ODE functions. """ # Plot the results using the values stored in the solution variable, "sol" # Plot Cp using the "0" element from the solution plt.figure(1) plt.plot(self.time, self.rv.pdf(self.time), color = 'blue', label = 'Input Function') plt.plot(self.time, self.sol[:,0],color="green", label = 'Cp') # Plot Cisf using the "1" element from the solution plt.plot(self.time, self.sol[:,1],color="purple", label = 'Cisf') plt.xlabel('Time [s]') plt.ylabel('Concentration [mM]') plt.legend(loc = 'best') plt.grid() # Plot mass of tracer using the "2" element from the solution plt.figure(2) plt.plot(self.time, self.Mass_plasma,color="red", label = 'Plasma') # Plot mass of tracer in tissue using the "3" element from the solution plt.plot(self.time, self.Mass_isf,color="black", label = 'Interstitial Space') plt.plot(self.time, self.Q, color="blue", label = 'Total mass') plt.xlabel('Time [s]') plt.ylabel('Mass [mg]') plt.legend(loc = 'best') plt.grid() print('Cp at 10 sec is ' + str(self.sol[100,0])) print('Cisf at 10 sec is ' + str(self.sol[100,1]))
def main(self)
-
Main function to solve ODEs
Expand source code
def main(self): """Main function to solve ODEs""" # Store the initial values in a list init = [self.Cp0, self.Cisf0] # Solve the odes with odeint self.sol = odeint(self.derivs, init, self.time) self.Mass_plasma = self.Vp * self.sol[:,0] #mass of tracer in plasma self.Mass_isf = self.Visf * self.sol[:,1] #mass of tracer in isf #Tp = Vp/(flow + PS) # mean transit time #E = 1 - np.exp(-PS/flow) #extraction fraction self.Q = self.Mass_plasma + self.Mass_isf