Next Generation Matrix.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 37,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "{-b*(-a*d_a*d_i*l*m*s_i + a*d_a*d_i*l*s_i + a*d_a*d_i*m*s_i - a*d_a*d_i*s_i + a*d_a*m*p - a*d_a*m + d_i*l*m*n*s_i - d_i*l*n*s_i - m*n*p)/(a*n): 1,\n",
       " 0: 8}"
      ]
     },
     "execution_count": 37,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Python program to inverse\n",
    "# a matrix using numpy\n",
    "  \n",
    "# Import required package\n",
    "import numpy as np\n",
    "import matplotlib as mpl\n",
    "import matplotlib.pyplot as plt\n",
    "from sympy import *\n",
    "\n",
    "# Finding the dominant eigenvalue for the non age split model \n",
    "b,m = symbols(\"b,m\")\n",
    "d_i, d_a, d_v = symbols(\"d_i,d_a,d_v\")\n",
    "s_i = symbols(\"s_i\")\n",
    "k = symbols(\"k\")\n",
    "p,q,l = symbols(\"p,q,l\")\n",
    "a,n = symbols(\"a,n\")\n",
    "\n",
    "T = Matrix([[0,0,0, (b * m), 0, (b * d_i * m), (b * d_a * m), (b * d_v * d_a * m), (b * d_i * d_a *m)], \n",
    "            [0, 0, 0, 0, 0, 0, 0, 0, 0],\n",
    "           [0, 0, 0, (b * s_i * (1-m)), 0, (b * s_i * d_i * (1-m)), (b * s_i * d_a * (1-m)), (b * s_i * d_v * d_a *(1-m)), (b * s_i * d_i * d_a *(1-m))],\n",
    "           [0, 0, 0, 0, 0, 0, 0, 0, 0], [0, 0, 0, 0, 0, 0, 0, 0, 0], [0, 0, 0, 0, 0, 0, 0, 0, 0],\n",
    "           [0, 0, 0, 0, 0, 0, 0, 0, 0], [0, 0, 0, 0, 0, 0, 0, 0, 0], [0, 0, 0, 0, 0, 0, 0, 0, 0]])\n",
    "T\n",
    "\n",
    "K = Matrix ([[-k, 0, 0, 0, 0, 0, 0, 0, 0], [0, -k, 0, 0, 0, 0, 0, 0, 0], [0, 0, -k, 0, 0, 0, 0, 0, 0],\n",
    "            [(p*k), 0, 0, -a, 0, 0, 0, 0, 0], [0, (q*k), 0, 0, -a, 0, 0, 0, 0], [0, 0, (l*k), 0, 0, -a, 0, 0, 0],\n",
    "            [((1-p)*k), 0, 0, 0, 0, 0, -n, 0, 0], [0, ((1-q)*k), 0, 0, 0, 0, 0, -n, 0], [0, 0, ((1-l)*k), 0, 0, 0, 0, 0, -n],])\n",
    "NGM = -T*K**-1\n",
    "NGM\n",
    "NGM.eigenvals()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 61,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/latex": [
       "$\\displaystyle 3.01846153846154$"
      ],
      "text/plain": [
       "3.01846153846154"
      ]
     },
     "execution_count": 61,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "#Calculating the scaling for the non age split model\n",
    "\n",
    "import numpy as np\n",
    "import matplotlib as mpl\n",
    "import matplotlib.pyplot as plt\n",
    "from sympy import *\n",
    "from sympy import symbols, exp\n",
    "b = symbols('b')\n",
    "m = symbols('m')\n",
    "d_i = symbols('d_i')\n",
    "d_a = symbols('d_a')\n",
    "d_v = symbols('d_v')\n",
    "s_i = symbols('s_i')\n",
    "k = symbols('k')\n",
    "p = symbols('p')\n",
    "q = symbols('q')\n",
    "l = symbols('l')\n",
    "a = symbols('a')\n",
    "n = symbols('n')\n",
    "\n",
    "expr = - b * ((-a * d_a * d_i * l * m * s_i) + (a * d_a * d_i * l * s_i) + (a * d_a * d_i * m * s_i) - (a * d_a * d_i * s_i) + (a * d_a * m * p) - (a * d_a * m) + (d_i * l * m * n * s_i) - (d_i * l * n * s_i) - (m * n * p))/(a * n)\n",
    "expr.subs(b,1.5).subs(m,0.6).subs(d_i, 1).subs(d_a,0.2).subs(d_v,1).subs(s_i,0.3).subs(k,0.71).subs(p,0.75).subs(q,0.25).subs(l,0.2).subs(a,0.26).subs(n,0.26)\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 66,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Therefore, with the current parameter values, the system has an R_0 value of 3.02. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/latex": [
       "$\\displaystyle 1.97206153846154$"
      ],
      "text/plain": [
       "1.97206153846154"
      ]
     },
     "execution_count": 1,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "import matplotlib.pyplot as plt\n",
    "from sympy import *\n",
    "from sympy import symbols, exp\n",
    "b_11, b_12, b_13, b_14 = symbols('b_11, b_12, b_13, b_14')\n",
    "b_21, b_22, b_23, b_24 = symbols('b_21, b_22, b_23, b_24')\n",
    "b_31, b_32, b_33, b_34 = symbols('b_31, b_32, b_33, b_34')\n",
    "b_41, b_42, b_43, b_44 = symbols('b_41, b_42, b_43, b_44')\n",
    "m = symbols('m')\n",
    "d_i = symbols('d_i')\n",
    "d_a = symbols('d_a')\n",
    "d_v = symbols('d_v')\n",
    "s_i = symbols('s_i')\n",
    "k = symbols('k')\n",
    "p = symbols('p')\n",
    "q = symbols('q')\n",
    "l = symbols('l')\n",
    "a = symbols('a')\n",
    "n = symbols('n')\n",
    "\n",
    "expr = - b_11 * ((-a * d_a * d_i * l * m * s_i) + (a * d_a * d_i * l * s_i) + (a * d_a * d_i * m * s_i) - (a * d_a * d_i * s_i) + (a * d_a * m * p) - (a * d_a * m) + (d_i * l * m * n * s_i) - (d_i * l * n * s_i) - (m * n * p))/(a * n)\n",
    "- b_12 * ((-a * d_a * d_i * l * m * s_i) + (a * d_a * d_i * l * s_i) + (a * d_a * d_i * m * s_i) - (a * d_a * d_i * s_i) + (a * d_a * m * p) - (a * d_a * m) + (d_i * l * m * n * s_i) - (d_i * l * n * s_i) - (m * n * p))/(a * n)\n",
    "- b_13 * ((-a * d_a * d_i * l * m * s_i) + (a * d_a * d_i * l * s_i) + (a * d_a * d_i * m * s_i) - (a * d_a * d_i * s_i) + (a * d_a * m * p) - (a * d_a * m) + (d_i * l * m * n * s_i) - (d_i * l * n * s_i) - (m * n * p))/(a * n)\n",
    "- b_14 * ((-a * d_a * d_i * l * m * s_i) + (a * d_a * d_i * l * s_i) + (a * d_a * d_i * m * s_i) - (a * d_a * d_i * s_i) + (a * d_a * m * p) - (a * d_a * m) + (d_i * l * m * n * s_i) - (d_i * l * n * s_i) - (m * n * p))/(a * n)\n",
    "- b_21 * ((-a * d_a * d_i * l * m * s_i) + (a * d_a * d_i * l * s_i) + (a * d_a * d_i * m * s_i) - (a * d_a * d_i * s_i) + (a * d_a * m * p) - (a * d_a * m) + (d_i * l * m * n * s_i) - (d_i * l * n * s_i) - (m * n * p))/(a * n)\n",
    "- b_22 * ((-a * d_a * d_i * l * m * s_i) + (a * d_a * d_i * l * s_i) + (a * d_a * d_i * m * s_i) - (a * d_a * d_i * s_i) + (a * d_a * m * p) - (a * d_a * m) + (d_i * l * m * n * s_i) - (d_i * l * n * s_i) - (m * n * p))/(a * n)\n",
    "- b_23 * ((-a * d_a * d_i * l * m * s_i) + (a * d_a * d_i * l * s_i) + (a * d_a * d_i * m * s_i) - (a * d_a * d_i * s_i) + (a * d_a * m * p) - (a * d_a * m) + (d_i * l * m * n * s_i) - (d_i * l * n * s_i) - (m * n * p))/(a * n)\n",
    "- b_24 * ((-a * d_a * d_i * l * m * s_i) + (a * d_a * d_i * l * s_i) + (a * d_a * d_i * m * s_i) - (a * d_a * d_i * s_i) + (a * d_a * m * p) - (a * d_a * m) + (d_i * l * m * n * s_i) - (d_i * l * n * s_i) - (m * n * p))/(a * n)\n",
    "- b_31 * ((-a * d_a * d_i * l * m * s_i) + (a * d_a * d_i * l * s_i) + (a * d_a * d_i * m * s_i) - (a * d_a * d_i * s_i) + (a * d_a * m * p) - (a * d_a * m) + (d_i * l * m * n * s_i) - (d_i * l * n * s_i) - (m * n * p))/(a * n)\n",
    "- b_32 * ((-a * d_a * d_i * l * m * s_i) + (a * d_a * d_i * l * s_i) + (a * d_a * d_i * m * s_i) - (a * d_a * d_i * s_i) + (a * d_a * m * p) - (a * d_a * m) + (d_i * l * m * n * s_i) - (d_i * l * n * s_i) - (m * n * p))/(a * n)\n",
    "- b_33 * ((-a * d_a * d_i * l * m * s_i) + (a * d_a * d_i * l * s_i) + (a * d_a * d_i * m * s_i) - (a * d_a * d_i * s_i) + (a * d_a * m * p) - (a * d_a * m) + (d_i * l * m * n * s_i) - (d_i * l * n * s_i) - (m * n * p))/(a * n)\n",
    "- b_34 * ((-a * d_a * d_i * l * m * s_i) + (a * d_a * d_i * l * s_i) + (a * d_a * d_i * m * s_i) - (a * d_a * d_i * s_i) + (a * d_a * m * p) - (a * d_a * m) + (d_i * l * m * n * s_i) - (d_i * l * n * s_i) - (m * n * p))/(a * n)\n",
    "- b_41 * ((-a * d_a * d_i * l * m * s_i) + (a * d_a * d_i * l * s_i) + (a * d_a * d_i * m * s_i) - (a * d_a * d_i * s_i) + (a * d_a * m * p) - (a * d_a * m) + (d_i * l * m * n * s_i) - (d_i * l * n * s_i) - (m * n * p))/(a * n)\n",
    "- b_42 * ((-a * d_a * d_i * l * m * s_i) + (a * d_a * d_i * l * s_i) + (a * d_a * d_i * m * s_i) - (a * d_a * d_i * s_i) + (a * d_a * m * p) - (a * d_a * m) + (d_i * l * m * n * s_i) - (d_i * l * n * s_i) - (m * n * p))/(a * n)\n",
    "- b_43 * ((-a * d_a * d_i * l * m * s_i) + (a * d_a * d_i * l * s_i) + (a * d_a * d_i * m * s_i) - (a * d_a * d_i * s_i) + (a * d_a * m * p) - (a * d_a * m) + (d_i * l * m * n * s_i) - (d_i * l * n * s_i) - (m * n * p))/(a * n)\n",
    "- b_44 * ((-a * d_a * d_i * l * m * s_i) + (a * d_a * d_i * l * s_i) + (a * d_a * d_i * m * s_i) - (a * d_a * d_i * s_i) + (a * d_a * m * p) - (a * d_a * m) + (d_i * l * m * n * s_i) - (d_i * l * n * s_i) - (m * n * p))/(a * n)\n",
    "\n",
    "expr.subs(b_11,0.98).subs(b_12,0.25).subs(b_13,1.1).subs(b_14,0.096).subs(b_21,0.22).subs(b_22,1.5).subs(b_23,1.1).subs(b_24,0.15).subs(b_31,0.21).subs(b_32,0.24).subs(b_33,1.5).subs(b_34,0.2).subs(b_41,0.069).subs(b_42,0.12).subs(b_43,0.74).subs(b_44,0.36).subs(m,0.6).subs(d_i, 1).subs(d_a,0.2).subs(d_v,1).subs(s_i,0.3).subs(k,0.71).subs(p,0.75).subs(q,0.25).subs(l,0.2).subs(a,0.26).subs(n,0.26)\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Therefore, with the current parameter values, the age split system has an R_0 of 1.97."
   ]
  }
 ],
 "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.8.5"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}