diff --git a/CHANGELOG.md b/CHANGELOG.md index 3b72411..6913b8a 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -5,6 +5,14 @@ All notable changes to this project will be documented in this file. The format is based on [Keep a Changelog](https://keepachangelog.com/en/1.1.0/), and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0.html). Dates formatted as YYYY-MM-DD as per [ISO standard](https://www.iso.org/iso-8601-date-and-time-format.html). +## Unreleased + +## Added + +* `10_multiple_arrival_processes.ipynb`: simulating multiple arrivals classes each with their own distribution +* `11_blocking.ipynb`: simulating the blocking of one resource while queuing for another. +* `distributions.py`: module containing some distributions to reduce code in notebooks. + ## [v0.2.0 - 11/02/2024](https://github.com/pythonhealthdatascience/intro-open-sim/releases/tag/v0.2.0) [![DOI](https://zenodo.org/badge/DOI/10.5281/zenodo.14849934.svg)](https://doi.org/10.5281/zenodo.14849934) Simulation workshop 2025 release. Additional notebooks to respond to reviewers requests. diff --git a/content/11_blocking.ipynb b/content/11_blocking.ipynb new file mode 100644 index 0000000..2343aab --- /dev/null +++ b/content/11_blocking.ipynb @@ -0,0 +1,796 @@ +{ + "cells": [ + { + "cell_type": "markdown", + "id": "0be7dabf-cb34-4faf-abb1-e2c8e735beda", + "metadata": {}, + "source": [ + "# Sequential Resource Holding and Blocking in SimPy\n", + "\n", + "In this notebook we will learn how to code simpy logic to mimic a process holding a resource while queuing for another resource. This simulates a blocking scenario in a process.\n", + "\n", + "We will work with a hypothetical stroke pathway where patients undergo acute treatment followed by transfer at a different hospital to undergo rehabilitation. A patient must remain in an acute stroke bed until a rehabilitation bed is free.\n", + "\n", + "> In this example we will not concern ourselves with a warm-up period or initial conditions.\n", + "\n", + "![model image](img/bed_blocking_image.png \"Bed blocking example\")\n" + ] + }, + { + "cell_type": "markdown", + "id": "0d9383eb-420c-49f8-b178-f2fe9e6b2a90", + "metadata": {}, + "source": [ + "## 1. Imports " + ] + }, + { + "cell_type": "code", + "execution_count": 1, + "id": "c1cee9f9-8696-4b13-94ff-bee2a2a2e5f8", + "metadata": {}, + "outputs": [], + "source": [ + "import numpy as np\n", + "import itertools\n", + "import simpy" + ] + }, + { + "cell_type": "code", + "execution_count": 2, + "id": "ea3d507f-9e6d-4ff0-8b90-f9c63c8a8bdf", + "metadata": {}, + "outputs": [], + "source": [ + "# to reduce code these classes can be found in distribution.py\n", + "from distributions import Exponential, Lognormal" + ] + }, + { + "cell_type": "markdown", + "id": "c422046d-488a-4743-8ad4-97e9f3dab420", + "metadata": {}, + "source": [ + "## 2. Constants" + ] + }, + { + "cell_type": "code", + "execution_count": 3, + "id": "1ecf0429-f03f-4ad2-abb4-46692a74e559", + "metadata": {}, + "outputs": [], + "source": [ + "# default mean inter-arrival times(exp)\n", + "# time is in days\n", + "IAT_STROKES = 1.5\n", + "\n", + "# resources\n", + "N_ACUTE_BEDS = 9\n", + "N_REHAB_BEDS = 15\n", + "\n", + "# Acute LoS (Lognormal)\n", + "ACUTE_LOS_MEAN = 7.0\n", + "ACUTE_LOC_STD = 1.0\n", + "\n", + "# Rehab LoS (Lognormal)\n", + "REHAB_LOS_MEAN = 30.0\n", + "REHAB_LOC_STD = 5.0\n", + "\n", + "# sampling settings\n", + "N_STREAMS = 3\n", + "DEFAULT_RND_SET = 0\n", + "\n", + "# Boolean switch to simulation results as the model runs\n", + "TRACE = False\n", + "\n", + "# run variables (units = days)\n", + "RUN_LENGTH = 100" + ] + }, + { + "cell_type": "markdown", + "id": "5f2a4ad9-6d5e-480d-850f-84d4882a738b", + "metadata": {}, + "source": [ + "## 2. Helper classes and functions" + ] + }, + { + "cell_type": "code", + "execution_count": 4, + "id": "52c9271f-1d05-454d-a199-8768bdf5b6e8", + "metadata": {}, + "outputs": [], + "source": [ + "def trace(msg):\n", + " \"\"\"\n", + " Turing printing of events on and off.\n", + "\n", + " Params:\n", + " -------\n", + " msg: str\n", + " string to print to screen.\n", + " \"\"\"\n", + " if TRACE:\n", + " print(msg)" + ] + }, + { + "cell_type": "markdown", + "id": "5a8c050c-4bb6-408f-a805-3a4aaab56916", + "metadata": {}, + "source": [ + "## 3. Experiment class" + ] + }, + { + "cell_type": "code", + "execution_count": 5, + "id": "576ae9b4-b21b-4ed0-9b13-e5898d423173", + "metadata": {}, + "outputs": [], + "source": [ + "class Experiment:\n", + " \"\"\"\n", + " Encapsulates the concept of an experiment 🧪 for the stroke pathway\n", + " bed blocking simulator. Manages parameters, PRNG streams and results.\n", + " \"\"\"\n", + " def __init__(\n", + " self,\n", + " random_number_set=DEFAULT_RND_SET,\n", + " n_streams=N_STREAMS,\n", + " iat_strokes=IAT_STROKES,\n", + " acute_los_mean=ACUTE_LOS_MEAN,\n", + " acute_los_std=ACUTE_LOC_STD,\n", + " rehab_los_mean=REHAB_LOS_MEAN,\n", + " rehab_los_std=REHAB_LOC_STD,\n", + " n_acute_beds=N_ACUTE_BEDS,\n", + " n_rehab_beds=N_REHAB_BEDS\n", + " ):\n", + " \"\"\"\n", + " The init method sets up our defaults.\n", + " \"\"\"\n", + " # sampling\n", + " self.random_number_set = random_number_set\n", + " self.n_streams = n_streams\n", + "\n", + " # store parameters for the run of the model\n", + " self.iat_strokes = iat_strokes\n", + " self.acute_los_mean = acute_los_mean\n", + " self.acute_los_std = acute_los_std\n", + " self.rehab_los_mean = rehab_los_mean\n", + " self.rehab_los_std = rehab_los_std\n", + "\n", + " # place holder for resources\n", + " self.acute_ward = None\n", + " self.rehab_unit = None\n", + " self.n_acute_beds = n_acute_beds\n", + " self.n_rehab_beds = n_rehab_beds\n", + " \n", + " # initialise results to zero\n", + " self.init_results_variables()\n", + "\n", + " # initialise sampling objects\n", + " self.init_sampling()\n", + "\n", + " def set_random_no_set(self, random_number_set):\n", + " \"\"\"\n", + " Controls the random sampling\n", + " Parameters:\n", + " ----------\n", + " random_number_set: int\n", + " Used to control the set of pseudo random numbers used by\n", + " the distributions in the simulation.\n", + " \"\"\"\n", + " self.random_number_set = random_number_set\n", + " self.init_sampling()\n", + "\n", + " def init_sampling(self):\n", + " \"\"\"\n", + " Create the distributions used by the model and initialise\n", + " the random seeds of each.\n", + " \"\"\"\n", + " # produce n non-overlapping streams\n", + " seed_sequence = np.random.SeedSequence(self.random_number_set)\n", + " self.seeds = seed_sequence.spawn(self.n_streams)\n", + "\n", + " # create distributions\n", + "\n", + " # inter-arrival time distributions\n", + " self.arrival_strokes = Exponential(\n", + " self.iat_strokes, random_seed=self.seeds[0]\n", + " )\n", + "\n", + " self.acute_los = Lognormal(\n", + " self.acute_los_mean, self.acute_los_std, random_seed=self.seeds[1]\n", + " )\n", + "\n", + " self.rehab_los = Lognormal(\n", + " self.rehab_los_mean, self.rehab_los_std, random_seed=self.seeds[2]\n", + " )\n", + "\n", + " def init_results_variables(self):\n", + " \"\"\"\n", + " Initialise all of the experiment variables used in results\n", + " collection. This method is called at the start of each run\n", + " of the model\n", + " \"\"\"\n", + " # variable used to store results of experiment\n", + " self.results = {}\n", + " self.results[\"n_arrivals\"] = 0\n", + " self.results[\"waiting_acute\"] = []\n", + " self.results[\"bed_blocking_times\"] = []" + ] + }, + { + "cell_type": "markdown", + "id": "94f0f9c5-22cb-493a-9f1f-4e2a8325beaa", + "metadata": {}, + "source": [ + "## 4. Pathway process logic\n", + "\n", + "The key things to recognise are \n", + "\n", + "* We request the bed from the acute stroke unit as usual using a `with` context manager\n", + "* We request the rehab bed within the acute bed `with` context manager. This means we do not release the acute bed while the patient waits for rehab.\n", + "* As we do not use a `with` context manager for rehab there is no teardown and we need to manually release the rehab bed." + ] + }, + { + "cell_type": "code", + "execution_count": 6, + "id": "bf2fae82-665f-4990-bc6c-4e43eb68bbeb", + "metadata": {}, + "outputs": [], + "source": [ + "def patient_pathway(patient_id, env, args):\n", + " \"\"\"Process a patient through the acute ward and rehab unit.\n", + " Simpy generator function.\n", + " \n", + " Parameters:\n", + " -----------\n", + " patient_id: int\n", + " A unique id representing the patient in the process\n", + "\n", + " env: simpy.Environment\n", + " The simulation environment\n", + "\n", + " args: Experiment\n", + " Container class for the simulation parameters/results.\n", + " \"\"\"\n", + " arrival_time = env.now\n", + "\n", + " with args.acute_ward.request() as acute_bed_request:\n", + " yield acute_bed_request\n", + " \n", + " acute_admit_time = env.now\n", + " wait_for_acute = acute_admit_time - arrival_time\n", + " args.results['waiting_acute'].append(wait_for_acute)\n", + " \n", + " trace(f\"{env.now:.2f}: Patient {patient_id} admitted to acute ward.\" \\\n", + " + f\"(waited {wait_for_acute:.2f} days)\")\n", + " \n", + " # Simulate acute care treatment\n", + " acute_care_duration = args.acute_los.sample()\n", + " yield env.timeout(acute_care_duration)\n", + " \n", + " # Patient is now medically ready for rehabilitation\n", + " medically_ready_time = env.now\n", + " trace(f\"{env.now:.2f}: Patient {patient_id} medically ready for rehab\")\n", + " \n", + " # Request a rehab bed but don't release the acute bed immediately\n", + " # Note we are still within the \"with\" context manager for the acute bed\n", + " # This is where bed blocking occurs. We wait here until the rehab bed\n", + " # is available. Make sure the indentation is correct or you will release\n", + " rehab_bed = args.rehab_unit.request()\n", + " yield rehab_bed\n", + " \n", + " # Now we have a rehab bed, we can transfer the patient\n", + " transfer_time = env.now\n", + " bed_blocking_duration = transfer_time - medically_ready_time\n", + " args.results['bed_blocking_times'].append(bed_blocking_duration)\n", + " \n", + " trace(f\"{env.now:.2f}: Patient {patient_id} transferred to rehab. \"\\\n", + " + f\"(blocked acute bed for {bed_blocking_duration:.2f} days)\")\n", + " \n", + " # Acute bed is now released\n", + " # Note the indentation! We are now outside of the with context manager.\n", + " # This automatically releases the simpy resource.\n", + " \n", + " # Simulate rehabilitation stay\n", + " rehab_duration = args.rehab_los.sample()\n", + " yield env.timeout(rehab_duration)\n", + " \n", + " # Patient completes rehabilitation and is discharged\n", + " discharge_time = env.now\n", + "\n", + " # Note: we need to explicitly call release on the rehab resource.\n", + " args.rehab_unit.release(rehab_bed)\n", + "\n", + " trace(f\"{env.now:.2f}: Patient {patient_id} discharged from Rehab.\")" + ] + }, + { + "cell_type": "markdown", + "id": "de8990c2-a330-4c02-ac77-26c30d3e0a41", + "metadata": {}, + "source": [ + "## 4. Arrivals generator\n", + "\n", + "This is a standard arrivals generator. We create stroke arrivals according to their distribution." + ] + }, + { + "cell_type": "code", + "execution_count": 7, + "id": "b3e686ce-5371-4471-a052-b9d43309bc85", + "metadata": {}, + "outputs": [], + "source": [ + "def stroke_arrivals_generator(env, args):\n", + " \"\"\"\n", + " Arrival process for strokes.\n", + "\n", + " Parameters:\n", + " ------\n", + " env: simpy.Environment\n", + " The simpy environment for the simulation\n", + "\n", + " args: Experiment\n", + " The settings and input parameters for the simulation.\n", + " \"\"\"\n", + " # use itertools as it provides an infinite loop\n", + " # with a counter variable that we can use for unique Ids\n", + " for patient_id in itertools.count(start=1):\n", + "\n", + " # the sample distribution is defined by the experiment.\n", + " inter_arrival_time = args.arrival_strokes.sample()\n", + " yield env.timeout(inter_arrival_time)\n", + "\n", + " args.results[\"n_arrivals\"] = patient_id\n", + " \n", + " trace(f\"{env.now:.2f}: Stroke arrival.\")\n", + "\n", + " # patient enters pathway\n", + " env.process(patient_pathway(patient_id, env, args))" + ] + }, + { + "cell_type": "markdown", + "id": "6058571e-9fdb-4961-be27-8a3b8c2fe26e", + "metadata": {}, + "source": [ + "## 5. Single run function" + ] + }, + { + "cell_type": "code", + "execution_count": 8, + "id": "0d0ea6cf-7d95-4d2c-9690-fcdbdae35d84", + "metadata": {}, + "outputs": [], + "source": [ + "def single_run(\n", + " experiment, \n", + " rep=0,\n", + " run_length=RUN_LENGTH\n", + "):\n", + " \"\"\"\n", + " Perform a single run of the model and return the results\n", + "\n", + " Parameters:\n", + " -----------\n", + "\n", + " experiment: Experiment\n", + " The experiment/paramaters to use with model\n", + "\n", + " rep: int\n", + " The replication number.\n", + "\n", + " rc_period: float, optional (default=RUN_LENGTH)\n", + " The run length of the model\n", + " \"\"\"\n", + "\n", + " # reset all results variables to zero and empty\n", + " experiment.init_results_variables()\n", + "\n", + " # set random number set to the replication no.\n", + " # this controls sampling for the run.\n", + " experiment.set_random_no_set(rep)\n", + "\n", + " # environment is (re)created inside single run\n", + " env = simpy.Environment()\n", + "\n", + " # simpy resources\n", + " experiment.acute_ward = simpy.Resource(env, experiment.n_acute_beds)\n", + " experiment.rehab_unit = simpy.Resource(env, experiment.n_rehab_beds)\n", + "\n", + " # we pass all arrival generators to simpy \n", + " env.process(stroke_arrivals_generator(env, experiment))\n", + "\n", + " # run model\n", + " env.run(until=run_length)\n", + "\n", + " # quick stats\n", + " results = {}\n", + " results['mean_acute_wait'] = np.array(\n", + " experiment.results[\"waiting_acute\"]\n", + " ).mean()\n", + " results['mean_bed_blocking'] = np.array(\n", + " experiment.results[\"bed_blocking_times\"]\n", + " ).mean()\n", + "\n", + " # return single run results\n", + " return results" + ] + }, + { + "cell_type": "code", + "execution_count": 9, + "id": "caf52390-5455-4fa1-bb22-60b5b91ad8d0", + "metadata": {}, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "4.94: Stroke arrival.\n", + "4.94: Patient 1 admitted to acute ward.(waited 0.00 days)\n", + "6.08: Stroke arrival.\n", + "6.08: Patient 2 admitted to acute ward.(waited 0.00 days)\n", + "7.96: Stroke arrival.\n", + "7.96: Patient 3 admitted to acute ward.(waited 0.00 days)\n", + "8.30: Stroke arrival.\n", + "8.30: Patient 4 admitted to acute ward.(waited 0.00 days)\n", + "8.63: Stroke arrival.\n", + "8.63: Patient 5 admitted to acute ward.(waited 0.00 days)\n", + "11.33: Stroke arrival.\n", + "11.33: Patient 6 admitted to acute ward.(waited 0.00 days)\n", + "11.37: Patient 2 medically ready for rehab\n", + "11.37: Patient 2 transferred to rehab. (blocked acute bed for 0.00 days)\n", + "11.43: Stroke arrival.\n", + "11.43: Patient 7 admitted to acute ward.(waited 0.00 days)\n", + "12.18: Patient 3 medically ready for rehab\n", + "12.18: Patient 3 transferred to rehab. (blocked acute bed for 0.00 days)\n", + "12.71: Patient 1 medically ready for rehab\n", + "12.71: Patient 1 transferred to rehab. (blocked acute bed for 0.00 days)\n", + "12.78: Stroke arrival.\n", + "12.78: Patient 8 admitted to acute ward.(waited 0.00 days)\n", + "13.23: Stroke arrival.\n", + "13.23: Patient 9 admitted to acute ward.(waited 0.00 days)\n", + "16.22: Patient 4 medically ready for rehab\n", + "16.22: Patient 4 transferred to rehab. (blocked acute bed for 0.00 days)\n", + "16.68: Stroke arrival.\n", + "16.68: Patient 10 admitted to acute ward.(waited 0.00 days)\n", + "16.96: Patient 5 medically ready for rehab\n", + "16.96: Patient 5 transferred to rehab. (blocked acute bed for 0.00 days)\n", + "17.93: Patient 6 medically ready for rehab\n", + "17.93: Patient 6 transferred to rehab. (blocked acute bed for 0.00 days)\n", + "18.68: Patient 7 medically ready for rehab\n", + "18.68: Patient 7 transferred to rehab. (blocked acute bed for 0.00 days)\n", + "18.76: Patient 8 medically ready for rehab\n", + "18.76: Patient 8 transferred to rehab. (blocked acute bed for 0.00 days)\n", + "19.46: Stroke arrival.\n", + "19.46: Patient 11 admitted to acute ward.(waited 0.00 days)\n", + "19.49: Stroke arrival.\n", + "19.49: Patient 12 admitted to acute ward.(waited 0.00 days)\n", + "19.74: Patient 9 medically ready for rehab\n", + "19.74: Patient 9 transferred to rehab. (blocked acute bed for 0.00 days)\n", + "24.39: Patient 10 medically ready for rehab\n", + "24.39: Patient 10 transferred to rehab. (blocked acute bed for 0.00 days)\n", + "27.41: Patient 11 medically ready for rehab\n", + "27.41: Patient 11 transferred to rehab. (blocked acute bed for 0.00 days)\n", + "27.70: Patient 12 medically ready for rehab\n", + "27.70: Patient 12 transferred to rehab. (blocked acute bed for 0.00 days)\n", + "27.88: Stroke arrival.\n", + "27.88: Patient 13 admitted to acute ward.(waited 0.00 days)\n", + "29.08: Stroke arrival.\n", + "29.08: Patient 14 admitted to acute ward.(waited 0.00 days)\n", + "29.60: Stroke arrival.\n", + "29.60: Patient 15 admitted to acute ward.(waited 0.00 days)\n", + "32.56: Stroke arrival.\n", + "32.56: Patient 16 admitted to acute ward.(waited 0.00 days)\n", + "33.67: Stroke arrival.\n", + "33.67: Patient 17 admitted to acute ward.(waited 0.00 days)\n", + "33.68: Stroke arrival.\n", + "33.68: Patient 18 admitted to acute ward.(waited 0.00 days)\n", + "34.92: Stroke arrival.\n", + "34.92: Patient 19 admitted to acute ward.(waited 0.00 days)\n", + "35.09: Patient 13 medically ready for rehab\n", + "35.09: Patient 13 transferred to rehab. (blocked acute bed for 0.00 days)\n", + "35.85: Stroke arrival.\n", + "35.85: Patient 20 admitted to acute ward.(waited 0.00 days)\n", + "36.02: Patient 15 medically ready for rehab\n", + "36.02: Patient 15 transferred to rehab. (blocked acute bed for 0.00 days)\n", + "36.62: Patient 14 medically ready for rehab\n", + "36.62: Patient 14 transferred to rehab. (blocked acute bed for 0.00 days)\n", + "36.76: Stroke arrival.\n", + "36.76: Patient 21 admitted to acute ward.(waited 0.00 days)\n", + "38.54: Patient 3 discharged from Rehab.\n", + "38.56: Patient 5 discharged from Rehab.\n", + "39.04: Patient 16 medically ready for rehab\n", + "39.04: Patient 16 transferred to rehab. (blocked acute bed for 0.00 days)\n", + "39.97: Patient 18 medically ready for rehab\n", + "39.97: Patient 18 transferred to rehab. (blocked acute bed for 0.00 days)\n", + "40.33: Patient 17 medically ready for rehab\n", + "41.02: Stroke arrival.\n", + "41.02: Patient 22 admitted to acute ward.(waited 0.00 days)\n", + "41.33: Patient 19 medically ready for rehab\n", + "41.37: Stroke arrival.\n", + "41.37: Patient 23 admitted to acute ward.(waited 0.00 days)\n", + "41.50: Patient 9 discharged from Rehab.\n", + "41.50: Patient 17 transferred to rehab. (blocked acute bed for 1.17 days)\n", + "42.66: Patient 20 medically ready for rehab\n", + "42.95: Patient 21 medically ready for rehab\n", + "45.13: Stroke arrival.\n", + "45.13: Patient 24 admitted to acute ward.(waited 0.00 days)\n", + "45.62: Stroke arrival.\n", + "45.62: Patient 25 admitted to acute ward.(waited 0.00 days)\n", + "45.95: Patient 2 discharged from Rehab.\n", + "45.95: Patient 19 transferred to rehab. (blocked acute bed for 4.62 days)\n", + "47.94: Stroke arrival.\n", + "47.94: Patient 26 admitted to acute ward.(waited 0.00 days)\n", + "49.03: Patient 23 medically ready for rehab\n", + "49.52: Stroke arrival.\n", + "49.52: Patient 27 admitted to acute ward.(waited 0.00 days)\n", + "49.94: Patient 22 medically ready for rehab\n", + "50.11: Stroke arrival.\n", + "50.11: Patient 28 admitted to acute ward.(waited 0.00 days)\n", + "50.81: Patient 7 discharged from Rehab.\n", + "50.81: Patient 20 transferred to rehab. (blocked acute bed for 8.16 days)\n", + "50.88: Patient 4 discharged from Rehab.\n", + "50.88: Patient 21 transferred to rehab. (blocked acute bed for 7.93 days)\n", + "51.13: Stroke arrival.\n", + "51.13: Patient 29 admitted to acute ward.(waited 0.00 days)\n", + "51.17: Patient 8 discharged from Rehab.\n", + "51.17: Patient 23 transferred to rehab. (blocked acute bed for 2.14 days)\n", + "51.20: Patient 1 discharged from Rehab.\n", + "51.20: Patient 22 transferred to rehab. (blocked acute bed for 1.27 days)\n", + "51.23: Patient 10 discharged from Rehab.\n", + "52.55: Patient 6 discharged from Rehab.\n", + "53.81: Patient 24 medically ready for rehab\n", + "53.81: Patient 24 transferred to rehab. (blocked acute bed for 0.00 days)\n", + "54.05: Stroke arrival.\n", + "54.05: Patient 30 admitted to acute ward.(waited 0.00 days)\n", + "54.30: Patient 26 medically ready for rehab\n", + "54.30: Patient 26 transferred to rehab. (blocked acute bed for 0.00 days)\n", + "54.82: Patient 25 medically ready for rehab\n", + "56.20: Patient 28 medically ready for rehab\n", + "56.79: Patient 27 medically ready for rehab\n", + "57.60: Patient 11 discharged from Rehab.\n", + "57.60: Patient 25 transferred to rehab. (blocked acute bed for 2.78 days)\n", + "57.62: Stroke arrival.\n", + "57.62: Patient 31 admitted to acute ward.(waited 0.00 days)\n", + "58.41: Patient 29 medically ready for rehab\n", + "58.87: Stroke arrival.\n", + "58.87: Patient 32 admitted to acute ward.(waited 0.00 days)\n", + "60.28: Patient 30 medically ready for rehab\n", + "62.00: Stroke arrival.\n", + "62.00: Patient 33 admitted to acute ward.(waited 0.00 days)\n", + "63.04: Stroke arrival.\n", + "63.04: Patient 34 admitted to acute ward.(waited 0.00 days)\n", + "63.74: Stroke arrival.\n", + "63.74: Patient 35 admitted to acute ward.(waited 0.00 days)\n", + "64.60: Patient 31 medically ready for rehab\n", + "65.55: Stroke arrival.\n", + "65.68: Stroke arrival.\n", + "66.18: Stroke arrival.\n", + "66.23: Patient 12 discharged from Rehab.\n", + "66.23: Patient 28 transferred to rehab. (blocked acute bed for 10.02 days)\n", + "66.23: Patient 36 admitted to acute ward.(waited 0.68 days)\n", + "66.40: Patient 32 medically ready for rehab\n", + "67.69: Stroke arrival.\n", + "68.04: Stroke arrival.\n", + "68.77: Patient 13 discharged from Rehab.\n", + "68.77: Patient 27 transferred to rehab. (blocked acute bed for 11.98 days)\n", + "68.77: Patient 37 admitted to acute ward.(waited 3.09 days)\n", + "68.84: Patient 34 medically ready for rehab\n", + "68.98: Patient 33 medically ready for rehab\n", + "69.69: Patient 16 discharged from Rehab.\n", + "69.69: Patient 29 transferred to rehab. (blocked acute bed for 11.29 days)\n", + "69.69: Patient 38 admitted to acute ward.(waited 3.51 days)\n", + "69.83: Stroke arrival.\n", + "70.54: Patient 35 medically ready for rehab\n", + "70.90: Stroke arrival.\n", + "71.72: Patient 19 discharged from Rehab.\n", + "71.72: Patient 30 transferred to rehab. (blocked acute bed for 11.44 days)\n", + "71.72: Patient 39 admitted to acute ward.(waited 4.03 days)\n", + "72.76: Patient 22 discharged from Rehab.\n", + "72.76: Patient 31 transferred to rehab. (blocked acute bed for 8.15 days)\n", + "72.76: Patient 40 admitted to acute ward.(waited 4.71 days)\n", + "72.79: Patient 36 medically ready for rehab\n", + "73.62: Patient 15 discharged from Rehab.\n", + "73.62: Patient 32 transferred to rehab. (blocked acute bed for 7.22 days)\n", + "73.62: Patient 41 admitted to acute ward.(waited 3.79 days)\n", + "73.65: Patient 18 discharged from Rehab.\n", + "73.65: Patient 34 transferred to rehab. (blocked acute bed for 4.81 days)\n", + "73.65: Patient 42 admitted to acute ward.(waited 2.74 days)\n", + "73.78: Patient 17 discharged from Rehab.\n", + "73.78: Patient 33 transferred to rehab. (blocked acute bed for 4.80 days)\n", + "75.99: Patient 37 medically ready for rehab\n", + "76.13: Patient 38 medically ready for rehab\n", + "76.53: Patient 39 medically ready for rehab\n", + "77.23: Stroke arrival.\n", + "77.23: Patient 43 admitted to acute ward.(waited 0.00 days)\n", + "77.29: Stroke arrival.\n", + "77.37: Patient 14 discharged from Rehab.\n", + "77.37: Patient 35 transferred to rehab. (blocked acute bed for 6.82 days)\n", + "77.37: Patient 44 admitted to acute ward.(waited 0.07 days)\n", + "77.49: Patient 20 discharged from Rehab.\n", + "77.49: Patient 36 transferred to rehab. (blocked acute bed for 4.70 days)\n", + "78.96: Patient 40 medically ready for rehab\n", + "79.28: Patient 21 discharged from Rehab.\n", + "79.28: Patient 37 transferred to rehab. (blocked acute bed for 3.29 days)\n", + "80.48: Patient 41 medically ready for rehab\n", + "80.65: Patient 42 medically ready for rehab\n", + "81.76: Patient 26 discharged from Rehab.\n", + "81.76: Patient 38 transferred to rehab. (blocked acute bed for 5.63 days)\n", + "81.86: Stroke arrival.\n", + "81.86: Patient 45 admitted to acute ward.(waited 0.00 days)\n", + "82.18: Stroke arrival.\n", + "82.18: Patient 46 admitted to acute ward.(waited 0.00 days)\n", + "82.30: Stroke arrival.\n", + "82.30: Patient 47 admitted to acute ward.(waited 0.00 days)\n", + "82.43: Patient 23 discharged from Rehab.\n", + "82.43: Patient 39 transferred to rehab. (blocked acute bed for 5.91 days)\n", + "83.69: Stroke arrival.\n", + "83.69: Patient 48 admitted to acute ward.(waited 0.00 days)\n", + "84.18: Patient 44 medically ready for rehab\n", + "85.15: Patient 24 discharged from Rehab.\n", + "85.15: Patient 40 transferred to rehab. (blocked acute bed for 6.19 days)\n", + "85.18: Patient 43 medically ready for rehab\n", + "86.77: Stroke arrival.\n", + "86.77: Patient 49 admitted to acute ward.(waited 0.00 days)\n", + "87.16: Patient 46 medically ready for rehab\n", + "87.85: Stroke arrival.\n", + "88.27: Stroke arrival.\n", + "89.09: Stroke arrival.\n", + "89.10: Patient 45 medically ready for rehab\n", + "89.38: Patient 47 medically ready for rehab\n", + "89.46: Stroke arrival.\n", + "90.62: Patient 48 medically ready for rehab\n", + "92.33: Stroke arrival.\n", + "92.74: Patient 25 discharged from Rehab.\n", + "92.74: Patient 41 transferred to rehab. (blocked acute bed for 12.26 days)\n", + "92.74: Patient 50 admitted to acute ward.(waited 4.88 days)\n", + "92.92: Patient 27 discharged from Rehab.\n", + "92.92: Patient 42 transferred to rehab. (blocked acute bed for 12.28 days)\n", + "92.92: Patient 51 admitted to acute ward.(waited 4.66 days)\n", + "92.99: Stroke arrival.\n", + "93.48: Stroke arrival.\n", + "93.86: Patient 28 discharged from Rehab.\n", + "93.86: Patient 44 transferred to rehab. (blocked acute bed for 9.68 days)\n", + "93.86: Patient 52 admitted to acute ward.(waited 4.76 days)\n", + "94.18: Stroke arrival.\n", + "94.27: Patient 49 medically ready for rehab\n", + "94.96: Patient 30 discharged from Rehab.\n", + "94.96: Patient 43 transferred to rehab. (blocked acute bed for 9.78 days)\n", + "94.96: Patient 53 admitted to acute ward.(waited 5.51 days)\n", + "95.92: Patient 29 discharged from Rehab.\n", + "95.92: Patient 46 transferred to rehab. (blocked acute bed for 8.76 days)\n", + "95.92: Patient 54 admitted to acute ward.(waited 3.59 days)\n", + "96.06: Stroke arrival.\n", + "97.43: Stroke arrival.\n", + "97.56: Stroke arrival.\n", + "98.97: Stroke arrival.\n", + "99.31: Patient 33 discharged from Rehab.\n", + "99.31: Patient 45 transferred to rehab. (blocked acute bed for 10.21 days)\n", + "99.31: Patient 55 admitted to acute ward.(waited 6.32 days)\n", + "99.48: Patient 34 discharged from Rehab.\n", + "99.48: Patient 47 transferred to rehab. (blocked acute bed for 10.10 days)\n", + "99.48: Patient 56 admitted to acute ward.(waited 5.99 days)\n", + "99.69: Patient 51 medically ready for rehab\n", + "99.97: Patient 52 medically ready for rehab\n" + ] + }, + { + "data": { + "text/plain": [ + "{'mean_acute_wait': 1.041965345852034, 'mean_bed_blocking': 4.3271192681703115}" + ] + }, + "execution_count": 9, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "TRACE = True\n", + "experiment = Experiment()\n", + "results = single_run(experiment)\n", + "results" + ] + }, + { + "cell_type": "markdown", + "id": "a48ffebd-5af0-4354-89bc-7de77ee60e8b", + "metadata": {}, + "source": [ + "## Quick tests\n", + "\n", + "### Remove rehab bed blocking" + ] + }, + { + "cell_type": "code", + "execution_count": 10, + "id": "8228ab0c-cc99-48e2-a2c9-c9dcce8d854f", + "metadata": {}, + "outputs": [ + { + "data": { + "text/plain": [ + "{'mean_acute_wait': 0.0, 'mean_bed_blocking': 0.0}" + ] + }, + "execution_count": 10, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "M = 1_000_000\n", + "TRACE = False\n", + "experiment = Experiment(n_rehab_beds=M)\n", + "results = single_run(experiment)\n", + "results" + ] + }, + { + "cell_type": "markdown", + "id": "f8e498cb-3999-4c69-8d41-e436cbebe6d8", + "metadata": {}, + "source": [ + "### Remove all waiting time for acute unit" + ] + }, + { + "cell_type": "code", + "execution_count": 11, + "id": "00d319d8-2a15-46f1-b1e9-a8efb7a69b2a", + "metadata": {}, + "outputs": [ + { + "data": { + "text/plain": [ + "{'mean_acute_wait': 0.0, 'mean_bed_blocking': 4.80862454209246}" + ] + }, + "execution_count": 11, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "M = 1_000_000\n", + "TRACE = False\n", + "experiment = Experiment(n_acute_beds=M)\n", + "results = single_run(experiment)\n", + "results" + ] + } + ], + "metadata": { + "kernelspec": { + "display_name": "Python 3 (ipykernel)", + "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.11.11" + } + }, + "nbformat": 4, + "nbformat_minor": 5 +} diff --git a/content/distributions.py b/content/distributions.py new file mode 100644 index 0000000..f8f42e0 --- /dev/null +++ b/content/distributions.py @@ -0,0 +1,143 @@ +""" +Distribution and utility classes avilable to these examples. +""" + +import numpy as np +from typing import Optional, Union, Tuple +from numpy.random import SeedSequence +from numpy.typing import NDArray +import math + +class Lognormal: + """ + Lognormal distribution implementation. + + A continuous probability distribution where the logarithm of a random + variable is normally distributed. It is useful for modeling variables that + are the product of many small independent factors. + + Notes: taken from my in development package sim-tools + https://github.com/TomMonks/sim-tools + """ + + def __init__( + self, + mean: float, + stdev: float, + random_seed: Optional[Union[int, SeedSequence]] = None, + ): + """ + Initialize a lognormal distribution. + + Parameters + ---------- + mean : float + Mean of the lognormal distribution. + + stdev : float + Standard deviation of the lognormal distribution. + + random_seed : Optional[Union[int, SeedSequence]], default=None + A random seed or SeedSequence to reproduce samples. If None, a + unique sample sequence is generated. + """ + self.rng = np.random.default_rng(random_seed) + mu, sigma = self.normal_moments_from_lognormal(mean, stdev**2) + self.mu = mu + self.sigma = sigma + self.mean = mean + self.stdev = stdev + + def __repr__(self): + return f"Lognormal(mean={self.mean}, stdev={self.stdev})" + + def normal_moments_from_lognormal( + self, m: float, v: float + ) -> Tuple[float, float]: + """ + Calculate mu and sigma of the normal distribution underlying + a lognormal with mean m and variance v. + + Parameters + ---------- + m : float + Mean of lognormal distribution. + v : float + Variance of lognormal distribution. + + Returns + ------- + Tuple[float, float] + The mu and sigma parameters of the underlying normal distribution. + + Notes + ----- + Formula source: + https://blogs.sas.com/content/iml/2014/06/04/simulate-lognormal-data- + with-specified-mean-and-variance.html + """ + phi = math.sqrt(v + m**2) + mu = math.log(m**2 / phi) + sigma = math.sqrt(math.log(phi**2 / m**2)) + return mu, sigma + + def sample( + self, size: Optional[Union[int, Tuple[int, ...]]] = None + ) -> Union[float, NDArray[np.float64]]: + """ + Generate random samples from the lognormal distribution. + + Parameters + ---------- + size : Optional[Union[int, Tuple[int, ...]]], default=None + The number/shape of samples to generate: + - If None: returns a single sample as a float + - If int: returns a 1-D array with that many samples + - If tuple of ints: returns an array with that shape + + Returns + ------- + Union[float, NDArray[np.float64]] + Random samples from the lognormal distribution: + - A single float when size is None + - A numpy array of floats with shape determined by size parameter + """ + return self.rng.lognormal(self.mu, self.sigma, size=size) + +class Exponential: + """ + Convenience class for the exponential distribution. + packages up distribution parameters, seed and random generator. + """ + + def __init__(self, mean, random_seed=None): + """ + Constructor + + Params: + ------ + mean: float + The mean of the exponential distribution + + random_seed: int| SeedSequence, optional (default=None) + A random seed to reproduce samples. If set to none then a unique + sample is created. + """ + self.rand = np.random.default_rng(seed=random_seed) + self.mean = mean + + def sample(self, size=None): + """ + Generate a sample from the exponential distribution + + Params: + ------- + size: int, optional (default=None) + the number of samples to return. If size=None then a single + sample is returned. + + Returns: + ------- + float or np.ndarray (if size >=1) + """ + return self.rand.exponential(self.mean, size=size) \ No newline at end of file diff --git a/content/img/bed_blocking_image.png b/content/img/bed_blocking_image.png new file mode 100644 index 0000000..65fcfae Binary files /dev/null and b/content/img/bed_blocking_image.png differ