{ "cells": [ { "cell_type": "code", "execution_count": 1, "metadata": { "collapsed": false }, "outputs": [ ], "source": [ "# Configure Jupyter so figures appear in the notebook\n", "%matplotlib inline\n", "\n", "# Configure Jupyter to display the assigned value after an assignment\n", "%config InteractiveShell.ast_node_interactivity='last_expr_or_assign'\n", "\n", "# import functions from the modsim.py module\n", "from modsim import *" ] }, { "cell_type": "code", "execution_count": 0, "metadata": { "collapsed": false }, "outputs": [ ], "source": [ ] }, { "cell_type": "code", "execution_count": 2, "metadata": { "collapsed": false }, "outputs": [ ], "source": [ "def make_system(gamma,mu,beta,tau,roe,alpha,pi,sigma,delta):\n", " \"\"\"Make a system object for the SIR model.\n", " \n", " beta: contact rate in days\n", " gamma: recovery rate in days\n", " \n", " returns: System object\n", " \"\"\"\n", " init = State(R=200, L=0, E=0, V=4e-7)\n", "\n", " t_0 = 0\n", " t_end = 120\n", "\n", " return System(init=init, t_0=t_0, t_end=t_end,\n", " beta=beta, gamma=gamma,mu=mu,tau=tau,roe=roe,alpha=alpha,sigma=sigma, pi=pi, delta = delta, dt=.0541)" ] }, { "cell_type": "code", "execution_count": 3, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/html": [ "
\n", "\n", "\n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", "
values
initR 2.000000e+02\n", "L 0.000000e+00\n", "E 0.000...
t_00
t_end120
beta0.00027
gamma1.36
mu0.00136
tau0.2
roe0.1
alpha0.036
sigma2
pi100
delta0.33
dt0.0541
\n", "
" ], "text/plain": [ "init R 2.000000e+02\n", "L 0.000000e+00\n", "E 0.000...\n", "t_0 0\n", "t_end 120\n", "beta 0.00027\n", "gamma 1.36\n", "mu 0.00136\n", "tau 0.2\n", "roe 0.1\n", "alpha 0.036\n", "sigma 2\n", "pi 100\n", "delta 0.33\n", "dt 0.0541\n", "dtype: object" ] }, "execution_count": 3, "metadata": { }, "output_type": "execute_result" } ], "source": [ "system = make_system(1.36,0.00136,.00027,0.2,0.1,0.036,100,2,.33)" ] }, { "cell_type": "code", "execution_count": 4, "metadata": { "collapsed": false }, "outputs": [ ], "source": [ "def update_func(state, t, system):\n", " R, L, E, V = state\n", " diffR = R\n", " diffL = L\n", " diffE = E\n", " diffV = V\n", " diffR = system.gamma*system.tau-system.mu*R-system.beta*R*V\n", "a diffE = (1-system.roe)*system.beta* R * V + system.alpha * L - system.delta * E\n", " diffV = system.pi * E - system.sigma * V\n", " R += diffR * system.dt\n", " L += diffL * system.dt\n", " E += diffE * system.dt\n", " V += diffV * system.dt\n", " return State(R=R, L=L, E=E, V=V)" ] }, { "cell_type": "code", "execution_count": 5, "metadata": { "collapsed": false }, "outputs": [ ], "source": [ "def run_simulation(system, update_func):\n", " \"\"\"Runs a simulation of the system.\n", " \n", " system: System object\n", " update_func: function that updates state\n", " \n", " returns: TimeFrame\n", " \"\"\"\n", " init = system.init\n", " t_0, t_end, dt = system.t_0, system.t_end, system.dt\n", " \n", " frame = TimeFrame(columns=init.index)\n", " frame.row[t_0] = init\n", " ts = linrange(t_0, t_end, system.dt)\n", " \n", " for t in ts:\n", " frame.row[t+dt] = update_func(frame.row[t], t, system)\n", " \n", " return frame" ] }, { "cell_type": "code", "execution_count": 0, "metadata": { "collapsed": false }, "outputs": [ ], "source": [ ] }, { "cell_type": "code", "execution_count": 6, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/html": [ "
\n", "\n", "\n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", "
RLEV
0.0000200.0000000.000000e+000.000000e+004.000000e-07
0.0541200.0000001.168560e-101.051704e-093.567200e-07
0.1082200.0000002.208320e-101.971065e-093.238126e-07
0.1623200.0000003.149843e-102.787693e-092.994396e-07
0.2164200.0000004.018259e-103.525843e-092.821216e-07
...............
119.777418.0081686.612353e-012.387806e-011.200246e+01
119.831518.0184016.602145e-012.386470e-011.199560e+01
119.885618.0286336.591958e-012.385137e-011.198875e+01
119.939718.0388656.581792e-012.383808e-011.198193e+01
119.993818.0490966.571646e-012.382483e-011.197512e+01
\n", "

2219 rows × 4 columns

\n", "
" ], "text/plain": [ " R L E V\n", "0.0000 200.000000 0.000000e+00 0.000000e+00 4.000000e-07\n", "0.0541 200.000000 1.168560e-10 1.051704e-09 3.567200e-07\n", "0.1082 200.000000 2.208320e-10 1.971065e-09 3.238126e-07\n", "0.1623 200.000000 3.149843e-10 2.787693e-09 2.994396e-07\n", "0.2164 200.000000 4.018259e-10 3.525843e-09 2.821216e-07\n", "... ... ... ... ...\n", "119.7774 18.008168 6.612353e-01 2.387806e-01 1.200246e+01\n", "119.8315 18.018401 6.602145e-01 2.386470e-01 1.199560e+01\n", "119.8856 18.028633 6.591958e-01 2.385137e-01 1.198875e+01\n", "119.9397 18.038865 6.581792e-01 2.383808e-01 1.198193e+01\n", "119.9938 18.049096 6.571646e-01 2.382483e-01 1.197512e+01\n", "\n", "[2219 rows x 4 columns]" ] }, "execution_count": 6, "metadata": { }, "output_type": "execute_result" } ], "source": [ "results = run_simulation(system, update_func)" ] }, { "cell_type": "code", "execution_count": 0, "metadata": { "collapsed": false }, "outputs": [ ], "source": [ "plot(results.V)\n", "decorate(title = 'Free Virions',\n", " xlabel = 'Days',\n", " ylabel = 'Virions')" ] }, { "cell_type": "code", "execution_count": 0, "metadata": { "collapsed": false }, "outputs": [ ], "source": [ "plot(results.R)\n", "decorate(title = 'Uninfected CD4',\n", " xlabel = 'Days',\n", " ylabel = 'CD4 Lymphosytes')" ] }, { "cell_type": "code", "execution_count": 0, "metadata": { "collapsed": false }, "outputs": [ ], "source": [ "plot(results.E)\n", "decorate(title = 'Infected CD4',\n", " xlabel = 'Days',\n", " ylabel = 'Infected CD4',\n", " yscale = 'log')" ] }, { "cell_type": "code", "execution_count": 0, "metadata": { "collapsed": false }, "outputs": [ ], "source": [ "plot(results.L)\n", "\n", "decorate(title = 'Free Virions',\n", " xlabel = 'Days',\n", " ylabel = 'Latently Infected CD4')" ] }, { "cell_type": "code", "execution_count": 0, "metadata": { "collapsed": false }, "outputs": [ ], "source": [ ] } ], "metadata": { "kernelspec": { "display_name": "Python 3", "language": "python", "name": "python3" }, "language_info": { "codemirror_mode": { "name": "ipython", "version": 3 }, "file_extension": ".py", "mimetype": "text/x-python", "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", "version": "3.7.3" } }, "nbformat": 4, "nbformat_minor": 0 }