he_model_oo.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "6e84e122",
   "metadata": {},
   "source": [
    "# Simple health economic DES model with Python"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "8672b03d",
   "metadata": {},
   "source": [
    "A very simple model to serve as an exploration of creating health economic Discrete-Event Simulation (DES) models in Python using SimPy library. The clinical context of the model is as follows:\n",
    "* During a cycle of treatment, patients can either die (p = 0.15) or experience a full cycle of treatment without any other events occurring.\n",
    "* Patients can receive up to five cycles of treatment.\n",
    "* A full cycle causes a longer delay or timeout compared to dying during a cycle. These timeouts are drawn from a Gamma distribution.\n",
    "* After surviving the maximum number of treatment cycles, patients will enter a follow up phase, in which patients can also die.\n",
    "\n",
    "The health economic context is as follows:\n",
    "* Each cycle of treatment incurs 5000 euro initially and 250 euro per day.\n",
    "* Followup costs 3500 euro per patient as a lumpsum.\n",
    "* Health-related quality of life (QoL) in the treatment phase is 0.7 as expressed over the duration of one year.\n",
    "* QoL in the followup phase is 0.8 as expressed over the duration of one year.\n",
    "\n",
    "Furthermore, all patients are simulated at `t = 0`, so the interarrival time = 0. Note that the model currently only represents one comparator. So for an actual cost-effectiveness analysis, another comparator must be added, but this is very straightforward."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "91142c9b",
   "metadata": {},
   "source": [
    "Three classes are created for the model: 1) a class `g` contains all constants, 2) a `Patient` class in which attributes of patients are set and 3) a `Model` class containing the model structure and Patient Generator method. The model is run by creating an instance of `Model` and subsequently executing the `run()` method of that `Model` instance. To accomodate the need to conduct multiple runs, this is performed in a for loop."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "id": "a312a585-bfdb-4851-ba3c-c79df22435a7",
   "metadata": {},
   "outputs": [],
   "source": [
    "# Remove the comment in the line below to install simpy if you're viewing this notebook in Google Colab. This cell needs to be run only once.\n",
    "#%pip install simpy"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "id": "486e7e65-2878-4c07-a61e-7194046ff6a6",
   "metadata": {},
   "outputs": [],
   "source": [
    "# Importing libraries and setting the random seed for reproducibility.\n",
    "from random import seed\n",
    "import random\n",
    "import simpy\n",
    "import pandas as pd\n",
    "import numpy as np\n",
    "seed(123)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "id": "cefc3edd-06ca-43b7-ab84-1839328f6d79",
   "metadata": {},
   "outputs": [],
   "source": [
    "class g:\n",
    "    max_cycles = 5 # Maximum number of treatment cycles per patient\n",
    "    prob_death = 0.15 #Probability of dying during one cycle per patient\n",
    "    n_patients = 10000 # Number of patients to simulate\n",
    "    c_treatment_init = 5000\n",
    "    c_treatment_daily = 250\n",
    "    c_followup = 3500\n",
    "    u_treatment = 0.7\n",
    "    u_followup = 0.8\n",
    "    days_per_year = 365.2422\n",
    "    sim_duration = simpy.core.Infinity\n",
    "    number_of_runs = 1"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "id": "567bce14-8c75-4569-bf24-8152d7ce2fff",
   "metadata": {},
   "outputs": [],
   "source": [
    "class Patient:\n",
    "    def __init__(self, patient_id):\n",
    "        \n",
    "        \"\"\"Initializing attributes of Patients. \n",
    "        Of course these can be expanded to more accurately reflect patient heterogeneity.\"\"\"\n",
    "        self.patient_id = patient_id\n",
    "        self.state = 'Alive'\n",
    "        self.treatment_cycles = 0\n",
    "        self.cost = 0\n",
    "        self.utility = 0"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "id": "84f94195-bc9a-4aac-b81e-632957ac4690",
   "metadata": {},
   "outputs": [],
   "source": [
    "class Model:\n",
    "    def __init__(self, run_number):\n",
    "        self.env = simpy.Environment()\n",
    "        self.patient_counter = 0\n",
    "        self.run_number = run_number\n",
    "\n",
    "    def generate_patients(self):\n",
    "        \"\"\"The method that generates patients.\n",
    "        \"\"\"\n",
    "        yield self.env.timeout(0)  # SimPy processes cannot be 'empty'\n",
    "        # Main generator loop that terminates when enough patients are simulated or\n",
    "        # when until is reached\n",
    "\n",
    "        self.run_number += 1\n",
    "        while self.env.now < g.sim_duration and self.patient_counter < g.n_patients:\n",
    "            self.patient_counter += 1\n",
    "            # Create a new instance of the Patient class\n",
    "            pat = Patient(patient_id=self.patient_counter)\n",
    "\n",
    "            # Use the SimPy environment and the enter_treatment and enter_followup methods\n",
    "            # with this patient\n",
    "            self.env.process(self.set_care_pathway(pat))\n",
    "\n",
    "    def run(self):\n",
    "        self.env.process(self.generate_patients())\n",
    "        self.env.run(until=g.sim_duration)\n",
    "\n",
    "        #################### START SECTION: MODEL STRUCTURE ####################\n",
    "\n",
    "    def set_care_pathway(self, patient):\n",
    "        \"\"\" Method that models the treatment phase.\n",
    "        \"\"\"\n",
    "\n",
    "        #### TREATMENT PHASE ####\n",
    "        while patient.treatment_cycles < g.max_cycles and patient.state == 'Alive':\n",
    "            patient.treatment_cycles += 1\n",
    "\n",
    "            # First, the event that occurs during a cycle is determined.\n",
    "            rand = random.uniform(0, 1)\n",
    "            if rand < g.prob_death:\n",
    "                ##### EVENT: DEATH\n",
    "\n",
    "                # CHANGE PATIENT'S STATE TO DEAD\n",
    "                patient.state = 'Dead'\n",
    "\n",
    "                # SAMPLE A TIME-TO-EVENT\n",
    "                time_to_death = np.random.gamma(1.5, 3, 1)\n",
    "\n",
    "                # TIMEOUT EQUAL TO THE TIME-TO-EVENT\n",
    "                yield self.env.timeout(time_to_death)\n",
    "\n",
    "                # INCREMENT ACCUMULATED COSTS AND UTILITY\n",
    "                patient.cost = self.increment_cost(patient, time_to_death)\n",
    "                patient.utility = self.increment_utility(time_to_death, g.u_treatment)\n",
    "\n",
    "            else:\n",
    "                ##### EVENT: FULL CYCLE\n",
    "                time_to_full_cycle = np.random.gamma(3, 10, 1)\n",
    "                yield self.env.timeout(time_to_full_cycle)\n",
    "                patient.cost = self.increment_cost(patient, time_to_full_cycle)\n",
    "                patient.utility = self.increment_utility(time_to_full_cycle, g.u_treatment)\n",
    "\n",
    "            # SAVE DATA AT THE END OF EACH CYCLE\n",
    "            self.save_data(patient, 'treatment')\n",
    "\n",
    "        #### FOLLOWUP PHASE ####\n",
    "        if patient.state == 'Alive':\n",
    "            time_in_folllowup = np.random.gamma(2, 15, 1)\n",
    "            yield self.env.timeout(time_in_folllowup)\n",
    "            patient.utility = self.increment_utility(time_in_folllowup, g.u_followup)\n",
    "            patient.cost = g.c_followup\n",
    "            self.save_data(patient, 'followup')\n",
    "\n",
    "            #################### END SECTION: MODEL STRUCTURE ####################\n",
    "\n",
    "    #################### START SECTION: HELPER METHODS ####################\n",
    "\n",
    "    def save_data(self, patient, phase):\n",
    "        \"\"\"Append a list of outcomes of interest as specified here to a list that\n",
    "        is created outside the Patient class. This method should be called whenever it\n",
    "        is appropriate to save data. E.g., each treatment cycle.\n",
    "        Note, appending a list to a list and converting the final list once to a dataframe is much more efficient\n",
    "        than appending directly to a pd dataframe.\"\"\"\n",
    "        output_list.append(\n",
    "            [patient.patient_id, patient.state, patient.treatment_cycles, phase, patient.cost, patient.utility,\n",
    "             self.run_number, self.env.now])\n",
    "\n",
    "    @staticmethod\n",
    "    def increment_cost(patient, duration):\n",
    "        \"\"\"Helper method to increment the cost attribute of Patient\"\"\"\n",
    "        cost_increment = int(duration) * g.c_treatment_daily\n",
    "        if patient.treatment_cycles == 1:\n",
    "            cost_increment += g.c_treatment_init\n",
    "        return cost_increment\n",
    "\n",
    "    @staticmethod\n",
    "    def increment_utility(duration, utility):\n",
    "        \"\"\"Helper method to increment the utility attribute of Patient\"\"\"\n",
    "        utility_increment = duration * (utility / g.days_per_year)\n",
    "        return utility_increment\n",
    "\n",
    "    #################### END SECTION: HELPER METHODS ####################"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "id": "295a0cdf-8690-45a3-be64-8590b5ceae6c",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Run 1 of 1\n",
      "\n"
     ]
    },
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>patient_id</th>\n",
       "      <th>state</th>\n",
       "      <th>treatment_cycle</th>\n",
       "      <th>phase</th>\n",
       "      <th>cost</th>\n",
       "      <th>utility</th>\n",
       "      <th>run_number</th>\n",
       "      <th>simulation_time</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>6349</th>\n",
       "      <td>7</td>\n",
       "      <td>Alive</td>\n",
       "      <td>1</td>\n",
       "      <td>treatment</td>\n",
       "      <td>11250</td>\n",
       "      <td>0.049277</td>\n",
       "      <td>1</td>\n",
       "      <td>25.711225</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>8846</th>\n",
       "      <td>7</td>\n",
       "      <td>Alive</td>\n",
       "      <td>2</td>\n",
       "      <td>treatment</td>\n",
       "      <td>2000</td>\n",
       "      <td>0.015548</td>\n",
       "      <td>1</td>\n",
       "      <td>33.823813</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>9387</th>\n",
       "      <td>7</td>\n",
       "      <td>Dead</td>\n",
       "      <td>3</td>\n",
       "      <td>treatment</td>\n",
       "      <td>250</td>\n",
       "      <td>0.003434</td>\n",
       "      <td>1</td>\n",
       "      <td>35.615664</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "      patient_id  state  treatment_cycle      phase   cost   utility  \\\n",
       "6349           7  Alive                1  treatment  11250  0.049277   \n",
       "8846           7  Alive                2  treatment   2000  0.015548   \n",
       "9387           7   Dead                3  treatment    250  0.003434   \n",
       "\n",
       "      run_number  simulation_time  \n",
       "6349           1        25.711225  \n",
       "8846           1        33.823813  \n",
       "9387           1        35.615664  "
      ]
     },
     "execution_count": 15,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Empty list to store outcomes of interest in. Used by the Patient.save_data() method.            \n",
    "output_list = []\n",
    "\n",
    "# Running the model for the required amount of runs and printing progress.\n",
    "for run in range(g.number_of_runs):\n",
    "    print(\"Run \", run+1, \" of \", g.number_of_runs, sep=\"\")\n",
    "    my_model = Model(run)\n",
    "    my_model.run()\n",
    "    print()\n",
    "\n",
    "output_df = pd.DataFrame(output_list, columns = ['patient_id', 'state', 'treatment_cycle', 'phase', 'cost', 'utility', 'run_number', 'simulation_time'])\n",
    "output_df[['utility', 'simulation_time']] = output_df[['utility', 'simulation_time']].astype(float)\n",
    "#print(output_df)\n",
    "\n",
    "output_df.loc[output_df['patient_id'] == 7]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "id": "07a71ee0-ab0f-4ad5-87a1-382a0f60f188",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>patient_id</th>\n",
       "      <th>state</th>\n",
       "      <th>treatment_cycles_rec</th>\n",
       "      <th>simulation_time</th>\n",
       "      <th>cost</th>\n",
       "      <th>utility</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>1874</th>\n",
       "      <td>1</td>\n",
       "      <td>Dead</td>\n",
       "      <td>1</td>\n",
       "      <td>8.673985</td>\n",
       "      <td>7000</td>\n",
       "      <td>0.016624</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2942</th>\n",
       "      <td>2</td>\n",
       "      <td>Dead</td>\n",
       "      <td>1</td>\n",
       "      <td>13.658946</td>\n",
       "      <td>8250</td>\n",
       "      <td>0.026178</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>29203</th>\n",
       "      <td>3</td>\n",
       "      <td>Alive</td>\n",
       "      <td>5</td>\n",
       "      <td>117.065124</td>\n",
       "      <td>32750</td>\n",
       "      <td>0.229181</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1103</th>\n",
       "      <td>4</td>\n",
       "      <td>Dead</td>\n",
       "      <td>1</td>\n",
       "      <td>4.700954</td>\n",
       "      <td>6000</td>\n",
       "      <td>0.009010</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>35665</th>\n",
       "      <td>5</td>\n",
       "      <td>Alive</td>\n",
       "      <td>5</td>\n",
       "      <td>158.807567</td>\n",
       "      <td>46500</td>\n",
       "      <td>0.305424</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>291</th>\n",
       "      <td>6</td>\n",
       "      <td>Dead</td>\n",
       "      <td>1</td>\n",
       "      <td>1.394590</td>\n",
       "      <td>5250</td>\n",
       "      <td>0.002673</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>9387</th>\n",
       "      <td>7</td>\n",
       "      <td>Dead</td>\n",
       "      <td>3</td>\n",
       "      <td>35.615664</td>\n",
       "      <td>13500</td>\n",
       "      <td>0.068259</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>11433</th>\n",
       "      <td>8</td>\n",
       "      <td>Dead</td>\n",
       "      <td>2</td>\n",
       "      <td>42.664439</td>\n",
       "      <td>15250</td>\n",
       "      <td>0.081768</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>39408</th>\n",
       "      <td>9</td>\n",
       "      <td>Alive</td>\n",
       "      <td>5</td>\n",
       "      <td>203.330209</td>\n",
       "      <td>51250</td>\n",
       "      <td>0.397874</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4383</th>\n",
       "      <td>10</td>\n",
       "      <td>Dead</td>\n",
       "      <td>2</td>\n",
       "      <td>18.966729</td>\n",
       "      <td>9500</td>\n",
       "      <td>0.036350</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "       patient_id  state  treatment_cycles_rec  simulation_time   cost  \\\n",
       "1874            1   Dead                     1         8.673985   7000   \n",
       "2942            2   Dead                     1        13.658946   8250   \n",
       "29203           3  Alive                     5       117.065124  32750   \n",
       "1103            4   Dead                     1         4.700954   6000   \n",
       "35665           5  Alive                     5       158.807567  46500   \n",
       "291             6   Dead                     1         1.394590   5250   \n",
       "9387            7   Dead                     3        35.615664  13500   \n",
       "11433           8   Dead                     2        42.664439  15250   \n",
       "39408           9  Alive                     5       203.330209  51250   \n",
       "4383           10   Dead                     2        18.966729   9500   \n",
       "\n",
       "        utility  \n",
       "1874   0.016624  \n",
       "2942   0.026178  \n",
       "29203  0.229181  \n",
       "1103   0.009010  \n",
       "35665  0.305424  \n",
       "291    0.002673  \n",
       "9387   0.068259  \n",
       "11433  0.081768  \n",
       "39408  0.397874  \n",
       "4383   0.036350  "
      ]
     },
     "execution_count": 16,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Transform the dataframe so that one row equals one patient where:\n",
    "# Costs and utility are summed, simulation_time equals the latest simulation_time,\n",
    "# State equals the final value of state, and treatment cycle equals the maximum value\n",
    "\n",
    "#summary_df = output_df.loc[output_df.groupby('patient_id')['simulation_time'].idxmax()]\n",
    "summary_df = pd.DataFrame()\n",
    "summary_df[['patient_id', 'state', 'treatment_cycles_rec', 'simulation_time']] = output_df[['patient_id', 'state', 'treatment_cycle', 'simulation_time']].loc[output_df.groupby(['patient_id', 'run_number'])['simulation_time'].idxmax()]\n",
    "\n",
    "# Function to sum costs and utility for each patient\n",
    "def calculate_on_group(x):\n",
    "    return pd.Series(x.sum(), index=x.index)\n",
    "\n",
    "summary_df['cost'] = output_df.groupby(['patient_id','run_number'])['cost'].apply(calculate_on_group)\n",
    "summary_df['utility'] = output_df.groupby(['patient_id','run_number'])['utility'].apply(calculate_on_group)\n",
    "\n",
    "summary_df.head(n = 10)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "id": "62e26b26-fa62-4925-bcce-826fc87188ff",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>cost</th>\n",
       "      <th>utility</th>\n",
       "      <th>treatment_cycles_rec</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>count</th>\n",
       "      <td>10000.000000</td>\n",
       "      <td>10000.000000</td>\n",
       "      <td>10000.000000</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>mean</th>\n",
       "      <td>30131.300000</td>\n",
       "      <td>0.212764</td>\n",
       "      <td>3.684300</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>std</th>\n",
       "      <td>17306.918075</td>\n",
       "      <td>0.151946</td>\n",
       "      <td>1.554051</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>min</th>\n",
       "      <td>5000.000000</td>\n",
       "      <td>0.000029</td>\n",
       "      <td>1.000000</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>25%</th>\n",
       "      <td>14250.000000</td>\n",
       "      <td>0.072431</td>\n",
       "      <td>2.000000</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>50%</th>\n",
       "      <td>31250.000000</td>\n",
       "      <td>0.212564</td>\n",
       "      <td>5.000000</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>75%</th>\n",
       "      <td>43500.000000</td>\n",
       "      <td>0.332453</td>\n",
       "      <td>5.000000</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>max</th>\n",
       "      <td>92250.000000</td>\n",
       "      <td>0.688399</td>\n",
       "      <td>5.000000</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "               cost       utility  treatment_cycles_rec\n",
       "count  10000.000000  10000.000000          10000.000000\n",
       "mean   30131.300000      0.212764              3.684300\n",
       "std    17306.918075      0.151946              1.554051\n",
       "min     5000.000000      0.000029              1.000000\n",
       "25%    14250.000000      0.072431              2.000000\n",
       "50%    31250.000000      0.212564              5.000000\n",
       "75%    43500.000000      0.332453              5.000000\n",
       "max    92250.000000      0.688399              5.000000"
      ]
     },
     "execution_count": 17,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "summary_df[['cost', 'utility', 'treatment_cycles_rec']].describe()"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "2d43beff-9ca3-40ef-8196-af4e6b3a623e",
   "metadata": {},
   "source": [
    "# Checking distributions"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "id": "b5605461-f293-4a63-bdf9-f5346b964dc0",
   "metadata": {},
   "outputs": [],
   "source": [
    "# # Just some code to check the values of distributions for the model\n",
    "# import matplotlib.pyplot as plt\n",
    "# import scipy.special as sps \n",
    "\n",
    "# shape = 1.5\n",
    "# scale = 3\n",
    "# x = np.random.gamma(shape, scale, 10000)\n",
    " \n",
    "# count, bins, ignored = plt.hist(x, 50, density=True)\n",
    "# y = bins**(shape-1)*(np.exp(-bins/scale) /  \n",
    "#                      (sps.gamma(shape)*scale**shape))\n",
    "# plt.plot(bins, y, linewidth=2, color='r')  \n",
    "# # plt.show()"
   ]
  }
 ],
 "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.8"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}