diff --git a/docs/pk_two_comp.html b/docs/pk_two_comp.html
new file mode 100644
index 0000000000000000000000000000000000000000..a7c4dbbe156bf2bd1d449fc445c47edd96923673
--- /dev/null
+++ b/docs/pk_two_comp.html
@@ -0,0 +1,1506 @@
+<!doctype html>
+<html lang="en">
+<head>
+<meta charset="utf-8">
+<meta name="viewport" content="width=device-width, initial-scale=1, minimum-scale=1" />
+<meta name="generator" content="pdoc 0.7.4" />
+<title>pk_two_comp API documentation</title>
+<meta name="description" content="The pk2Comp object is a two compartment PK model
+that outputs graphs of concentration of tracer over time." />
+<link href='https://cdnjs.cloudflare.com/ajax/libs/normalize/8.0.0/normalize.min.css' rel='stylesheet'>
+<link href='https://cdnjs.cloudflare.com/ajax/libs/10up-sanitize.css/8.0.0/sanitize.min.css' rel='stylesheet'>
+<link href="https://cdnjs.cloudflare.com/ajax/libs/highlight.js/9.12.0/styles/github.min.css" rel="stylesheet">
+<style>.flex{display:flex !important}body{line-height:1.5em}#content{padding:20px}#sidebar{padding:30px;overflow:hidden}.http-server-breadcrumbs{font-size:130%;margin:0 0 15px 0}#footer{font-size:.75em;padding:5px 30px;border-top:1px solid #ddd;text-align:right}#footer p{margin:0 0 0 1em;display:inline-block}#footer p:last-child{margin-right:30px}h1,h2,h3,h4,h5{font-weight:300}h1{font-size:2.5em;line-height:1.1em}h2{font-size:1.75em;margin:1em 0 .50em 0}h3{font-size:1.4em;margin:25px 0 10px 0}h4{margin:0;font-size:105%}a{color:#058;text-decoration:none;transition:color .3s ease-in-out}a:hover{color:#e82}.title code{font-weight:bold}h2[id^="header-"]{margin-top:2em}.ident{color:#900}pre code{background:#f8f8f8;font-size:.8em;line-height:1.4em}code{background:#f2f2f1;padding:1px 4px;overflow-wrap:break-word}h1 code{background:transparent}pre{background:#f8f8f8;border:0;border-top:1px solid #ccc;border-bottom:1px solid #ccc;margin:1em 0;padding:1ex}#http-server-module-list{display:flex;flex-flow:column}#http-server-module-list div{display:flex}#http-server-module-list dt{min-width:10%}#http-server-module-list p{margin-top:0}.toc ul,#index{list-style-type:none;margin:0;padding:0}#index code{background:transparent}#index h3{border-bottom:1px solid #ddd}#index ul{padding:0}#index h4{font-weight:bold}#index h4 + ul{margin-bottom:.6em}@media (min-width:200ex){#index .two-column{column-count:2}}@media (min-width:300ex){#index .two-column{column-count:3}}dl{margin-bottom:2em}dl dl:last-child{margin-bottom:4em}dd{margin:0 0 1em 3em}#header-classes + dl > dd{margin-bottom:3em}dd dd{margin-left:2em}dd p{margin:10px 0}.name{background:#eee;font-weight:bold;font-size:.85em;padding:5px 10px;display:inline-block;min-width:40%}.name:hover{background:#e0e0e0}.name > span:first-child{white-space:nowrap}.name.class > span:nth-child(2){margin-left:.4em}.inherited{color:#999;border-left:5px solid #eee;padding-left:1em}.inheritance em{font-style:normal;font-weight:bold}.desc h2{font-weight:400;font-size:1.25em}.desc h3{font-size:1em}.desc dt code{background:inherit}.source summary,.git-link-div{color:#666;text-align:right;font-weight:400;font-size:.8em;text-transform:uppercase}.source summary > *{white-space:nowrap;cursor:pointer}.git-link{color:inherit;margin-left:1em}.source pre{max-height:500px;overflow:auto;margin:0}.source pre code{font-size:12px;overflow:visible}.hlist{list-style:none}.hlist li{display:inline}.hlist li:after{content:',\2002'}.hlist li:last-child:after{content:none}.hlist .hlist{display:inline;padding-left:1em}img{max-width:100%}.admonition{padding:.1em .5em;margin-bottom:1em}.admonition-title{font-weight:bold}.admonition.note,.admonition.info,.admonition.important{background:#aef}.admonition.todo,.admonition.versionadded,.admonition.tip,.admonition.hint{background:#dfd}.admonition.warning,.admonition.versionchanged,.admonition.deprecated{background:#fd4}.admonition.error,.admonition.danger,.admonition.caution{background:lightpink}</style>
+<style media="screen and (min-width: 700px)">@media screen and (min-width:700px){#sidebar{width:30%}#content{width:70%;max-width:100ch;padding:3em 4em;border-left:1px solid #ddd}pre code{font-size:1em}.item .name{font-size:1em}main{display:flex;flex-direction:row-reverse;justify-content:flex-end}.toc ul ul,#index ul{padding-left:1.5em}.toc > ul > li{margin-top:.5em}}</style>
+<style media="print">@media print{#sidebar h1{page-break-before:always}.source{display:none}}@media print{*{background:transparent !important;color:#000 !important;box-shadow:none !important;text-shadow:none !important}a[href]:after{content:" (" attr(href) ")";font-size:90%}a[href][title]:after{content:none}abbr[title]:after{content:" (" attr(title) ")"}.ir a:after,a[href^="javascript:"]:after,a[href^="#"]:after{content:""}pre,blockquote{border:1px solid #999;page-break-inside:avoid}thead{display:table-header-group}tr,img{page-break-inside:avoid}img{max-width:100% !important}@page{margin:0.5cm}p,h2,h3{orphans:3;widows:3}h1,h2,h3,h4,h5,h6{page-break-after:avoid}}</style>
+</head>
+<body>
+<main>
+<article id="content">
+<header>
+<h1 class="title">Module <code>pk_two_comp</code></h1>
+</header>
+<section id="section-intro">
+<p>The pk2Comp object is a two compartment PK model
+that outputs graphs of concentration of tracer over time.</p>
+<details class="source">
+<summary>
+<span>Expand source code</span>
+</summary>
+<pre><code class="python">#!/usr/bin/env python
+# coding: utf-8
+
+# In[1]:
+
+
+"""The pk2Comp object is a two compartment PK model
+ that outputs graphs of concentration of tracer over time.
+"""
+#!/usr/bin/env python
+# coding: utf-8
+
+# In[1]:
+import pathlib
+import os
+import csv
+import re
+import math
+import numpy as np
+import matplotlib.pyplot as plt
+
+from scipy.integrate import solve_ivp
+from scipy.optimize import fmin
+
+class pk_two_comp:
+ """The pk2Comp object is a two compartment PK model
+ that outputs graphs of concentration of tracer over time.
+ """
+
+ def __init__(self, wd=pathlib.Path('Data').absolute(), filename='CTPERF005_stress.csv'):
+ """Initializes the model with default parameter values for flow, Vp, Visf, and PS.
+ Parameters
+ ----------
+ wd : path
+ Absolute path name to current directory. Defaults to ./Data
+ filename : String
+ Name of the data file you'd like to access
+
+ Attributes
+ -----------
+ time : double[]
+ list of all timepoints
+ aorta : double[]
+ concentration of tracer in aorta (input function)
+ myo : double[]
+ concentration of tracer in myocardial tissue (conc_isf)
+ 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.
+ ymax : int
+ Magnitude of Gamma-var peak.
+ tmax : double
+ Time of which highest peak in Gamma-var appears
+ alpha : double
+ Scale factor
+ delay : double
+ Delay to start Gamma-var curve.
+
+ """
+ # Subject Data
+ if os.path.basename(os.path.normpath(pathlib.Path().absolute())) != 'Data':
+ self.wd = pathlib.Path('Data').absolute()
+ else:
+ self.wd = wd
+
+ if not isinstance(filename, str):
+ raise ValueError('Filename must be a string')
+
+ self.filename = filename
+ self.time = []
+ self.aorta = []
+ self.myo = []
+
+ # Declare Variables for initial conditions
+ # Compartment variables to be fitted
+ self.flow = 1/60
+ self.visf = 0.15
+ self.baseline = 60
+
+ # Other Compartmental Modelvariables
+ self.perm_surf = 0.35
+ self.vol_plasma = 0.10
+
+ # Solved ode
+ self.sol = []
+
+ # Gamma variables
+ self.ymax = 250
+ self.tmax = 6.5
+ self.alpha = 2.5
+ self.delay = 0
+
+ self.deriv_sol = np.array([])
+ self.fit_myo = np.array([])
+
+ def get_data(self, filename):
+ """Imports data from all .csv files in directory.
+ Parameters
+ ----------
+ filename : string
+ Name of the data file you'd like to access
+ wd : str
+ wd is the working directory path
+
+ Attributes
+ ----------
+ time : double[]
+ list of all timepoints
+ aorta : double[]
+ concentration of tracer in aorta (input function)
+ myo : double[]
+ concentration of tracer in myocardial tissue (conc_isf)
+
+ Returns
+ -------
+ time : double[]
+ list of all timepoints
+ aorta : double[]
+ concentration of tracer in aorta (input function)
+ myo : double[]
+ concentration of tracer in myocardial tissue (conc_isf)
+ """
+ self.time = []
+ self.aorta = []
+ self.myo = []
+
+ os.chdir(self.wd)
+ # File not found error
+ if not os.path.isfile(filename):
+ raise ValueError(
+ "Input file does not exist: {0}. I'll quit now.".format(filename))
+
+ data = list(csv.reader(open(filename), delimiter='\t'))
+
+ for i in range(12):
+ self.time.append(
+ float(re.compile(r'\d+[.]+\d+|\d+').findall(data[i+1][0])[0]))
+ self.aorta.append(
+ float(re.compile(r'\d+[.]+\d+|\d+').findall(data[i+1][1])[0]))
+ self.myo.append(
+ float(re.compile(r'\d+[.]+\d+|\d+').findall(data[i+1][2])[0]))
+
+ return self.time, self.aorta, self.myo
+
+ # gamma_var distribution curve
+ def gamma_var(self, time=np.arange(0, 25), ymax=10, tmax=10, alpha=2, delay=0):
+ """Creates a gamma variate probability density function with given alpha,
+ location, and scale values.
+ Parameters
+ ----------
+ t : double[]
+ array of timepoints
+ ymax : double
+ peak y value of gamma distribution
+ tmax : double
+ location of 50th percentile of function
+ alpha : double
+ scale parameter
+ delay : double
+ time delay of which to start gamma distribution
+
+ Returns
+ -------
+ y : double[]
+ probability density function of your gamma variate.
+ """
+ # Following Madsen 1992 simplified parameterization for gamma variate
+ t = time
+ self.ymax = ymax
+ self.tmax = tmax
+ self.alpha = alpha
+ self.delay = delay
+
+ y = np.zeros(np.size(t)) # preallocate output
+
+ # For odeint, checks if t input is array or float
+ if isinstance(t, (list, np.ndarray)):
+ for i in range(np.size(y)):
+ if t[i] < delay:
+ y[i] = 0
+ else:
+ y[i] = round((ymax*tmax**(-alpha)*math.exp(alpha))*(t[i]-delay)
+ ** alpha*math.exp(-alpha*(t[i]-delay)/tmax), 3)
+ else:
+ y = round((ymax*tmax**(-alpha)*math.exp(alpha))*(t-delay)
+ ** alpha*math.exp(-alpha*(t-delay)/tmax), 3)
+ return y
+
+ # gamma_var_error
+ def input_mse(self, guess=[10, 10, 2, 5]):
+ """Calculates Mean squared error (MSE) between data and
+ gamma variate with given parameters values.
+ Parameters
+ ----------
+ param : ndarray[]
+ time : double[]
+ array of timepoints
+ ymax : double
+ peak y value of gamma distribution
+ tmax : double
+ location of 50th percentile of function
+ alpha : double
+ scale parameter
+ delay : double
+ time delay of which to start gamma distribution
+
+ Returns
+ -------
+ MSE : double
+ Mean squared error
+ """
+ if len(guess) < 1:
+ self.ymax = 10
+ self.tmax = 10
+ self.alpha = 2
+ self.delay = 5
+ elif len(guess) < 2:
+ self.ymax = guess[0]
+ self.tmax = 10
+ self.alpha = 2
+ self.delay = 5
+ elif len(guess) < 3:
+ self.ymax = guess[0]
+ self.tmax = guess[1]
+ self.alpha = 2
+ self.delay = 5
+ elif len(guess) < 4:
+ self.ymax = guess[0]
+ self.tmax = guess[1]
+ self.alpha = guess[2]
+ self.delay = 5
+ else:
+ # Mean squared error (MSE) between data and gamma variate with given parameters
+ self.ymax = guess[0]
+ self.tmax = guess[1]
+ self.alpha = guess[2]
+ self.delay = guess[3]
+
+ mse = 0
+
+ if self.tmax <= 0 or self.ymax <= 10 or self.delay < 0 or self.alpha < 0 \
+ or self.alpha > 1000 or self.tmax > 1000:
+ mse = 1000000 # just return a big number
+
+ else:
+ model_vals = self.gamma_var(
+ self.time, self.ymax, self.tmax, self.alpha, self.delay)
+
+ for i in range(len(self.aorta)):
+ mse = (self.aorta[i] - model_vals[i])**2 + mse
+ mse = mse / len(self.aorta)
+ return round(mse, 3)
+
+ def input_func_fit(self, initGuesses):
+ """Uses fmin algorithm (Nelder-Mead simplex algorithm) to
+ minimize loss function (MSE) of input function.
+ Parameters
+ ----------
+ initGuesses : ndarray[]
+ Array of initial guesses containing:
+ time : double[]
+ array of timepoints
+ ymax : double
+ peak y value of gamma distribution
+ tmax : double
+ location of 50th percentile of function
+ alpha : double
+ scale parameter
+ delay : double
+ time delay of which to start gamma distribution
+ Returns
+ -------
+ opt : double[]
+ optimized parameters
+ """
+ # Mean squared error (MSE) between data and gamma variate with given parameters
+ opt = fmin(self.input_mse, initGuesses, maxiter=1000)
+
+ self.ymax = opt[0]
+ self.tmax = opt[1]
+ self.alpha = opt[2]
+ self.delay = opt[3]
+
+ return opt.round(2)
+
+ # Derivative function
+ def derivs(self, time, curr_vals):
+ """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
+ -------
+ dconc_plasma_dt : double[]
+ contains the derivative of concentration in plasma with respect to time.
+ dconc_isf_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
+ conc_plasma, conc_isf = curr_vals
+
+ # Define value of input function conc_in
+ conc_in = self.gamma_var(time, self.ymax, self.tmax,
+ self.alpha, self.delay)
+
+ # Right-hand side of odes, which are used to computer the derivative
+ dconc_plasma_dt = (self.flow/self.vol_plasma)*(conc_in - conc_plasma) + \
+ (self.perm_surf/self.vol_plasma)*(conc_isf - conc_plasma)
+ dconc_isf_dt = (self.perm_surf/self.visf)*(conc_plasma - conc_isf)
+ return dconc_plasma_dt, dconc_isf_dt
+
+ def output_mse(self, guess):
+ """Calculates Mean squared error (MSE) between data and
+ output derivs with given parameters values.
+ Parameters
+ ----------
+ guess : ndarray[]
+
+ 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.
+ Returns
+ -------
+ MSE : double
+ Mean squared error
+ """
+ self.flow = guess[0]
+ self.visf = guess[1]
+ self.baseline = guess[2]
+
+ mse = 0
+
+ if self.flow <= 0 or self.flow >= 25 or self.visf > 100 \
+ or self.visf < 0 or self.baseline > 150 or self.baseline < 0:
+ mse = 100000 # just return a big number
+
+ else:
+ sol = solve_ivp(self.derivs, [0, 30], [0, 0], t_eval=self.time)
+ MBF = sol.y[0] + sol.y[1]
+
+ temp = np.asarray(self.myo) - self.baseline
+
+ for i in range(len(self.myo)):
+ mse = (temp[i] - MBF[i])**2 + mse
+
+ mse = mse / len(self.myo)
+ return mse
+
+ def output_func_fit(self, initGuesses):
+ """Uses fmin algorithm (Nelder-Mead simplex algorithm) to minimize
+ loss function (MSE) of output function.
+ Parameters
+ ----------
+ initGuesses : ndarray[]
+ Array of initial guesses containing:
+ 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.
+ Returns
+ -------
+ opt : double[]
+ optimized parameters
+ """
+ # Mean squared error (MSE) between data and gamma variate with given parameters
+ opt1 = fmin(self.output_mse, initGuesses, maxiter=10000)
+
+ self.flow = opt1[0].round(4)
+ self.visf = opt1[1].round(4)
+ self.baseline = opt1[2].round(4)
+
+ return opt1.round(4)
+
+ def main(self):
+
+ # Gets data from file
+ self.get_data(self.filename)
+
+ # Plots original data
+ plt.plot(self.time, self.aorta, 'bo')
+ plt.plot(self.time, self.myo, 'ro')
+
+ # Fit gamma_var input function and plots it
+ opt = self.input_func_fit([250, 7, 4, 0])
+ plt.plot(np.arange(0, 25, 0.01), self.gamma_var(np.arange(0, 25, 0.01),
+ opt[0], opt[1], opt[2], opt[3]), 'k-')
+
+ # Fit uptake function and plot it
+ opt2 = self.output_func_fit([.011418, .62, self.myo[0]])
+ self.deriv_sol = solve_ivp(self.derivs, [0, 30],
+ [0, 0], t_eval=self.time)
+ self.fit_myo = self.deriv_sol.y[0] + self.deriv_sol.y[1]
+ plt.plot(self.time, self.fit_myo + self.baseline, 'm-')
+ </code></pre>
+</details>
+</section>
+<section>
+</section>
+<section>
+</section>
+<section>
+</section>
+<section>
+<h2 class="section-title" id="header-classes">Classes</h2>
+<dl>
+<dt id="pk_two_comp.pk_two_comp"><code class="flex name class">
+<span>class <span class="ident">pk_two_comp</span></span>
+<span>(</span><span>wd=WindowsPath('C:/Users/Ethan/OneDrive - Michigan State University/MSU/Classwork/Computational Modeling/Models/Data'), filename='CTPERF005_stress.csv')</span>
+</code></dt>
+<dd>
+<section class="desc"><p>The pk2Comp object is a two compartment PK model
+that outputs graphs of concentration of tracer over time.</p>
+<p>Initializes the model with default parameter values for flow, Vp, Visf, and PS.
+Parameters</p>
+<hr>
+<dl>
+<dt><strong><code>wd</code></strong> : <code>path</code></dt>
+<dd>Absolute path name to current directory. Defaults to ./Data</dd>
+<dt><strong><code>filename</code></strong> : <code>String</code></dt>
+<dd>Name of the data file you'd like to access</dd>
+</dl>
+<h2 id="attributes">Attributes</h2>
+<dl>
+<dt><strong><code>time</code></strong> : <code>double</code>[]</dt>
+<dd>list of all timepoints</dd>
+<dt><strong><code>aorta</code></strong> : <code>double</code>[]</dt>
+<dd>concentration of tracer in aorta (input function)</dd>
+<dt><strong><code>myo</code></strong> : <code>double</code>[]</dt>
+<dd>concentration of tracer in myocardial tissue (conc_isf)</dd>
+<dt><strong><code>Flow</code></strong> : <code>double</code></dt>
+<dd>Flow is the flow of plasma through the blood vessel in mL/(mL*min). Defaults to 1/60.</dd>
+<dt><strong><code>Vp</code></strong> : <code>double</code></dt>
+<dd>Vp is the volume of plasma in mL. Defaults to 0.05.</dd>
+<dt><strong><code>Visf</code></strong> : <code>double</code></dt>
+<dd>Visf is the volume of interstitial fluid in mL. Defaults to 0.15.</dd>
+<dt><strong><code>PS</code></strong> : <code>double</code></dt>
+<dd>PS is the permeability-surface area constant in mL/(g*min). Defaults to 1/60.</dd>
+<dt><strong><code>ymax</code></strong> : <code>int</code></dt>
+<dd>Magnitude of Gamma-var peak.</dd>
+<dt><strong><code>tmax</code></strong> : <code>double</code></dt>
+<dd>Time of which highest peak in Gamma-var appears</dd>
+<dt><strong><code>alpha</code></strong> : <code>double</code></dt>
+<dd>Scale factor</dd>
+<dt><strong><code>delay</code></strong> : <code>double</code></dt>
+<dd>Delay to start Gamma-var curve.</dd>
+</dl></section>
+<details class="source">
+<summary>
+<span>Expand source code</span>
+</summary>
+<pre><code class="python">class pk_two_comp:
+ """The pk2Comp object is a two compartment PK model
+ that outputs graphs of concentration of tracer over time.
+ """
+
+ def __init__(self, wd=pathlib.Path('Data').absolute(), filename='CTPERF005_stress.csv'):
+ """Initializes the model with default parameter values for flow, Vp, Visf, and PS.
+ Parameters
+ ----------
+ wd : path
+ Absolute path name to current directory. Defaults to ./Data
+ filename : String
+ Name of the data file you'd like to access
+
+ Attributes
+ -----------
+ time : double[]
+ list of all timepoints
+ aorta : double[]
+ concentration of tracer in aorta (input function)
+ myo : double[]
+ concentration of tracer in myocardial tissue (conc_isf)
+ 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.
+ ymax : int
+ Magnitude of Gamma-var peak.
+ tmax : double
+ Time of which highest peak in Gamma-var appears
+ alpha : double
+ Scale factor
+ delay : double
+ Delay to start Gamma-var curve.
+
+ """
+ # Subject Data
+ if os.path.basename(os.path.normpath(pathlib.Path().absolute())) != 'Data':
+ self.wd = pathlib.Path('Data').absolute()
+ else:
+ self.wd = wd
+
+ if not isinstance(filename, str):
+ raise ValueError('Filename must be a string')
+
+ self.filename = filename
+ self.time = []
+ self.aorta = []
+ self.myo = []
+
+ # Declare Variables for initial conditions
+ # Compartment variables to be fitted
+ self.flow = 1/60
+ self.visf = 0.15
+ self.baseline = 60
+
+ # Other Compartmental Modelvariables
+ self.perm_surf = 0.35
+ self.vol_plasma = 0.10
+
+ # Solved ode
+ self.sol = []
+
+ # Gamma variables
+ self.ymax = 250
+ self.tmax = 6.5
+ self.alpha = 2.5
+ self.delay = 0
+
+ self.deriv_sol = np.array([])
+ self.fit_myo = np.array([])
+
+ def get_data(self, filename):
+ """Imports data from all .csv files in directory.
+ Parameters
+ ----------
+ filename : string
+ Name of the data file you'd like to access
+ wd : str
+ wd is the working directory path
+
+ Attributes
+ ----------
+ time : double[]
+ list of all timepoints
+ aorta : double[]
+ concentration of tracer in aorta (input function)
+ myo : double[]
+ concentration of tracer in myocardial tissue (conc_isf)
+
+ Returns
+ -------
+ time : double[]
+ list of all timepoints
+ aorta : double[]
+ concentration of tracer in aorta (input function)
+ myo : double[]
+ concentration of tracer in myocardial tissue (conc_isf)
+ """
+ self.time = []
+ self.aorta = []
+ self.myo = []
+
+ os.chdir(self.wd)
+ # File not found error
+ if not os.path.isfile(filename):
+ raise ValueError(
+ "Input file does not exist: {0}. I'll quit now.".format(filename))
+
+ data = list(csv.reader(open(filename), delimiter='\t'))
+
+ for i in range(12):
+ self.time.append(
+ float(re.compile(r'\d+[.]+\d+|\d+').findall(data[i+1][0])[0]))
+ self.aorta.append(
+ float(re.compile(r'\d+[.]+\d+|\d+').findall(data[i+1][1])[0]))
+ self.myo.append(
+ float(re.compile(r'\d+[.]+\d+|\d+').findall(data[i+1][2])[0]))
+
+ return self.time, self.aorta, self.myo
+
+ # gamma_var distribution curve
+ def gamma_var(self, time=np.arange(0, 25), ymax=10, tmax=10, alpha=2, delay=0):
+ """Creates a gamma variate probability density function with given alpha,
+ location, and scale values.
+ Parameters
+ ----------
+ t : double[]
+ array of timepoints
+ ymax : double
+ peak y value of gamma distribution
+ tmax : double
+ location of 50th percentile of function
+ alpha : double
+ scale parameter
+ delay : double
+ time delay of which to start gamma distribution
+
+ Returns
+ -------
+ y : double[]
+ probability density function of your gamma variate.
+ """
+ # Following Madsen 1992 simplified parameterization for gamma variate
+ t = time
+ self.ymax = ymax
+ self.tmax = tmax
+ self.alpha = alpha
+ self.delay = delay
+
+ y = np.zeros(np.size(t)) # preallocate output
+
+ # For odeint, checks if t input is array or float
+ if isinstance(t, (list, np.ndarray)):
+ for i in range(np.size(y)):
+ if t[i] < delay:
+ y[i] = 0
+ else:
+ y[i] = round((ymax*tmax**(-alpha)*math.exp(alpha))*(t[i]-delay)
+ ** alpha*math.exp(-alpha*(t[i]-delay)/tmax), 3)
+ else:
+ y = round((ymax*tmax**(-alpha)*math.exp(alpha))*(t-delay)
+ ** alpha*math.exp(-alpha*(t-delay)/tmax), 3)
+ return y
+
+ # gamma_var_error
+ def input_mse(self, guess=[10, 10, 2, 5]):
+ """Calculates Mean squared error (MSE) between data and
+ gamma variate with given parameters values.
+ Parameters
+ ----------
+ param : ndarray[]
+ time : double[]
+ array of timepoints
+ ymax : double
+ peak y value of gamma distribution
+ tmax : double
+ location of 50th percentile of function
+ alpha : double
+ scale parameter
+ delay : double
+ time delay of which to start gamma distribution
+
+ Returns
+ -------
+ MSE : double
+ Mean squared error
+ """
+ if len(guess) < 1:
+ self.ymax = 10
+ self.tmax = 10
+ self.alpha = 2
+ self.delay = 5
+ elif len(guess) < 2:
+ self.ymax = guess[0]
+ self.tmax = 10
+ self.alpha = 2
+ self.delay = 5
+ elif len(guess) < 3:
+ self.ymax = guess[0]
+ self.tmax = guess[1]
+ self.alpha = 2
+ self.delay = 5
+ elif len(guess) < 4:
+ self.ymax = guess[0]
+ self.tmax = guess[1]
+ self.alpha = guess[2]
+ self.delay = 5
+ else:
+ # Mean squared error (MSE) between data and gamma variate with given parameters
+ self.ymax = guess[0]
+ self.tmax = guess[1]
+ self.alpha = guess[2]
+ self.delay = guess[3]
+
+ mse = 0
+
+ if self.tmax <= 0 or self.ymax <= 10 or self.delay < 0 or self.alpha < 0 \
+ or self.alpha > 1000 or self.tmax > 1000:
+ mse = 1000000 # just return a big number
+
+ else:
+ model_vals = self.gamma_var(
+ self.time, self.ymax, self.tmax, self.alpha, self.delay)
+
+ for i in range(len(self.aorta)):
+ mse = (self.aorta[i] - model_vals[i])**2 + mse
+ mse = mse / len(self.aorta)
+ return round(mse, 3)
+
+ def input_func_fit(self, initGuesses):
+ """Uses fmin algorithm (Nelder-Mead simplex algorithm) to
+ minimize loss function (MSE) of input function.
+ Parameters
+ ----------
+ initGuesses : ndarray[]
+ Array of initial guesses containing:
+ time : double[]
+ array of timepoints
+ ymax : double
+ peak y value of gamma distribution
+ tmax : double
+ location of 50th percentile of function
+ alpha : double
+ scale parameter
+ delay : double
+ time delay of which to start gamma distribution
+ Returns
+ -------
+ opt : double[]
+ optimized parameters
+ """
+ # Mean squared error (MSE) between data and gamma variate with given parameters
+ opt = fmin(self.input_mse, initGuesses, maxiter=1000)
+
+ self.ymax = opt[0]
+ self.tmax = opt[1]
+ self.alpha = opt[2]
+ self.delay = opt[3]
+
+ return opt.round(2)
+
+ # Derivative function
+ def derivs(self, time, curr_vals):
+ """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
+ -------
+ dconc_plasma_dt : double[]
+ contains the derivative of concentration in plasma with respect to time.
+ dconc_isf_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
+ conc_plasma, conc_isf = curr_vals
+
+ # Define value of input function conc_in
+ conc_in = self.gamma_var(time, self.ymax, self.tmax,
+ self.alpha, self.delay)
+
+ # Right-hand side of odes, which are used to computer the derivative
+ dconc_plasma_dt = (self.flow/self.vol_plasma)*(conc_in - conc_plasma) + \
+ (self.perm_surf/self.vol_plasma)*(conc_isf - conc_plasma)
+ dconc_isf_dt = (self.perm_surf/self.visf)*(conc_plasma - conc_isf)
+ return dconc_plasma_dt, dconc_isf_dt
+
+ def output_mse(self, guess):
+ """Calculates Mean squared error (MSE) between data and
+ output derivs with given parameters values.
+ Parameters
+ ----------
+ guess : ndarray[]
+
+ 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.
+ Returns
+ -------
+ MSE : double
+ Mean squared error
+ """
+ self.flow = guess[0]
+ self.visf = guess[1]
+ self.baseline = guess[2]
+
+ mse = 0
+
+ if self.flow <= 0 or self.flow >= 25 or self.visf > 100 \
+ or self.visf < 0 or self.baseline > 150 or self.baseline < 0:
+ mse = 100000 # just return a big number
+
+ else:
+ sol = solve_ivp(self.derivs, [0, 30], [0, 0], t_eval=self.time)
+ MBF = sol.y[0] + sol.y[1]
+
+ temp = np.asarray(self.myo) - self.baseline
+
+ for i in range(len(self.myo)):
+ mse = (temp[i] - MBF[i])**2 + mse
+
+ mse = mse / len(self.myo)
+ return mse
+
+ def output_func_fit(self, initGuesses):
+ """Uses fmin algorithm (Nelder-Mead simplex algorithm) to minimize
+ loss function (MSE) of output function.
+ Parameters
+ ----------
+ initGuesses : ndarray[]
+ Array of initial guesses containing:
+ 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.
+ Returns
+ -------
+ opt : double[]
+ optimized parameters
+ """
+ # Mean squared error (MSE) between data and gamma variate with given parameters
+ opt1 = fmin(self.output_mse, initGuesses, maxiter=10000)
+
+ self.flow = opt1[0].round(4)
+ self.visf = opt1[1].round(4)
+ self.baseline = opt1[2].round(4)
+
+ return opt1.round(4)
+
+ def main(self):
+
+ # Gets data from file
+ self.get_data(self.filename)
+
+ # Plots original data
+ plt.plot(self.time, self.aorta, 'bo')
+ plt.plot(self.time, self.myo, 'ro')
+
+ # Fit gamma_var input function and plots it
+ opt = self.input_func_fit([250, 7, 4, 0])
+ plt.plot(np.arange(0, 25, 0.01), self.gamma_var(np.arange(0, 25, 0.01),
+ opt[0], opt[1], opt[2], opt[3]), 'k-')
+
+ # Fit uptake function and plot it
+ opt2 = self.output_func_fit([.011418, .62, self.myo[0]])
+ self.deriv_sol = solve_ivp(self.derivs, [0, 30],
+ [0, 0], t_eval=self.time)
+ self.fit_myo = self.deriv_sol.y[0] + self.deriv_sol.y[1]
+ plt.plot(self.time, self.fit_myo + self.baseline, 'm-')</code></pre>
+</details>
+<h3>Methods</h3>
+<dl>
+<dt id="pk_two_comp.pk_two_comp.derivs"><code class="name flex">
+<span>def <span class="ident">derivs</span></span>(<span>self, time, curr_vals)</span>
+</code></dt>
+<dd>
+<section class="desc"><p>Finds derivatives of ODEs.</p>
+<h2 id="parameters">Parameters</h2>
+<dl>
+<dt><strong><code>curr_vals</code></strong> : <code>double</code>[]</dt>
+<dd>curr_vals it he current values of the variables we wish to
+"update" from the curr_vals list.</dd>
+<dt><strong><code>time</code></strong> : <code>double</code>[]</dt>
+<dd>time is our time array from 0 to tmax with timestep dt.</dd>
+</dl>
+<h2 id="returns">Returns</h2>
+<dl>
+<dt><strong><code>dconc_plasma_dt</code></strong> : <code>double</code>[]</dt>
+<dd>contains the derivative of concentration in plasma with respect to time.</dd>
+<dt><strong><code>dconc_isf_dt</code></strong> : <code>double</code>[]</dt>
+<dd>contains the derivative of concentration in interstitial fluid with respect to time.</dd>
+</dl></section>
+<details class="source">
+<summary>
+<span>Expand source code</span>
+</summary>
+<pre><code class="python">def derivs(self, time, curr_vals):
+ """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
+ -------
+ dconc_plasma_dt : double[]
+ contains the derivative of concentration in plasma with respect to time.
+ dconc_isf_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
+ conc_plasma, conc_isf = curr_vals
+
+ # Define value of input function conc_in
+ conc_in = self.gamma_var(time, self.ymax, self.tmax,
+ self.alpha, self.delay)
+
+ # Right-hand side of odes, which are used to computer the derivative
+ dconc_plasma_dt = (self.flow/self.vol_plasma)*(conc_in - conc_plasma) + \
+ (self.perm_surf/self.vol_plasma)*(conc_isf - conc_plasma)
+ dconc_isf_dt = (self.perm_surf/self.visf)*(conc_plasma - conc_isf)
+ return dconc_plasma_dt, dconc_isf_dt</code></pre>
+</details>
+</dd>
+<dt id="pk_two_comp.pk_two_comp.gamma_var"><code class="name flex">
+<span>def <span class="ident">gamma_var</span></span>(<span>self, time=array([ 0,
+1,
+2,
+3,
+4,
+5,
+6,
+7,
+8,
+9, 10, 11, 12, 13, 14, 15, 16,
+17, 18, 19, 20, 21, 22, 23, 24]), ymax=10, tmax=10, alpha=2, delay=0)</span>
+</code></dt>
+<dd>
+<section class="desc"><p>Creates a gamma variate probability density function with given alpha,
+location, and scale values.
+Parameters</p>
+<hr>
+<dl>
+<dt><strong><code>t</code></strong> : <code>double</code>[]</dt>
+<dd>array of timepoints</dd>
+<dt><strong><code>ymax</code></strong> : <code>double</code></dt>
+<dd>peak y value of gamma distribution</dd>
+<dt><strong><code>tmax</code></strong> : <code>double</code></dt>
+<dd>location of 50th percentile of function</dd>
+<dt><strong><code>alpha</code></strong> : <code>double</code></dt>
+<dd>scale parameter</dd>
+<dt><strong><code>delay</code></strong> : <code>double</code></dt>
+<dd>time delay of which to start gamma distribution</dd>
+</dl>
+<h2 id="returns">Returns</h2>
+<dl>
+<dt><strong><code>y</code></strong> : <code>double</code>[]</dt>
+<dd>probability density function of your gamma variate.</dd>
+</dl></section>
+<details class="source">
+<summary>
+<span>Expand source code</span>
+</summary>
+<pre><code class="python">def gamma_var(self, time=np.arange(0, 25), ymax=10, tmax=10, alpha=2, delay=0):
+ """Creates a gamma variate probability density function with given alpha,
+ location, and scale values.
+ Parameters
+ ----------
+ t : double[]
+ array of timepoints
+ ymax : double
+ peak y value of gamma distribution
+ tmax : double
+ location of 50th percentile of function
+ alpha : double
+ scale parameter
+ delay : double
+ time delay of which to start gamma distribution
+
+ Returns
+ -------
+ y : double[]
+ probability density function of your gamma variate.
+ """
+ # Following Madsen 1992 simplified parameterization for gamma variate
+ t = time
+ self.ymax = ymax
+ self.tmax = tmax
+ self.alpha = alpha
+ self.delay = delay
+
+ y = np.zeros(np.size(t)) # preallocate output
+
+ # For odeint, checks if t input is array or float
+ if isinstance(t, (list, np.ndarray)):
+ for i in range(np.size(y)):
+ if t[i] < delay:
+ y[i] = 0
+ else:
+ y[i] = round((ymax*tmax**(-alpha)*math.exp(alpha))*(t[i]-delay)
+ ** alpha*math.exp(-alpha*(t[i]-delay)/tmax), 3)
+ else:
+ y = round((ymax*tmax**(-alpha)*math.exp(alpha))*(t-delay)
+ ** alpha*math.exp(-alpha*(t-delay)/tmax), 3)
+ return y</code></pre>
+</details>
+</dd>
+<dt id="pk_two_comp.pk_two_comp.get_data"><code class="name flex">
+<span>def <span class="ident">get_data</span></span>(<span>self, filename)</span>
+</code></dt>
+<dd>
+<section class="desc"><p>Imports data from all .csv files in directory.
+Parameters</p>
+<hr>
+<dl>
+<dt><strong><code>filename</code></strong> : <code>string</code></dt>
+<dd>Name of the data file you'd like to access</dd>
+<dt><strong><code>wd</code></strong> : <code>str</code></dt>
+<dd>wd is the working directory path</dd>
+</dl>
+<h2 id="attributes">Attributes</h2>
+<dl>
+<dt><strong><code>time</code></strong> : <code>double</code>[]</dt>
+<dd>list of all timepoints</dd>
+<dt><strong><code>aorta</code></strong> : <code>double</code>[]</dt>
+<dd>concentration of tracer in aorta (input function)</dd>
+<dt><strong><code>myo</code></strong> : <code>double</code>[]</dt>
+<dd>concentration of tracer in myocardial tissue (conc_isf)</dd>
+</dl>
+<h2 id="returns">Returns</h2>
+<dl>
+<dt><strong><code>time</code></strong> : <code>double</code>[]</dt>
+<dd>list of all timepoints</dd>
+<dt><strong><code>aorta</code></strong> : <code>double</code>[]</dt>
+<dd>concentration of tracer in aorta (input function)</dd>
+<dt><strong><code>myo</code></strong> : <code>double</code>[]</dt>
+<dd>concentration of tracer in myocardial tissue (conc_isf)</dd>
+</dl></section>
+<details class="source">
+<summary>
+<span>Expand source code</span>
+</summary>
+<pre><code class="python">def get_data(self, filename):
+ """Imports data from all .csv files in directory.
+ Parameters
+ ----------
+ filename : string
+ Name of the data file you'd like to access
+ wd : str
+ wd is the working directory path
+
+ Attributes
+ ----------
+ time : double[]
+ list of all timepoints
+ aorta : double[]
+ concentration of tracer in aorta (input function)
+ myo : double[]
+ concentration of tracer in myocardial tissue (conc_isf)
+
+ Returns
+ -------
+ time : double[]
+ list of all timepoints
+ aorta : double[]
+ concentration of tracer in aorta (input function)
+ myo : double[]
+ concentration of tracer in myocardial tissue (conc_isf)
+ """
+ self.time = []
+ self.aorta = []
+ self.myo = []
+
+ os.chdir(self.wd)
+ # File not found error
+ if not os.path.isfile(filename):
+ raise ValueError(
+ "Input file does not exist: {0}. I'll quit now.".format(filename))
+
+ data = list(csv.reader(open(filename), delimiter='\t'))
+
+ for i in range(12):
+ self.time.append(
+ float(re.compile(r'\d+[.]+\d+|\d+').findall(data[i+1][0])[0]))
+ self.aorta.append(
+ float(re.compile(r'\d+[.]+\d+|\d+').findall(data[i+1][1])[0]))
+ self.myo.append(
+ float(re.compile(r'\d+[.]+\d+|\d+').findall(data[i+1][2])[0]))
+
+ return self.time, self.aorta, self.myo</code></pre>
+</details>
+</dd>
+<dt id="pk_two_comp.pk_two_comp.input_func_fit"><code class="name flex">
+<span>def <span class="ident">input_func_fit</span></span>(<span>self, initGuesses)</span>
+</code></dt>
+<dd>
+<section class="desc"><p>Uses fmin algorithm (Nelder-Mead simplex algorithm) to
+minimize loss function (MSE) of input function.
+Parameters</p>
+<hr>
+<dl>
+<dt><strong><code>initGuesses</code></strong> : <code>ndarray</code>[]</dt>
+<dd>Array of initial guesses containing:
+time : double[]
+array of timepoints
+ymax : double
+peak y value of gamma distribution
+tmax : double
+location of 50th percentile of function
+alpha : double
+scale parameter
+delay : double
+time delay of which to start gamma distribution</dd>
+</dl>
+<h2 id="returns">Returns</h2>
+<dl>
+<dt><strong><code>opt</code></strong> : <code>double</code>[]</dt>
+<dd>optimized parameters</dd>
+</dl></section>
+<details class="source">
+<summary>
+<span>Expand source code</span>
+</summary>
+<pre><code class="python">def input_func_fit(self, initGuesses):
+ """Uses fmin algorithm (Nelder-Mead simplex algorithm) to
+ minimize loss function (MSE) of input function.
+ Parameters
+ ----------
+ initGuesses : ndarray[]
+ Array of initial guesses containing:
+ time : double[]
+ array of timepoints
+ ymax : double
+ peak y value of gamma distribution
+ tmax : double
+ location of 50th percentile of function
+ alpha : double
+ scale parameter
+ delay : double
+ time delay of which to start gamma distribution
+ Returns
+ -------
+ opt : double[]
+ optimized parameters
+ """
+ # Mean squared error (MSE) between data and gamma variate with given parameters
+ opt = fmin(self.input_mse, initGuesses, maxiter=1000)
+
+ self.ymax = opt[0]
+ self.tmax = opt[1]
+ self.alpha = opt[2]
+ self.delay = opt[3]
+
+ return opt.round(2)</code></pre>
+</details>
+</dd>
+<dt id="pk_two_comp.pk_two_comp.input_mse"><code class="name flex">
+<span>def <span class="ident">input_mse</span></span>(<span>self, guess=[10, 10, 2, 5])</span>
+</code></dt>
+<dd>
+<section class="desc"><p>Calculates Mean squared error (MSE) between data and
+gamma variate with given parameters values.
+Parameters</p>
+<hr>
+<dl>
+<dt><strong><code>param</code></strong> : <code>ndarray</code>[]</dt>
+<dd>time : double[]
+array of timepoints
+ymax : double
+peak y value of gamma distribution
+tmax : double
+location of 50th percentile of function
+alpha : double
+scale parameter
+delay : double
+time delay of which to start gamma distribution</dd>
+</dl>
+<h2 id="returns">Returns</h2>
+<dl>
+<dt><strong><code>MSE</code></strong> : <code>double</code></dt>
+<dd>Mean squared error</dd>
+</dl></section>
+<details class="source">
+<summary>
+<span>Expand source code</span>
+</summary>
+<pre><code class="python">def input_mse(self, guess=[10, 10, 2, 5]):
+ """Calculates Mean squared error (MSE) between data and
+ gamma variate with given parameters values.
+ Parameters
+ ----------
+ param : ndarray[]
+ time : double[]
+ array of timepoints
+ ymax : double
+ peak y value of gamma distribution
+ tmax : double
+ location of 50th percentile of function
+ alpha : double
+ scale parameter
+ delay : double
+ time delay of which to start gamma distribution
+
+ Returns
+ -------
+ MSE : double
+ Mean squared error
+ """
+ if len(guess) < 1:
+ self.ymax = 10
+ self.tmax = 10
+ self.alpha = 2
+ self.delay = 5
+ elif len(guess) < 2:
+ self.ymax = guess[0]
+ self.tmax = 10
+ self.alpha = 2
+ self.delay = 5
+ elif len(guess) < 3:
+ self.ymax = guess[0]
+ self.tmax = guess[1]
+ self.alpha = 2
+ self.delay = 5
+ elif len(guess) < 4:
+ self.ymax = guess[0]
+ self.tmax = guess[1]
+ self.alpha = guess[2]
+ self.delay = 5
+ else:
+ # Mean squared error (MSE) between data and gamma variate with given parameters
+ self.ymax = guess[0]
+ self.tmax = guess[1]
+ self.alpha = guess[2]
+ self.delay = guess[3]
+
+ mse = 0
+
+ if self.tmax <= 0 or self.ymax <= 10 or self.delay < 0 or self.alpha < 0 \
+ or self.alpha > 1000 or self.tmax > 1000:
+ mse = 1000000 # just return a big number
+
+ else:
+ model_vals = self.gamma_var(
+ self.time, self.ymax, self.tmax, self.alpha, self.delay)
+
+ for i in range(len(self.aorta)):
+ mse = (self.aorta[i] - model_vals[i])**2 + mse
+ mse = mse / len(self.aorta)
+ return round(mse, 3)</code></pre>
+</details>
+</dd>
+<dt id="pk_two_comp.pk_two_comp.main"><code class="name flex">
+<span>def <span class="ident">main</span></span>(<span>self)</span>
+</code></dt>
+<dd>
+<section class="desc"></section>
+<details class="source">
+<summary>
+<span>Expand source code</span>
+</summary>
+<pre><code class="python">def main(self):
+
+ # Gets data from file
+ self.get_data(self.filename)
+
+ # Plots original data
+ plt.plot(self.time, self.aorta, 'bo')
+ plt.plot(self.time, self.myo, 'ro')
+
+ # Fit gamma_var input function and plots it
+ opt = self.input_func_fit([250, 7, 4, 0])
+ plt.plot(np.arange(0, 25, 0.01), self.gamma_var(np.arange(0, 25, 0.01),
+ opt[0], opt[1], opt[2], opt[3]), 'k-')
+
+ # Fit uptake function and plot it
+ opt2 = self.output_func_fit([.011418, .62, self.myo[0]])
+ self.deriv_sol = solve_ivp(self.derivs, [0, 30],
+ [0, 0], t_eval=self.time)
+ self.fit_myo = self.deriv_sol.y[0] + self.deriv_sol.y[1]
+ plt.plot(self.time, self.fit_myo + self.baseline, 'm-')</code></pre>
+</details>
+</dd>
+<dt id="pk_two_comp.pk_two_comp.output_func_fit"><code class="name flex">
+<span>def <span class="ident">output_func_fit</span></span>(<span>self, initGuesses)</span>
+</code></dt>
+<dd>
+<section class="desc"><p>Uses fmin algorithm (Nelder-Mead simplex algorithm) to minimize
+loss function (MSE) of output function.
+Parameters</p>
+<hr>
+<dl>
+<dt><strong><code>initGuesses</code></strong> : <code>ndarray</code>[]</dt>
+<dd>
+<p>Array of initial guesses containing:
+flow : double
+Flow is the flow of plasma through the blood vessel in mL/(mL*min).
+Defaults to 1/60.</p>
+<p>Vp : double
+Vp is the volume of plasma in mL. Defaults to 0.05.</p>
+<p>Visf : double
+Visf is the volume of interstitial fluid in mL. Defaults to 0.15.</p>
+<p>PS : double
+PS is the permeability-surface area constant in mL/(g*min). Defaults to 1/60.</p>
+</dd>
+</dl>
+<h2 id="returns">Returns</h2>
+<dl>
+<dt><strong><code>opt</code></strong> : <code>double</code>[]</dt>
+<dd>optimized parameters</dd>
+</dl></section>
+<details class="source">
+<summary>
+<span>Expand source code</span>
+</summary>
+<pre><code class="python">def output_func_fit(self, initGuesses):
+ """Uses fmin algorithm (Nelder-Mead simplex algorithm) to minimize
+ loss function (MSE) of output function.
+ Parameters
+ ----------
+ initGuesses : ndarray[]
+ Array of initial guesses containing:
+ 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.
+ Returns
+ -------
+ opt : double[]
+ optimized parameters
+ """
+ # Mean squared error (MSE) between data and gamma variate with given parameters
+ opt1 = fmin(self.output_mse, initGuesses, maxiter=10000)
+
+ self.flow = opt1[0].round(4)
+ self.visf = opt1[1].round(4)
+ self.baseline = opt1[2].round(4)
+
+ return opt1.round(4)</code></pre>
+</details>
+</dd>
+<dt id="pk_two_comp.pk_two_comp.output_mse"><code class="name flex">
+<span>def <span class="ident">output_mse</span></span>(<span>self, guess)</span>
+</code></dt>
+<dd>
+<section class="desc"><p>Calculates Mean squared error (MSE) between data and
+output derivs with given parameters values.
+Parameters</p>
+<hr>
+<dl>
+<dt><strong><code>guess</code></strong> : <code>ndarray</code>[]</dt>
+<dd>
+<p>Flow : double
+Flow is the flow of plasma through the blood vessel in mL/(mL*min).
+Defaults to 1/60.</p>
+<p>Vp : double
+Vp is the volume of plasma in mL. Defaults to 0.05.</p>
+<p>Visf : double
+Visf is the volume of interstitial fluid in mL. Defaults to 0.15.</p>
+<p>PS : double
+PS is the permeability-surface area constant in mL/(g*min). Defaults to 1/60.</p>
+</dd>
+</dl>
+<h2 id="returns">Returns</h2>
+<dl>
+<dt><strong><code>MSE</code></strong> : <code>double</code></dt>
+<dd>Mean squared error</dd>
+</dl></section>
+<details class="source">
+<summary>
+<span>Expand source code</span>
+</summary>
+<pre><code class="python">def output_mse(self, guess):
+ """Calculates Mean squared error (MSE) between data and
+ output derivs with given parameters values.
+ Parameters
+ ----------
+ guess : ndarray[]
+
+ 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.
+ Returns
+ -------
+ MSE : double
+ Mean squared error
+ """
+ self.flow = guess[0]
+ self.visf = guess[1]
+ self.baseline = guess[2]
+
+ mse = 0
+
+ if self.flow <= 0 or self.flow >= 25 or self.visf > 100 \
+ or self.visf < 0 or self.baseline > 150 or self.baseline < 0:
+ mse = 100000 # just return a big number
+
+ else:
+ sol = solve_ivp(self.derivs, [0, 30], [0, 0], t_eval=self.time)
+ MBF = sol.y[0] + sol.y[1]
+
+ temp = np.asarray(self.myo) - self.baseline
+
+ for i in range(len(self.myo)):
+ mse = (temp[i] - MBF[i])**2 + mse
+
+ mse = mse / len(self.myo)
+ return mse</code></pre>
+</details>
+</dd>
+</dl>
+</dd>
+</dl>
+</section>
+</article>
+<nav id="sidebar">
+<h1>Index</h1>
+<div class="toc">
+<ul></ul>
+</div>
+<ul id="index">
+<li><h3><a href="#header-classes">Classes</a></h3>
+<ul>
+<li>
+<h4><code><a title="pk_two_comp.pk_two_comp" href="#pk_two_comp.pk_two_comp">pk_two_comp</a></code></h4>
+<ul class="two-column">
+<li><code><a title="pk_two_comp.pk_two_comp.derivs" href="#pk_two_comp.pk_two_comp.derivs">derivs</a></code></li>
+<li><code><a title="pk_two_comp.pk_two_comp.gamma_var" href="#pk_two_comp.pk_two_comp.gamma_var">gamma_var</a></code></li>
+<li><code><a title="pk_two_comp.pk_two_comp.get_data" href="#pk_two_comp.pk_two_comp.get_data">get_data</a></code></li>
+<li><code><a title="pk_two_comp.pk_two_comp.input_func_fit" href="#pk_two_comp.pk_two_comp.input_func_fit">input_func_fit</a></code></li>
+<li><code><a title="pk_two_comp.pk_two_comp.input_mse" href="#pk_two_comp.pk_two_comp.input_mse">input_mse</a></code></li>
+<li><code><a title="pk_two_comp.pk_two_comp.main" href="#pk_two_comp.pk_two_comp.main">main</a></code></li>
+<li><code><a title="pk_two_comp.pk_two_comp.output_func_fit" href="#pk_two_comp.pk_two_comp.output_func_fit">output_func_fit</a></code></li>
+<li><code><a title="pk_two_comp.pk_two_comp.output_mse" href="#pk_two_comp.pk_two_comp.output_mse">output_mse</a></code></li>
+</ul>
+</li>
+</ul>
+</li>
+</ul>
+</nav>
+</main>
+<footer id="footer">
+<p>Generated by <a href="https://pdoc3.github.io/pdoc"><cite>pdoc</cite> 0.7.4</a>.</p>
+</footer>
+<script src="https://cdnjs.cloudflare.com/ajax/libs/highlight.js/9.12.0/highlight.min.js"></script>
+<script>hljs.initHighlightingOnLoad()</script>
+</body>
+</html>
\ No newline at end of file