Module pkOptimizer

Expand source code
#!/usr/bin/env python
# coding: utf-8

# In[9]:


from scipy.stats import gamma
from scipy.integrate import odeint 
from scipy.optimize import minimize
from scipy.optimize import curve_fit

import os
import csv
import re
import math as math
import numpy as np
import matplotlib.pyplot as plt
#%matplotlib inline

class pkOptimizer:
    """The pkOptimizer object is an optimizer for parameters in pk models."""
    
    def __init__ (self, wd, Flow = 1/60, Vp = 0.05, Visf = 0.15, PS = 1/60):
        """Initializes the model with initial guess parameter values for flow, Vp, Visf, and PS.
        Parameters
        ----------      
        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.    
        """
        
    def getData(self, wd):
        """Imports data from all .csv files in directory.
        Parameters
        ----------  
        wd : str
            wd is the working directory path
            
        Attributes
        ----------
        t : double[]
            list of all timepoints
        aorta : double[]
            concentration of tracer in aorta (input function)
        myo : double[]
            concentration of tracer in myocardial tissue (Cisf)
        
        Returns
        -------
        t : double[]
            list of all timepoints
        aorta : double[]
            concentration of tracer in aorta (input function)
        myo : double[]
            concentration of tracer in myocardial tissue (Cisf)
        """
    
        os.chdir(wd)
        #os.chdir(r"C:\Users\Ethan\OneDrive - Michigan State University\MSU\Classwork\Computational Modeling\Models\Data")
        #create directory of all csv files,
        data = list(csv.reader(open('CTPERF005_stress.csv'), delimiter = '\t'))

        t = []
        aorta = []
        myo = []
        
        for i in range(12):
            t.append(float(re.compile('\d+[.]+\d+|\d+').findall(data[i+1][0])[0]))
            aorta.append(float(re.compile('\d+[.]+\d+|\d+').findall(data[i+1][1])[0]))
            myo.append(float(re.compile('\d+[.]+\d+|\d+').findall(data[i+1][2])[0]))

        return t, aorta, myo

    def gammaFunc(self, time, a, l, s):
        """Creates a gamma variate probability density function with given alpha, location, and scale values.
        Parameters
        ----------  
        time : double[]
            array of timepoints
        a : double
            alpha value of gamma PDF
        l : double
            location of 50th percentile of function
        s : double
            scale parameter 
            
        Returns
        -------
        rv.pdf(time)
            probability density function of your gamma variate.
        """
        rv = gamma(a, loc = l, scale = s) #input function
        return rv.pdf(time)
    
    def curveFit(self, t, aorta, myo, model):
        """Takes in data and fits gamma curve to aorta and Cisf from model to myo. Returns parameters for best fit.
        
        Parameters
        ----------  
        t : double[]
            list of all timepoints
        aorta : double[]
            concentration of tracer in aorta (input function)
        myo : double[]
            concentration of tracer in myocardial tissue (Cisf)
        model : pkModel object
            a pk model, either 1Comp or 2Comp
        
        Returns
        -------
        Flow : double
            Flow is the flow of plasma through the blood vessel in mL/(mL*min).
        
        Vp : double
            Vp is the volume of plasma in mL.
            
        Visf : double
            Visf is the volume of interstitial fluid in mL.
        
        PS : double
            PS is the permeability-surface area constant in mL/(g*min).
        """
    
    def getPlot(self):
        """Plots the original data to the fitted curve."""
        plt.plot(t, aorta, 'bo', label='data')
        #plt.plot(t, y, 'b-', label='data')
        popt, pcov = curve_fit(gammaFunc, t, aorta, p0 = [2, 8, 10000], method = 'trf')

        print(f'alpha = {popt[0]}, loc = {popt[1]}, scale = {popt[2]}')

        plt.plot(t, gammaFunc(t, *popt), 'r-', label='fit: a=%5.3f, b=%5.3f, c=%5.3f' % tuple(popt))
        plt.plot(time, gammaFunc(time, .1313, 8.533, 10000), 'b-')


# In[ ]:

Classes

class pkOptimizer (wd, Flow=0.016666666666666666, Vp=0.05, Visf=0.15, PS=0.016666666666666666)

The pkOptimizer object is an optimizer for parameters in pk models.

Initializes the model with initial guess parameter values for flow, Vp, Visf, and PS. Parameters


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 pkOptimizer:
    """The pkOptimizer object is an optimizer for parameters in pk models."""
    
    def __init__ (self, wd, Flow = 1/60, Vp = 0.05, Visf = 0.15, PS = 1/60):
        """Initializes the model with initial guess parameter values for flow, Vp, Visf, and PS.
        Parameters
        ----------      
        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.    
        """
        
    def getData(self, wd):
        """Imports data from all .csv files in directory.
        Parameters
        ----------  
        wd : str
            wd is the working directory path
            
        Attributes
        ----------
        t : double[]
            list of all timepoints
        aorta : double[]
            concentration of tracer in aorta (input function)
        myo : double[]
            concentration of tracer in myocardial tissue (Cisf)
        
        Returns
        -------
        t : double[]
            list of all timepoints
        aorta : double[]
            concentration of tracer in aorta (input function)
        myo : double[]
            concentration of tracer in myocardial tissue (Cisf)
        """
    
        os.chdir(wd)
        #os.chdir(r"C:\Users\Ethan\OneDrive - Michigan State University\MSU\Classwork\Computational Modeling\Models\Data")
        #create directory of all csv files,
        data = list(csv.reader(open('CTPERF005_stress.csv'), delimiter = '\t'))

        t = []
        aorta = []
        myo = []
        
        for i in range(12):
            t.append(float(re.compile('\d+[.]+\d+|\d+').findall(data[i+1][0])[0]))
            aorta.append(float(re.compile('\d+[.]+\d+|\d+').findall(data[i+1][1])[0]))
            myo.append(float(re.compile('\d+[.]+\d+|\d+').findall(data[i+1][2])[0]))

        return t, aorta, myo

    def gammaFunc(self, time, a, l, s):
        """Creates a gamma variate probability density function with given alpha, location, and scale values.
        Parameters
        ----------  
        time : double[]
            array of timepoints
        a : double
            alpha value of gamma PDF
        l : double
            location of 50th percentile of function
        s : double
            scale parameter 
            
        Returns
        -------
        rv.pdf(time)
            probability density function of your gamma variate.
        """
        rv = gamma(a, loc = l, scale = s) #input function
        return rv.pdf(time)
    
    def curveFit(self, t, aorta, myo, model):
        """Takes in data and fits gamma curve to aorta and Cisf from model to myo. Returns parameters for best fit.
        
        Parameters
        ----------  
        t : double[]
            list of all timepoints
        aorta : double[]
            concentration of tracer in aorta (input function)
        myo : double[]
            concentration of tracer in myocardial tissue (Cisf)
        model : pkModel object
            a pk model, either 1Comp or 2Comp
        
        Returns
        -------
        Flow : double
            Flow is the flow of plasma through the blood vessel in mL/(mL*min).
        
        Vp : double
            Vp is the volume of plasma in mL.
            
        Visf : double
            Visf is the volume of interstitial fluid in mL.
        
        PS : double
            PS is the permeability-surface area constant in mL/(g*min).
        """
    
    def getPlot(self):
        """Plots the original data to the fitted curve."""
        plt.plot(t, aorta, 'bo', label='data')
        #plt.plot(t, y, 'b-', label='data')
        popt, pcov = curve_fit(gammaFunc, t, aorta, p0 = [2, 8, 10000], method = 'trf')

        print(f'alpha = {popt[0]}, loc = {popt[1]}, scale = {popt[2]}')

        plt.plot(t, gammaFunc(t, *popt), 'r-', label='fit: a=%5.3f, b=%5.3f, c=%5.3f' % tuple(popt))
        plt.plot(time, gammaFunc(time, .1313, 8.533, 10000), 'b-')

Methods

def curveFit(self, t, aorta, myo, model)

Takes in data and fits gamma curve to aorta and Cisf from model to myo. Returns parameters for best fit.

Parameters

t : double[] list of all timepoints aorta : double[] concentration of tracer in aorta (input function) myo : double[] concentration of tracer in myocardial tissue (Cisf) model : pkModel object a pk model, either 1Comp or 2Comp

Returns

Flow : double
Flow is the flow of plasma through the blood vessel in mL/(mL*min).
Vp : double
Vp is the volume of plasma in mL.
Visf : double
Visf is the volume of interstitial fluid in mL.
PS : double
PS is the permeability-surface area constant in mL/(g*min).
Expand source code
def curveFit(self, t, aorta, myo, model):
    """Takes in data and fits gamma curve to aorta and Cisf from model to myo. Returns parameters for best fit.
    
    Parameters
    ----------  
    t : double[]
        list of all timepoints
    aorta : double[]
        concentration of tracer in aorta (input function)
    myo : double[]
        concentration of tracer in myocardial tissue (Cisf)
    model : pkModel object
        a pk model, either 1Comp or 2Comp
    
    Returns
    -------
    Flow : double
        Flow is the flow of plasma through the blood vessel in mL/(mL*min).
    
    Vp : double
        Vp is the volume of plasma in mL.
        
    Visf : double
        Visf is the volume of interstitial fluid in mL.
    
    PS : double
        PS is the permeability-surface area constant in mL/(g*min).
    """
def gammaFunc(self, time, a, l, s)

Creates a gamma variate probability density function with given alpha, location, and scale values. Parameters


time : double[] array of timepoints a : double alpha value of gamma PDF l : double location of 50th percentile of function s : double scale parameter

Returns

rv.pdf(time)
probability density function of your gamma variate.
Expand source code
def gammaFunc(self, time, a, l, s):
    """Creates a gamma variate probability density function with given alpha, location, and scale values.
    Parameters
    ----------  
    time : double[]
        array of timepoints
    a : double
        alpha value of gamma PDF
    l : double
        location of 50th percentile of function
    s : double
        scale parameter 
        
    Returns
    -------
    rv.pdf(time)
        probability density function of your gamma variate.
    """
    rv = gamma(a, loc = l, scale = s) #input function
    return rv.pdf(time)
def getData(self, wd)

Imports data from all .csv files in directory. Parameters


wd : str wd is the working directory path

Attributes

t : double[]
list of all timepoints
aorta : double[]
concentration of tracer in aorta (input function)
myo : double[]
concentration of tracer in myocardial tissue (Cisf)

Returns

t : double[]
list of all timepoints
aorta : double[]
concentration of tracer in aorta (input function)
myo : double[]
concentration of tracer in myocardial tissue (Cisf)
Expand source code
def getData(self, wd):
    """Imports data from all .csv files in directory.
    Parameters
    ----------  
    wd : str
        wd is the working directory path
        
    Attributes
    ----------
    t : double[]
        list of all timepoints
    aorta : double[]
        concentration of tracer in aorta (input function)
    myo : double[]
        concentration of tracer in myocardial tissue (Cisf)
    
    Returns
    -------
    t : double[]
        list of all timepoints
    aorta : double[]
        concentration of tracer in aorta (input function)
    myo : double[]
        concentration of tracer in myocardial tissue (Cisf)
    """

    os.chdir(wd)
    #os.chdir(r"C:\Users\Ethan\OneDrive - Michigan State University\MSU\Classwork\Computational Modeling\Models\Data")
    #create directory of all csv files,
    data = list(csv.reader(open('CTPERF005_stress.csv'), delimiter = '\t'))

    t = []
    aorta = []
    myo = []
    
    for i in range(12):
        t.append(float(re.compile('\d+[.]+\d+|\d+').findall(data[i+1][0])[0]))
        aorta.append(float(re.compile('\d+[.]+\d+|\d+').findall(data[i+1][1])[0]))
        myo.append(float(re.compile('\d+[.]+\d+|\d+').findall(data[i+1][2])[0]))

    return t, aorta, myo
def getPlot(self)

Plots the original data to the fitted curve.

Expand source code
def getPlot(self):
    """Plots the original data to the fitted curve."""
    plt.plot(t, aorta, 'bo', label='data')
    #plt.plot(t, y, 'b-', label='data')
    popt, pcov = curve_fit(gammaFunc, t, aorta, p0 = [2, 8, 10000], method = 'trf')

    print(f'alpha = {popt[0]}, loc = {popt[1]}, scale = {popt[2]}')

    plt.plot(t, gammaFunc(t, *popt), 'r-', label='fit: a=%5.3f, b=%5.3f, c=%5.3f' % tuple(popt))
    plt.plot(time, gammaFunc(time, .1313, 8.533, 10000), 'b-')