diff --git a/CHANGES.md b/CHANGES.md index bc791e9..90c0f52 100644 --- a/CHANGES.md +++ b/CHANGES.md @@ -7,7 +7,24 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0 Consistent identifier (represents all versions, resolves to latest): [![DOI](https://zenodo.org/badge/DOI/10.5281/zenodo.10026326.svg)](https://doi.org/10.5281/zenodo.10026326) -## [v2.1.0]() - 2024-05-30 +## [v2.2.0]() UNRELEASED + +### Added + +* All model classes and functions now have python type hints +* `treat_sim.datasets` module with `load_nelson_arrivals`, `load_alternative_arrivals` and `valid_arrival_profile` functions +* `tests/test_datasets.py` contains functional and dirty tests for loading and using internal arrival profile datasets. + +### Changed + +* `Scenario` defaults to the time dependent arrival profile given in Nelson (2013), but also accepts `arrival_profile` a `pandas.DataFrame` parameter for scenario analysis. +* Default arrival profile is sourced from local package rather than GitHub URL. + +### Fixed + +* MODEL: thinning alg: `np.Inf` -> `np.inf` for compatibility with `numpy>=2` + +## [v2.1.0](https://github.com/pythonhealthdatascience/stars-treat-sim/releases/tag/v2.1.0) - 2024-05-30 - [![DOI](https://zenodo.org/badge/DOI/10.5281/zenodo.11396022.svg)](https://doi.org/10.5281/zenodo.11396022) ### Changes diff --git a/README.md b/README.md index b77983d..bae7e9d 100644 --- a/README.md +++ b/README.md @@ -1,26 +1,25 @@ + +# πŸ’« Towards Sharing Tools, and Artifacts, for Reusable Simulation (STARS): a minimal model example + [![Binder](https://mybinder.org/badge_logo.svg)](https://mybinder.org/v2/gh/pythonhealthdatascience/stars-treat-sim/HEAD) [![License: MIT](https://img.shields.io/badge/License-MIT-yellow.svg)](https://opensource.org/licenses/MIT) [![Python](https://img.shields.io/pypi/pyversions/treat_sim)](https://pypi.org/project/treat_sim/) [![DOI](https://zenodo.org/badge/DOI/10.5281/zenodo.10026326.svg)](https://doi.org/10.5281/zenodo.10026326) [![PyPI version fury.io](https://badge.fury.io/py/treat-sim.svg)](https://pypi.org/project/treat-sim/) - - [](https://hub.docker.com/r/tommonks01/treat_sim) - -# Towards Sharing Tools, and Artifacts, for Reusable Simulation: a minimal model example - ## Overview -The materials and methods in this repository support work towards developing the S.T.A.R.S healthcare framework (**S**haring **T**ools and **A**rtifacts for **R**eusable **S**imulations in healthcare). The code and written materials here demonstrate the application of S.T.A.R.S' version 1 to sharing a `SimPy` discrete-event simulation model and associated research artifacts. +The materials and methods in this repository support work towards developing the STARShealthcare framework (**S**haring **T**ools and **A**rtifacts for **R**eusable **S**imulations in healthcare). The code and written materials here demonstrate the application of STARS version 1 to sharing a `SimPy` discrete-event simulation model and associated research artifacts. * All artifacts in this repository are linked to study researchers via ORCIDs; * Model code is made available under an MIT license; -* Python dependencies are managed through `conda`; -* Documentation of the model is enhanced using a Jupyter notebook. -* The python code itself can be viewed and executed in Jupyter notebooks via [Binder](https://mybinder.org); +* Python dependencies are managed through `mamba`; +* Documentation of the model is enhanced using a simple Jupyter notebook. +* The python model itself can be viewed and executed in Jupyter notebooks via [Binder](https://mybinder.org); * The materials are deposited and made citable using Zenodo; * The model is sharable with other researchers and the NHS without the need to install software. +* A full suite of automated tests are provided with the model. ## Author ORCIDs @@ -35,7 +34,7 @@ This code is part of independent research supported by the National Institute fo ### Install from PyPI -If you do not wish to you the code or would like to use the model as part of your own work you can install the model as a python package. +If you do not wish to view the code or would like to use the model as part of your own work you can install the model as a python package. ```bash pip install treat-sim @@ -59,7 +58,7 @@ git clone https://github.com/pythonhealthdatascience/stars-treat-sim #### Installing dependencies -[![Python 3.8+](https://img.shields.io/badge/python-3.8+-blue.svg)](https://www.python.org/downloads/release/python-380/) +[![Python](https://img.shields.io/pypi/pyversions/treat_sim)](https://pypi.org/project/treat_sim/) All dependencies can be found in [`binder/environment.yml`]() and are pulled from conda-forge. To run the code locally, we recommend installing [miniforge](https://github.com/conda-forge/miniforge); @@ -111,6 +110,23 @@ if __name__ == '__main__': print(results) ``` + +The model can be run with different time dependent arrival profiles. By default the model runs with the arrival profile taken from Nelson (2013). The `datasets` module provides access to an alternative example dataset where arrivals are slightly skewed towards the end of the working day. + +```python + +from treat_sim.model import Scenario, multiple_replications +from treat_sim.datasets import load_alternative_arrivals + +if __name__ == '__main__': + + # set the arrival profile to later in the day + scenario1 = Scenario(arrival_porfile=load_alternative_arrivals()) + + alternative_results = multiple_replications(scenario1).describe().round(2).T + print(alternative_results) +``` + #### Testing the model > See our [online documentation](https://pythonhealthdatascience.github.io/stars-simpy-example-docs/content/02_model_code/05_testing.html) for an overview of testing @@ -142,12 +158,15 @@ pytest --cov=treat_sim tests/ β”œβ”€β”€ pyproject.toml β”œβ”€β”€ README.md β”œβ”€β”€ tests +β”‚Β Β  └── test_datasets.ipynb β”‚Β Β  └── test_model.ipynb └── treat_sim β”œβ”€β”€ data β”‚Β Β  └── ed_arrivals.csv - β”œβ”€β”€ distributions.py + β”‚Β Β  └── ed_arrivals_scenario1.csv β”œβ”€β”€ __init__.py + β”œβ”€β”€ datasets.py + β”œβ”€β”€ distributions.py └── model.py ``` @@ -161,8 +180,9 @@ pytest --cov=treat_sim tests/ * `tests/` - contains automated testing code * `treat_sim/` - contains packaged version of the model. * `data/` - directory containing data file used by package. - * `distributions.py` - distribution classes. * `__init__.py` - required as part of package - contains author and version. + * `datasets.py` - functions to load example dataset for parameterising the model. + * `distributions.py` - distribution classes. * `model.py` - example SimPy model. @@ -182,7 +202,7 @@ Monks, T., Harper, A., & Heather, A. (2024). Towards Sharing Tools, and Artifact month = May, year = 2024, publisher = {Zenodo}, - version = {v2.0.0}, + version = {v2.2.0}, doi = {10.5281//zenodo.10026326.}, url = {https://doi.org/10.5281//zenodo.10026326} } diff --git a/notebooks/test_package.ipynb b/notebooks/test_package.ipynb index ede3a82..a2a70b7 100644 --- a/notebooks/test_package.ipynb +++ b/notebooks/test_package.ipynb @@ -19,7 +19,7 @@ { "data": { "text/plain": [ - "'1.0.0'" + "'2.2.0'" ] }, "execution_count": 1, @@ -42,7 +42,7 @@ }, { "cell_type": "code", - "execution_count": null, + "execution_count": 8, "id": "cf77eab0-8038-4e89-b0a8-abae02c61046", "metadata": {}, "outputs": [], @@ -54,10 +54,33 @@ }, { "cell_type": "code", - "execution_count": null, + "execution_count": 14, "id": "a4427a5c-a048-4e78-81ad-9a72a62b982f", "metadata": {}, - "outputs": [], + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "Scenario Analysis\n", + "No. Scenario: 6\n", + "Replications: 100\n", + "Running base => done.\n", + "\n", + "Running triage+1 => done.\n", + "\n", + "Running exam+1 => done.\n", + "\n", + "Running treat+1 => done.\n", + "\n", + "Running swap_exam_treat => done.\n", + "\n", + "Running short_exam => done.\n", + "\n", + "Scenario analysis complete.\n" + ] + } + ], "source": [ "results = run_scenario_analysis(get_scenarios(), DEFAULT_RESULTS_COLLECTION_PERIOD,\n", " n_reps=50)" @@ -65,10 +88,231 @@ }, { "cell_type": "code", - "execution_count": null, + "execution_count": 10, "id": "42d1d046-102e-4edc-8f22-98ea40429e27", "metadata": {}, - "outputs": [], + "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", + " \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", + "
basetriage+1exam+1treat+1swap_exam_treatshort_exam
00_arrivals227.720000227.720000227.720000227.720000227.720000227.720000
01a_triage_wait35.2430191.32517635.24301935.2430191.3251761.325176
01b_triage_util0.6070390.3035190.6070390.6070390.3035190.303519
02a_registration_wait105.572644135.915119105.572644105.572644135.915119135.915119
02b_registration_util0.8406130.8439540.8406130.8406130.8439540.843954
03a_examination_wait25.55252526.6888660.14656525.552525148.86596566.003847
03b_examination_util0.8504480.8535670.6678560.8504480.8936290.868998
04a_treatment_wait(non_trauma)136.660635138.575405151.8082042.28096233.568594111.337103
04b_treatment_util(non_trauma)0.8672380.8693240.8676910.6293460.8363840.871039
05_total_time(non-trauma)234.335783232.133275223.063821192.191909280.729440244.797808
06a_trauma_wait151.680557175.135269151.680557151.680557175.135269175.135269
06b_trauma_util0.8269410.8409770.8269410.8269410.8409770.840977
07a_treatment_wait(trauma)14.30504414.99728714.30504414.30504414.99728714.997287
07b_treatment_util(trauma)0.5010100.5069080.5010100.5010100.5069080.506908
08_total_time(trauma)292.280487282.783517292.280487292.280487282.783517282.783517
09_throughput162.160000162.920000165.880000195.300000138.540000156.600000
\n", + "
" + ], + "text/plain": [ + " base triage+1 exam+1 \\\n", + "00_arrivals 227.720000 227.720000 227.720000 \n", + "01a_triage_wait 35.243019 1.325176 35.243019 \n", + "01b_triage_util 0.607039 0.303519 0.607039 \n", + "02a_registration_wait 105.572644 135.915119 105.572644 \n", + "02b_registration_util 0.840613 0.843954 0.840613 \n", + "03a_examination_wait 25.552525 26.688866 0.146565 \n", + "03b_examination_util 0.850448 0.853567 0.667856 \n", + "04a_treatment_wait(non_trauma) 136.660635 138.575405 151.808204 \n", + "04b_treatment_util(non_trauma) 0.867238 0.869324 0.867691 \n", + "05_total_time(non-trauma) 234.335783 232.133275 223.063821 \n", + "06a_trauma_wait 151.680557 175.135269 151.680557 \n", + "06b_trauma_util 0.826941 0.840977 0.826941 \n", + "07a_treatment_wait(trauma) 14.305044 14.997287 14.305044 \n", + "07b_treatment_util(trauma) 0.501010 0.506908 0.501010 \n", + "08_total_time(trauma) 292.280487 282.783517 292.280487 \n", + "09_throughput 162.160000 162.920000 165.880000 \n", + "\n", + " treat+1 swap_exam_treat short_exam \n", + "00_arrivals 227.720000 227.720000 227.720000 \n", + "01a_triage_wait 35.243019 1.325176 1.325176 \n", + "01b_triage_util 0.607039 0.303519 0.303519 \n", + "02a_registration_wait 105.572644 135.915119 135.915119 \n", + "02b_registration_util 0.840613 0.843954 0.843954 \n", + "03a_examination_wait 25.552525 148.865965 66.003847 \n", + "03b_examination_util 0.850448 0.893629 0.868998 \n", + "04a_treatment_wait(non_trauma) 2.280962 33.568594 111.337103 \n", + "04b_treatment_util(non_trauma) 0.629346 0.836384 0.871039 \n", + "05_total_time(non-trauma) 192.191909 280.729440 244.797808 \n", + "06a_trauma_wait 151.680557 175.135269 175.135269 \n", + "06b_trauma_util 0.826941 0.840977 0.840977 \n", + "07a_treatment_wait(trauma) 14.305044 14.997287 14.997287 \n", + "07b_treatment_util(trauma) 0.501010 0.506908 0.506908 \n", + "08_total_time(trauma) 292.280487 282.783517 282.783517 \n", + "09_throughput 195.300000 138.540000 156.600000 " + ] + }, + "execution_count": 10, + "metadata": {}, + "output_type": "execute_result" + } + ], "source": [ "scenario_summary_frame(results)" ] @@ -83,17 +327,18 @@ }, { "cell_type": "code", - "execution_count": null, + "execution_count": 15, "id": "cda7e07f-7ab7-45a9-9559-62073c3a4fc0", "metadata": {}, "outputs": [], "source": [ - "from treat_sim.model import Scenario, multiple_replications" + "from treat_sim.model import Scenario, multiple_replications\n", + "from treat_sim.datasets import load_alternative_arrivals" ] }, { "cell_type": "code", - "execution_count": null, + "execution_count": 12, "id": "cdb93109-1054-4493-b4c6-743d70be8a8c", "metadata": {}, "outputs": [], @@ -103,13 +348,559 @@ }, { "cell_type": "code", - "execution_count": null, + "execution_count": 13, "id": "f1b2bdd9-4015-49db-9ca7-92c32c70935e", "metadata": {}, - "outputs": [], + "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", + " \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", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + "
countmeanstdmin25%50%75%max
00_arrivals5.0231.805.07227.00229.00230.00233.00240.00
01a_triage_wait5.030.5115.3917.6824.2824.8028.6657.12
01b_triage_util5.00.610.030.570.610.610.620.64
02a_registration_wait5.0106.1311.7290.00103.24103.61112.24121.54
02b_registration_util5.00.860.010.840.850.850.870.87
03a_examination_wait5.024.847.0714.6921.3725.5131.0931.56
03b_examination_util5.00.860.010.850.860.860.870.88
04a_treatment_wait(non_trauma)5.0134.8727.0894.02120.25150.35152.48157.26
04b_treatment_util(non_trauma)5.00.900.020.870.890.890.910.93
05_total_time(non-trauma)5.0232.0014.87208.36233.30233.88234.76249.69
06a_trauma_wait5.0209.3763.74133.81149.37236.51250.71276.42
06b_trauma_util5.00.900.080.830.860.870.881.03
07a_treatment_wait(trauma)5.011.084.893.729.2312.2513.7016.52
07b_treatment_util(trauma)5.00.500.130.370.420.460.610.66
08_total_time(trauma)5.0350.1363.75281.06301.52346.70380.83440.52
09_throughput5.0167.004.64161.00164.00167.00171.00172.00
\n", + "
" + ], + "text/plain": [ + " count mean std min 25% 50% \\\n", + "00_arrivals 5.0 231.80 5.07 227.00 229.00 230.00 \n", + "01a_triage_wait 5.0 30.51 15.39 17.68 24.28 24.80 \n", + "01b_triage_util 5.0 0.61 0.03 0.57 0.61 0.61 \n", + "02a_registration_wait 5.0 106.13 11.72 90.00 103.24 103.61 \n", + "02b_registration_util 5.0 0.86 0.01 0.84 0.85 0.85 \n", + "03a_examination_wait 5.0 24.84 7.07 14.69 21.37 25.51 \n", + "03b_examination_util 5.0 0.86 0.01 0.85 0.86 0.86 \n", + "04a_treatment_wait(non_trauma) 5.0 134.87 27.08 94.02 120.25 150.35 \n", + "04b_treatment_util(non_trauma) 5.0 0.90 0.02 0.87 0.89 0.89 \n", + "05_total_time(non-trauma) 5.0 232.00 14.87 208.36 233.30 233.88 \n", + "06a_trauma_wait 5.0 209.37 63.74 133.81 149.37 236.51 \n", + "06b_trauma_util 5.0 0.90 0.08 0.83 0.86 0.87 \n", + "07a_treatment_wait(trauma) 5.0 11.08 4.89 3.72 9.23 12.25 \n", + "07b_treatment_util(trauma) 5.0 0.50 0.13 0.37 0.42 0.46 \n", + "08_total_time(trauma) 5.0 350.13 63.75 281.06 301.52 346.70 \n", + "09_throughput 5.0 167.00 4.64 161.00 164.00 167.00 \n", + "\n", + " 75% max \n", + "00_arrivals 233.00 240.00 \n", + "01a_triage_wait 28.66 57.12 \n", + "01b_triage_util 0.62 0.64 \n", + "02a_registration_wait 112.24 121.54 \n", + "02b_registration_util 0.87 0.87 \n", + "03a_examination_wait 31.09 31.56 \n", + "03b_examination_util 0.87 0.88 \n", + "04a_treatment_wait(non_trauma) 152.48 157.26 \n", + "04b_treatment_util(non_trauma) 0.91 0.93 \n", + "05_total_time(non-trauma) 234.76 249.69 \n", + "06a_trauma_wait 250.71 276.42 \n", + "06b_trauma_util 0.88 1.03 \n", + "07a_treatment_wait(trauma) 13.70 16.52 \n", + "07b_treatment_util(trauma) 0.61 0.66 \n", + "08_total_time(trauma) 380.83 440.52 \n", + "09_throughput 171.00 172.00 " + ] + }, + "execution_count": 13, + "metadata": {}, + "output_type": "execute_result" + } + ], "source": [ "multiple_replications(base_case).describe().round(2).T" ] + }, + { + "cell_type": "markdown", + "id": "92b095d3-47a3-4e42-9756-a94d6042ab1f", + "metadata": {}, + "source": [ + "## Test creation of scenario with a different arrival profile" + ] + }, + { + "cell_type": "code", + "execution_count": 16, + "id": "81573b67", + "metadata": {}, + "outputs": [], + "source": [ + "later_arrivals_scenario = Scenario(arrival_profile = load_alternative_arrivals())" + ] + }, + { + "cell_type": "code", + "execution_count": 17, + "id": "bd8b3071-0924-46eb-9e3c-80b417aa1553", + "metadata": {}, + "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", + " \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", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + "
countmeanstdmin25%50%75%max
00_arrivals5.0232.807.66222.00229.00233.00240.00240.00
01a_triage_wait5.031.7921.4218.9619.0322.8728.6769.45
01b_triage_util5.00.610.030.560.610.610.640.64
02a_registration_wait5.0105.519.5591.42100.32110.16110.54115.12
02b_registration_util5.00.840.020.810.830.850.850.86
03a_examination_wait5.024.436.8514.4721.4224.8130.7330.73
03b_examination_util5.00.840.020.820.840.850.860.86
04a_treatment_wait(non_trauma)5.0139.8432.3695.91117.15149.76163.87172.52
04b_treatment_util(non_trauma)5.00.870.030.810.870.880.890.90
05_total_time(non-trauma)5.0231.1215.00206.38227.45238.76241.00242.03
06a_trauma_wait5.0213.1068.45134.03147.32237.89260.12286.14
06b_trauma_util5.00.890.080.830.840.860.881.03
07a_treatment_wait(trauma)5.011.095.182.9910.5110.9614.3416.63
07b_treatment_util(trauma)5.00.490.120.370.420.440.580.66
08_total_time(trauma)5.0341.0064.56261.20301.66340.15373.03428.93
09_throughput5.0162.406.31157.00157.00160.00167.00171.00
\n", + "
" + ], + "text/plain": [ + " count mean std min 25% 50% \\\n", + "00_arrivals 5.0 232.80 7.66 222.00 229.00 233.00 \n", + "01a_triage_wait 5.0 31.79 21.42 18.96 19.03 22.87 \n", + "01b_triage_util 5.0 0.61 0.03 0.56 0.61 0.61 \n", + "02a_registration_wait 5.0 105.51 9.55 91.42 100.32 110.16 \n", + "02b_registration_util 5.0 0.84 0.02 0.81 0.83 0.85 \n", + "03a_examination_wait 5.0 24.43 6.85 14.47 21.42 24.81 \n", + "03b_examination_util 5.0 0.84 0.02 0.82 0.84 0.85 \n", + "04a_treatment_wait(non_trauma) 5.0 139.84 32.36 95.91 117.15 149.76 \n", + "04b_treatment_util(non_trauma) 5.0 0.87 0.03 0.81 0.87 0.88 \n", + "05_total_time(non-trauma) 5.0 231.12 15.00 206.38 227.45 238.76 \n", + "06a_trauma_wait 5.0 213.10 68.45 134.03 147.32 237.89 \n", + "06b_trauma_util 5.0 0.89 0.08 0.83 0.84 0.86 \n", + "07a_treatment_wait(trauma) 5.0 11.09 5.18 2.99 10.51 10.96 \n", + "07b_treatment_util(trauma) 5.0 0.49 0.12 0.37 0.42 0.44 \n", + "08_total_time(trauma) 5.0 341.00 64.56 261.20 301.66 340.15 \n", + "09_throughput 5.0 162.40 6.31 157.00 157.00 160.00 \n", + "\n", + " 75% max \n", + "00_arrivals 240.00 240.00 \n", + "01a_triage_wait 28.67 69.45 \n", + "01b_triage_util 0.64 0.64 \n", + "02a_registration_wait 110.54 115.12 \n", + "02b_registration_util 0.85 0.86 \n", + "03a_examination_wait 30.73 30.73 \n", + "03b_examination_util 0.86 0.86 \n", + "04a_treatment_wait(non_trauma) 163.87 172.52 \n", + "04b_treatment_util(non_trauma) 0.89 0.90 \n", + "05_total_time(non-trauma) 241.00 242.03 \n", + "06a_trauma_wait 260.12 286.14 \n", + "06b_trauma_util 0.88 1.03 \n", + "07a_treatment_wait(trauma) 14.34 16.63 \n", + "07b_treatment_util(trauma) 0.58 0.66 \n", + "08_total_time(trauma) 373.03 428.93 \n", + "09_throughput 167.00 171.00 " + ] + }, + "execution_count": 17, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "multiple_replications(later_arrivals_scenario).describe().round(2).T" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "8b1ff049-9492-45a1-8a17-11459a196cb1", + "metadata": {}, + "outputs": [], + "source": [] } ], "metadata": { @@ -128,7 +919,7 @@ "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", - "version": "3.8.12" + "version": "3.10.14" } }, "nbformat": 4, diff --git a/pyproject.toml b/pyproject.toml index 9bc4caf..b85e614 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -26,7 +26,7 @@ dependencies = [ "numpy>=1.19.2", "pandas>=1.2.3", "scipy>=1.6.1", - "simpy>=4.0.1", + "simpy>=4.0.1" ] [project.urls] @@ -40,4 +40,4 @@ path = "treat_sim/__init__.py" [tool.hatch.build.targets.sdist] include = [ "/treat_sim", -] +] \ No newline at end of file diff --git a/tests/test_datasets.py b/tests/test_datasets.py new file mode 100644 index 0000000..076fa64 --- /dev/null +++ b/tests/test_datasets.py @@ -0,0 +1,107 @@ +"""Datasets testing + +This module provides a set of tests to run against `treat_sim.datasets`. +These tests are either pass or fail and no interpretation is needed. + +There are two example arrivals profiles built into the package (very small size ~380 bytes). +They are loaded from the functions load_nelson_arrivals() and load_alternative_arrivals(). +The package also provides the function valid_arrival_profile() that will return True if all is +okay with a profile format or raise exceptions if there are problems. + +Tests are divided as follows + +1. Dirty tests +Tests that the `treat-sim.datasets` functions fail as expected when certain invalid values. + +2 Functional tests +These provide tests of the correct loading and use of internal datasets plus +#check that a scenario is created correctly when presented with a valid dataset. + +""" +import numpy as np +import pandas as pd + +import pytest + +from treat_sim.datasets import ( + load_nelson_arrivals, + load_alternative_arrivals, + valid_arrival_profile, +) + +from treat_sim.model import Scenario + +### 1. Dirty tests +# tests to check that `treat-sim` fails as +# expected when given certain values. + + +def nelson_arrivals_wrong_cols(): + """Create an arrival profile with incorrect col names + """ + df = load_nelson_arrivals() + df.columns = ["random_col1", "random_col2"] + return df + +def nelson_arrivals_wrong_rows(): + """Create an arrival profile with incorrect num rows. + """ + df1 = load_nelson_arrivals() + df2 = load_nelson_arrivals() + return pd.concat([df1,df2],ignore_index=True) + +@pytest.mark.parametrize( + "arrival_profile, exception_type", + [ + (["period", "arrival_rate"], TypeError), + (["random_str1", "random_str2"], TypeError), + (np.arange(10), TypeError), + (nelson_arrivals_wrong_cols(), ValueError), + (nelson_arrivals_wrong_rows(), ValueError), + + ], +) +def test_invalid_profile(arrival_profile, exception_type): + """tests that exceptions are thrown for various types + of problems with profile + wrong type + wrong columns + wrong number of rows""" + with pytest.raises(exception_type): + valid_arrival_profile(arrival_profile) + + +### 2 Functional tests +# These provide tests of the correct +# loading and use of internal datasets plus +# check that a scenario is created correctly. + +@pytest.mark.parametrize( + "arrival_profile", + [ + (load_nelson_arrivals()), + (load_alternative_arrivals()), + ], +) +def test_valid_profile(arrival_profile): + """ + Test that included datasets pass the + validation test. + """ + assert valid_arrival_profile(arrival_profile) + + +@pytest.mark.parametrize( + "arrival_profile", + [ + (load_nelson_arrivals()), + (load_alternative_arrivals()), + ], +) +def test_scenario_accepts_valid_profile(arrival_profile): + """ + Test that included datasets work with Scenario. + """ + test_scenario = Scenario(arrival_profile = arrival_profile) + + assert isinstance(test_scenario, Scenario) diff --git a/treat_sim/__init__.py b/treat_sim/__init__.py index 098ac8a..e9c887a 100644 --- a/treat_sim/__init__.py +++ b/treat_sim/__init__.py @@ -1,2 +1,2 @@ __author__ = 'Thomas Monks, Alison Harper and Amy Heather' -__version__ = '2.1.0' +__version__ = '2.2.0' diff --git a/treat_sim/data/ed_arrivals_scenario1.csv b/treat_sim/data/ed_arrivals_scenario1.csv new file mode 100644 index 0000000..c61af89 --- /dev/null +++ b/treat_sim/data/ed_arrivals_scenario1.csv @@ -0,0 +1,19 @@ +period,arrival_rate +6AM-7AM,1.36666666666667 +7AM-8AM,1.8 +8AM-9AM,7.83333333333333 +9AM-10AM,9.4333333333333 +10AM-11AM,13.8 +11AM-12PM,25.2666666666667 +12PM-1PM,30.4 +1PM-2PM,17.0666666666667 +2PM-3PM,17.4666666666667 +3PM-4PM,13.0333333333333 +4PM-5PM,12.6 +5PM-6PM,30.8666666666667 +6PM-7PM,19.0333333333333 +7PM-8PM,12.5 +8PM-9PM,6.3 +9PM-10PM,4.06666666666667 +10PM-11PM,2.2 +11PM-12AM,2.1 diff --git a/treat_sim/datasets.py b/treat_sim/datasets.py new file mode 100644 index 0000000..9881c60 --- /dev/null +++ b/treat_sim/datasets.py @@ -0,0 +1,76 @@ +""" +Datasets + +treat_sim bundles small datasets that can be used by the model for experimentation. +Datasets are returned as pandas dataframes. +""" + +from pathlib import Path + +import pandas as pd + +# path to Nelson arrival profile CSV file +DEFAULT_NSPP_PROFILE_FILE = "ed_arrivals.csv" +SCENARIO_NSPP_PROFILE_FILE = "ed_arrivals_scenario1.csv" + + +def load_nelson_arrivals() -> pd.DataFrame: + """Default time dependent arrival profile from Nelson 2013 + + Arrival rates between between 6am and 12am broken down into + 60 minute intervals. Duration of day is 1080 minutes (18 hours * 60 mins). + + Returns a pd.DataFrame with 2 columns: period and arrival_rate + + Returns: + ------- + pd.DataFrame + """ + path_to_file = Path(__file__).parent.joinpath("data", DEFAULT_NSPP_PROFILE_FILE) + return pd.read_csv(path_to_file) + + +def load_alternative_arrivals() -> pd.DataFrame: + """An example alternative arrival profile. + + In this scenario the Treatment Centre has number of arrivals overall + as the default scenario, but their is a shift in numbers towards later in the day + with a higher peak at 6pm + + Arrival rates between between 6am and 12am broken down into + 60 minute intervals. Duration of day is 1080 minutes (18 hours * 60 mins). + + Returns a pd.DataFrame with 2 columns: period and arrival_rate + + Returns: + ------- + pd.DataFrame + """ + path_to_file = Path(__file__).parent.joinpath("data", SCENARIO_NSPP_PROFILE_FILE) + return pd.read_csv(path_to_file) + + +def valid_arrival_profile(arrival_profile: pd.DataFrame) -> bool: + """ + Provides a simple check that a dataframe containing an arrival + profile is in a valid format. + + Raise an exception if invalid + """ + + if not isinstance(arrival_profile, pd.DataFrame): + raise TypeError( + "Invalid arrival profile. arrival_profile must a DataFrame in the correct format." + ) + + if not {"period", "arrival_rate"}.issubset(arrival_profile.columns): + raise ValueError( + "Invalid arrival profile. DataFrame must contain period and arrival_rate columns" + ) + + if arrival_profile.shape[0] != 18: + raise ValueError( + f"Invalid arrival profile. Profile should contain 18 period. But selected DataFrame contains {arrival_profile.shape[0]} rows." + ) + + return True diff --git a/treat_sim/distributions.py b/treat_sim/distributions.py index 4fe0a96..e504bfd 100644 --- a/treat_sim/distributions.py +++ b/treat_sim/distributions.py @@ -11,102 +11,140 @@ * Bernoulli * Normal * Uniform + +Code taken from our MIT licensed sim-tools package (`pip install sim-tools`) +https://github.com/TomMonks/sim-tools ''' +from abc import ABC, abstractmethod import numpy as np import math -class Exponential: - ''' +from typing import Optional + + +class Distribution(ABC): + """ + Distribution abstract class + All distributions derived from it. + """ + + def __init__(self, random_seed: Optional[int] = None): + self.rng = np.random.default_rng(random_seed) + + @abstractmethod + def sample(self, size: Optional[int] = None) -> float | np.ndarray: + """ + Generate a sample from the distribution + + Params: + ------- + size: int, optional (default=None) + the number of samples to return. If size=None then a single + sample is returned. + + Returns: + ------- + np.ndarray or scalar + """ + pass + + +class Exponential(Distribution): + """ Convenience class for the exponential distribution. packages up distribution parameters, seed and random generator. - ''' - def __init__(self, mean, random_seed=None): - ''' + """ + + def __init__(self, mean: float, random_seed: Optional[int] = None): + """ Constructor - + Params: ------ mean: float The mean of the exponential distribution - + random_seed: int, optional (default=None) A random seed to reproduce samples. If set to none then a unique sample is created. - ''' - self.rng = np.random.default_rng(seed=random_seed) + """ + super().__init__(random_seed) self.mean = mean - - def sample(self, size=None): - ''' + + def sample(self, size: Optional[int] = None) -> float | np.ndarray: + """ 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. - ''' + """ return self.rng.exponential(self.mean, size=size) - -class Bernoulli: - ''' + +class Bernoulli(Distribution): + """ Convenience class for the Bernoulli distribution. packages up distribution parameters, seed and random generator. - ''' - def __init__(self, p, random_seed=None): - ''' + """ + + def __init__(self, p: float, random_seed: Optional[int] = None): + """ Constructor - + Params: ------ p: float probability of drawing a 1 - + random_seed: int, optional (default=None) A random seed to reproduce samples. If set to none then a unique sample is created. - ''' - self.rng = np.random.default_rng(seed=random_seed) + """ + super().__init__(random_seed) self.p = p - - def sample(self, size=None): - ''' + + def sample(self, size: Optional[int] = None) -> float | np.ndarray: + """ 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. - ''' + """ return self.rng.binomial(n=1, p=self.p, size=size) -class Lognormal: + +class Lognormal(Distribution): """ Encapsulates a lognormal distirbution """ - def __init__(self, mean, stdev, random_seed=None): + + def __init__(self, mean: float, stdev: float, random_seed: Optional[int] = None): """ Params: ------- mean: float mean of the lognormal distribution - + stdev: float standard dev of the lognormal distribution - + random_seed: int, optional (default=None) Random seed to control sampling """ - self.rng = np.random.default_rng(seed=random_seed) + super().__init__(random_seed) mu, sigma = self.normal_moments_from_lognormal(mean, stdev**2) self.mu = mu self.sigma = sigma - + def normal_moments_from_lognormal(self, m, v): - ''' + """ Returns mu and sigma of normal distribution underlying a lognormal with mean m and variance v source: https://blogs.sas.com/content/iml/2014/06/04/simulate-lognormal @@ -118,31 +156,37 @@ def normal_moments_from_lognormal(self, m, v): mean of lognormal distribution v: float variance of lognormal distribution - + Returns: ------- (float, float) - ''' + """ phi = math.sqrt(v + m**2) - mu = math.log(m**2/phi) - sigma = math.sqrt(math.log(phi**2/m**2)) + mu = math.log(m**2 / phi) + sigma = math.sqrt(math.log(phi**2 / m**2)) return mu, sigma - - def sample(self): + + def sample(self, size: Optional[int] = None) -> float | np.ndarray: """ Sample from the normal distribution """ - return self.rng.lognormal(self.mu, self.sigma) + return self.rng.lognormal(self.mu, self.sigma, size=size) -class Normal: +class Normal(Distribution): ''' Convenience class for the normal distribution. packages up distribution parameters, seed and random generator. Use the minimum parameter to truncate the distribution ''' - def __init__(self, mean, sigma, minimum=None, random_seed=None): + def __init__( + self, + mean: float, + sigma: float, + minimum: Optional[float] = None, + random_seed: Optional[int] = None, + ): ''' Constructor @@ -153,6 +197,10 @@ def __init__(self, mean, sigma, minimum=None, random_seed=None): sigma: float The stdev of the normal distribution + + minimum: float + Truncate the normal distribution to a minimum + value. random_seed: int, optional (default=None) A random seed to reproduce samples. If set to none then a unique @@ -163,7 +211,7 @@ def __init__(self, mean, sigma, minimum=None, random_seed=None): self.sigma = sigma self.minimum = minimum - def sample(self, size=None): + def sample(self, size: Optional[int] = None) -> float | np.ndarray: ''' Generate a sample from the normal distribution @@ -185,40 +233,45 @@ def sample(self, size=None): samples[neg_idx] = self.minimum return samples - -class Uniform(): - ''' + + + +class Uniform(Distribution): + """ Convenience class for the Uniform distribution. packages up distribution parameters, seed and random generator. - ''' - def __init__(self, low, high, random_seed=None): - ''' + """ + + def __init__( + self, low: float, high: float, random_seed: Optional[int] = None + ) -> float | np.ndarray: + """ Constructor - + Params: ------ low: float lower range of the uniform - + high: float upper range of the uniform - + random_seed: int, 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) + """ + super().__init__(random_seed) self.low = low self.high = high - - def sample(self, size=None): - ''' + + def sample(self, size: Optional[int] = None) -> float | np.ndarray: + """ Generate a sample from the uniform distribution - + Params: ------- size: int, optional (default=None) the number of samples to return. If size=None then a single sample is returned. - ''' - return self.rand.uniform(low=self.low, high=self.high, size=size) + """ + return self.rng.uniform(low=self.low, high=self.high, size=size) diff --git a/treat_sim/model.py b/treat_sim/model.py index 9e5dc6b..c60663d 100644 --- a/treat_sim/model.py +++ b/treat_sim/model.py @@ -1,7 +1,7 @@ #!/usr/bin/env python # coding: utf-8 -''' +""" FirstTreatment: A health clinic based in the US. This example is based on exercise 13 from Nelson (2013) page 170. @@ -25,15 +25,18 @@ in a cubicle before being discharged. In this model treatment of trauma and non-trauma patients is modelled seperately -''' +""" import numpy as np import pandas as pd import itertools + import simpy -from treat_sim.distributions import ( - Exponential, Normal, Uniform, Bernoulli, Lognormal) +from typing import Optional, Union, List, Dict + +from treat_sim.distributions import Exponential, Normal, Uniform, Bernoulli, Lognormal +from treat_sim.datasets import load_nelson_arrivals, valid_arrival_profile # Constants and defaults for modelling **as-is** @@ -44,7 +47,7 @@ # registration parameters DEFAULT_REG_MEAN = 5.0 -DEFAULT_REG_VAR= 2.0 +DEFAULT_REG_VAR = 2.0 # examination parameters DEFAULT_EXAM_MEAN = 16.0 @@ -71,10 +74,8 @@ # Time dependent arrival rates data # The data for arrival rates varies between clinic opening at 6am and closure at -# 12am. - -NSPP_PATH = 'https://raw.githubusercontent.com/pythonhealthdatascience/' \ - + 'stars-treat-sim/main/treat_sim/data/ed_arrivals.csv' +# 12am +DEFAULT_NSPP_PROFILE_FILE = "auto" # Resource counts @@ -106,161 +107,183 @@ TRACE = False # list of metrics useful for external apps -RESULT_FIELDS = ['00_arrivals', - '01a_triage_wait', - '01b_triage_util', - '02a_registration_wait', - '02b_registration_util', - '03a_examination_wait', - '03b_examination_util', - '04a_treatment_wait(non_trauma)', - '04b_treatment_util(non_trauma)', - '05_total_time(non-trauma)', - '06a_trauma_wait', - '06b_trauma_util', - '07a_treatment_wait(trauma)', - '07b_treatment_util(trauma)', - '08_total_time(trauma)', - '09_throughput'] +RESULT_FIELDS = [ + "00_arrivals", + "01a_triage_wait", + "01b_triage_util", + "02a_registration_wait", + "02b_registration_util", + "03a_examination_wait", + "03b_examination_util", + "04a_treatment_wait(non_trauma)", + "04b_treatment_util(non_trauma)", + "05_total_time(non-trauma)", + "06a_trauma_wait", + "06b_trauma_util", + "07a_treatment_wait(trauma)", + "07b_treatment_util(trauma)", + "08_total_time(trauma)", + "09_throughput", +] # list of metrics useful for external apps -RESULT_LABELS = {'00_arrivals': 'Arrivals', - '01a_triage_wait': 'Triage Wait (mins)', - '01b_triage_util': 'Triage Utilisation', - '02a_registration_wait': 'Registration Waiting Time (mins)', - '02b_registration_util': 'Registration Utilisation', - '03a_examination_wait': 'Examination Waiting Time (mins)', - '03b_examination_util': 'Examination Utilisation', - '04a_treatment_wait(non_trauma)': 'Non-trauma cubicle waiting time (mins)', - '04b_treatment_util(non_trauma)': 'Non-trauma cubicle utilisation', - '05_total_time(non-trauma)': 'Total time (non-trauma)', - '06a_trauma_wait': 'Trauma stabilisation waiting time (mins)', - '06b_trauma_util': 'Trauma stabilisation utilisation', - '07a_treatment_wait(trauma)': 'Trauma cubicle waiting time (mins)', - '07b_treatment_util(trauma)': 'Trauma cubicle utilisation', - '08_total_time(trauma)': 'Total time (trauma)', - '09_throughput': 'throughput'} +RESULT_LABELS = { + "00_arrivals": "Arrivals", + "01a_triage_wait": "Triage Wait (mins)", + "01b_triage_util": "Triage Utilisation", + "02a_registration_wait": "Registration Waiting Time (mins)", + "02b_registration_util": "Registration Utilisation", + "03a_examination_wait": "Examination Waiting Time (mins)", + "03b_examination_util": "Examination Utilisation", + "04a_treatment_wait(non_trauma)": "Non-trauma cubicle waiting time (mins)", + "04b_treatment_util(non_trauma)": "Non-trauma cubicle utilisation", + "05_total_time(non-trauma)": "Total time (non-trauma)", + "06a_trauma_wait": "Trauma stabilisation waiting time (mins)", + "06b_trauma_util": "Trauma stabilisation utilisation", + "07a_treatment_wait(trauma)": "Trauma cubicle waiting time (mins)", + "07b_treatment_util(trauma)": "Trauma cubicle utilisation", + "08_total_time(trauma)": "Total time (trauma)", + "09_throughput": "throughput", +} # Utility functions -def trace(msg): - ''' + +def trace(msg: str) -> None: + """ Utility function for printing a trace as the simulation model executes. Set the TRACE constant to False, to turn tracing off. - + Params: ------- msg: str string to print to screen. - ''' + """ if TRACE: print(msg) # ## Model parameterisation + class Scenario: - ''' + """ Container class for scenario parameters/arguments - + Passed to a model and its process classes - ''' - def __init__(self, random_number_set=DEFAULT_RNG_SET, - n_triage=DEFAULT_N_TRIAGE, - n_reg=DEFAULT_N_REG, - n_exam=DEFAULT_N_EXAM, - n_trauma=DEFAULT_N_TRAUMA, - n_cubicles_1=DEFAULT_N_CUBICLES_1, - n_cubicles_2=DEFAULT_N_CUBICLES_2, - triage_mean=DEFAULT_TRIAGE_MEAN, - reg_mean=DEFAULT_REG_MEAN, - reg_var=DEFAULT_REG_VAR, - exam_mean=DEFAULT_EXAM_MEAN, - exam_var=DEFAULT_EXAM_VAR, - exam_min=DEFAULT_EXAM_MIN, - trauma_mean=DEFAULT_TRAUMA_MEAN, - trauma_treat_mean=DEFAULT_TRAUMA_TREAT_MEAN, - trauma_treat_var=DEFAULT_TRAUMA_TREAT_VAR, - non_trauma_treat_mean=DEFAULT_NON_TRAUMA_TREAT_MEAN, - non_trauma_treat_var=DEFAULT_NON_TRAUMA_TREAT_VAR, - non_trauma_treat_p=DEFAULT_NON_TRAUMA_TREAT_P, - prob_trauma=DEFAULT_PROB_TRAUMA): - ''' + """ + + def __init__( + self, + random_number_set: Optional[int | None] = DEFAULT_RNG_SET, + n_triage: Optional[int] = DEFAULT_N_TRIAGE, + n_reg: Optional[int] = DEFAULT_N_REG, + n_exam: Optional[int] = DEFAULT_N_EXAM, + n_trauma: Optional[int] = DEFAULT_N_TRAUMA, + n_cubicles_1: Optional[int] = DEFAULT_N_CUBICLES_1, + n_cubicles_2: Optional[int] = DEFAULT_N_CUBICLES_2, + triage_mean: Optional[int] = DEFAULT_TRIAGE_MEAN, + reg_mean: Optional[float] = DEFAULT_REG_MEAN, + reg_var: Optional[float] = DEFAULT_REG_VAR, + exam_mean: Optional[float] = DEFAULT_EXAM_MEAN, + exam_var: Optional[float] = DEFAULT_EXAM_VAR, + exam_min: Optional[float] = DEFAULT_EXAM_MIN, + trauma_mean: Optional[float] = DEFAULT_TRAUMA_MEAN, + trauma_treat_mean: Optional[float] = DEFAULT_TRAUMA_TREAT_MEAN, + trauma_treat_var: Optional[float] = DEFAULT_TRAUMA_TREAT_VAR, + non_trauma_treat_mean: Optional[float] = DEFAULT_NON_TRAUMA_TREAT_MEAN, + non_trauma_treat_var: Optional[float] = DEFAULT_NON_TRAUMA_TREAT_VAR, + non_trauma_treat_p: Optional[float] = DEFAULT_NON_TRAUMA_TREAT_P, + prob_trauma: Optional[float] = DEFAULT_PROB_TRAUMA, + arrival_profile: Optional[pd.DataFrame | None] = None, + ): + """ Create a scenario to parameterise the simulation model - + Parameters: ----------- random_number_set: int, optional (default=DEFAULT_RNG_SET) Set to control the initial seeds of each stream of pseudo random numbers used in the model. - + n_triage: int The number of triage cubicles - + n_reg: int The number of registration clerks - + n_exam: int The number of examination rooms - + n_trauma: int The number of trauma bays for stabilisation - + n_cubicles_1: int The number of non-trauma treatment cubicles - + n_cubicles_2: int The number of trauma treatment cubicles - + triage_mean: float Mean duration of the triage distribution (Exponential) - + reg_mean: float Mean duration of the registration distribution (Lognormal) - + reg_var: float Variance of the registration distribution (Lognormal) - + exam_mean: float Mean of the examination distribution (Normal) - + exam_var: float Variance of the examination distribution (Normal) exam_min: float The minimum value that an examination can take (Truncated Normal) - + trauma_mean: float Mean of the trauma stabilisation distribution (Exponential) - + trauma_treat_mean: float Mean of the trauma cubicle treatment distribution (Lognormal) - + trauma_treat_var: float Variance of the trauma cubicle treatment distribution (Lognormal) - + non_trauma_treat_mean: float Mean of the non trauma treatment distribution - + non_trauma_treat_var: float Variance of the non trauma treatment distribution - + non_trauma_treat_p: float Probability non trauma patient requires treatment - + prob_trauma: float probability that a new arrival is a trauma patient. - ''' + + arrival_profile: None | pd.DataFrame + pandas dataframe containing the arrival profile for the day. + This should have two columns 'period' and 'arrival_rate'. There + should be 18 rows (each representing an hour) + + See also: + --------- + datasets.load_nelson_arrivals() + datasets.load_alternative_arrivals() + datasets.valid_arrival_profile() + """ # sampling self.random_number_set = random_number_set - + # store parameters for sampling self.triage_mean = triage_mean self.reg_mean = reg_mean self.reg_var = reg_var - self.exam_mean= exam_mean + self.exam_mean = exam_mean self.exam_var = exam_var self.exam_min = exam_min self.trauma_mean = trauma_mean @@ -270,272 +293,291 @@ def __init__(self, random_number_set=DEFAULT_RNG_SET, self.non_trauma_treat_var = non_trauma_treat_var self.non_trauma_treat_p = non_trauma_treat_p self.prob_trauma = prob_trauma - + + if arrival_profile is None: + self.arrivals = load_nelson_arrivals() + else: + if valid_arrival_profile(arrival_profile): + self.arrivals = arrival_profile + self.init_sampling() - + # count of each type of resource - self.init_resourse_counts(n_triage, n_reg, n_exam, n_trauma, - n_cubicles_1, n_cubicles_2) - - def set_random_no_set(self, random_number_set): - ''' - Controls the random sampling + self.init_resourse_counts( + n_triage, n_reg, n_exam, n_trauma, n_cubicles_1, n_cubicles_2 + ) + + def set_random_no_set(self, random_number_set: int) -> None: + """ + Controls the random sampling Parameters: ---------- random_number_set: int Used to control the set of pseudo random numbers used by the distributions in the simulation. - ''' + """ self.random_number_set = random_number_set self.init_sampling() - def init_resourse_counts(self, n_triage, n_reg, n_exam, n_trauma, - n_cubicles_1, n_cubicles_2): - ''' + def init_resourse_counts( + self, + n_triage: int, + n_reg: int, + n_exam: int, + n_trauma: int, + n_cubicles_1: int, + n_cubicles_2: int, + ): + """ Init the counts of resources to default values... - ''' + """ self.n_triage = n_triage self.n_reg = n_reg self.n_exam = n_exam self.n_trauma = n_trauma - + # non-trauma (1), trauma (2) treatment cubicles self.n_cubicles_1 = n_cubicles_1 self.n_cubicles_2 = n_cubicles_2 def init_sampling(self): - ''' - Create the distributions used by the model and initialise + """ + Create the distributions used by the model and initialise the random seeds of each. - ''' + """ # MODIFICATION. Better method for producing n non-overlapping streams seed_sequence = np.random.SeedSequence(self.random_number_set) - + # Generate n high quality child seeds self.seeds = seed_sequence.spawn(N_STREAMS) - + # create distributions - + # Triage duration - self.triage_dist = Exponential(self.triage_mean, - random_seed=self.seeds[0]) - + self.triage_dist = Exponential(self.triage_mean, random_seed=self.seeds[0]) + # Registration duration (non-trauma only) - self.reg_dist = Lognormal(self.reg_mean, - np.sqrt(self.reg_var), - random_seed=self.seeds[1]) - + self.reg_dist = Lognormal( + self.reg_mean, np.sqrt(self.reg_var), random_seed=self.seeds[1] + ) + # Evaluation (non-trauma only) - self.exam_dist = Normal(self.exam_mean, - np.sqrt(self.exam_var), - minimum=self.exam_min, - random_seed=self.seeds[2]) - + self.exam_dist = Normal( + self.exam_mean, + np.sqrt(self.exam_var), + minimum=self.exam_min, + random_seed=self.seeds[2], + ) + # Trauma/stablisation duration (trauma only) - self.trauma_dist = Exponential(self.trauma_mean, - random_seed=self.seeds[3]) - + self.trauma_dist = Exponential(self.trauma_mean, random_seed=self.seeds[3]) + # Non-trauma treatment - self.nt_treat_dist = Lognormal(self.non_trauma_treat_mean, - np.sqrt(self.non_trauma_treat_var), - random_seed=self.seeds[4]) - + self.nt_treat_dist = Lognormal( + self.non_trauma_treat_mean, + np.sqrt(self.non_trauma_treat_var), + random_seed=self.seeds[4], + ) + # treatment of trauma patients - self.treat_dist = Lognormal(self.trauma_treat_mean, - np.sqrt(self.trauma_treat_var), - random_seed=self.seeds[5]) - + self.treat_dist = Lognormal( + self.trauma_treat_mean, + np.sqrt(self.trauma_treat_var), + random_seed=self.seeds[5], + ) + # probability of non-trauma patient requiring treatment - self.nt_p_treat_dist = Bernoulli(self.non_trauma_treat_p, - random_seed=self.seeds[6]) - - + self.nt_p_treat_dist = Bernoulli( + self.non_trauma_treat_p, random_seed=self.seeds[6] + ) + # probability of non-trauma versus trauma patient - self.p_trauma_dist = Bernoulli(self.prob_trauma, - random_seed=self.seeds[7]) - + self.p_trauma_dist = Bernoulli(self.prob_trauma, random_seed=self.seeds[7]) + # init sampling for non-stationary poisson process self.init_nspp() - + def init_nspp(self): - + # read arrival profile - self.arrivals = pd.read_csv(NSPP_PATH) - self.arrivals['mean_iat'] = 60 / self.arrivals['arrival_rate'] - + self.arrivals["mean_iat"] = 60.0 / self.arrivals["arrival_rate"] + # maximum arrival rate (smallest time between arrivals) - self.lambda_max = self.arrivals['arrival_rate'].max() - + self.lambda_max = self.arrivals["arrival_rate"].max() + # thinning exponential - self.arrival_dist = Exponential(60.0 / self.lambda_max, - random_seed=self.seeds[8]) - + self.arrival_dist = Exponential( + 60.0 / self.lambda_max, random_seed=self.seeds[8] + ) + # thinning uniform rng - self.thinning_rng = Uniform(low=0.0, high=1.0, - random_seed=self.seeds[9]) + self.thinning_rng = Uniform(low=0.0, high=1.0, random_seed=self.seeds[9]) # ## Patient Pathways Process Logic class TraumaPathway: - ''' + """ Encapsulates the process a patient with severe injuries or illness. - + These patients are signed into the ED and triaged as having severe injuries or illness. - - Patients are stabilised in resus (trauma) and then sent to Treatment. + + Patients are stabilised in resus (trauma) and then sent to Treatment. Following treatment they are discharged. - ''' - def __init__(self, identifier, env, args): - ''' + """ + + def __init__(self, identifier: int, env: simpy.Environment, args: Scenario) -> None: + """ Constructor method - + Params: ----- identifier: int a numeric identifier for the patient. - + env: simpy.Environment the simulation environment - + args: Scenario Container class for the simulation parameters - - ''' + + """ self.identifier = identifier self.env = env self.args = args - + # metrics self.arrival = -np.inf self.wait_triage = -np.inf self.wait_trauma = -np.inf self.wait_treat = -np.inf self.total_time = -np.inf - + self.triage_duration = -np.inf self.trauma_duration = -np.inf self.treat_duration = -np.inf - + def execute(self): - ''' + """ simulates the major treatment process for a patient - + 1. request and wait for sign-in/triage 2. trauma 3. treatment - ''' + """ # record the time of arrival and entered the triage queue self.arrival = self.env.now - # request sign-in/triage + # request sign-in/triage with self.args.triage.request() as req: yield req # record the waiting time for triage - self.wait_triage = self.env.now - self.arrival - - trace(f'patient {self.identifier} triaged to trauma ' - f'{self.env.now:.3f}') - + self.wait_triage = self.env.now - self.arrival + + trace(f"patient {self.identifier} triaged to trauma " f"{self.env.now:.3f}") + # sample triage duration. self.triage_duration = self.args.triage_dist.sample() yield self.env.timeout(self.triage_duration) self.triage_complete() - + # record the time that entered the trauma queue start_wait = self.env.now - - # request trauma room + + # request trauma room with self.args.trauma.request() as req: yield req - + # record the waiting time for trauma self.wait_trauma = self.env.now - start_wait - + # sample stablisation duration. self.trauma_duration = self.args.trauma_dist.sample() yield self.env.timeout(self.trauma_duration) - + self.trauma_complete() - + # record the time that entered the treatment queue start_wait = self.env.now - - # request treatment cubicle + + # request treatment cubicle with self.args.cubicle_2.request() as req: yield req - + # record the waiting time for trauma self.wait_treat = self.env.now - start_wait - trace(f'treatment of patient {self.identifier} at ' - f'{self.env.now:.3f}') - + trace(f"treatment of patient {self.identifier} at " f"{self.env.now:.3f}") + # sample treatment duration. self.treat_duration = self.args.treat_dist.sample() yield self.env.timeout(self.treat_duration) - + self.treatment_complete() - + # total time in system - self.total_time = self.env.now - self.arrival - - def triage_complete(self): - ''' + self.total_time = self.env.now - self.arrival + + def triage_complete(self) -> None: + """ Triage complete event - ''' - trace(f'triage {self.identifier} complete {self.env.now:.3f}; ' - f'waiting time was {self.wait_triage:.3f}') - - def trauma_complete(self): - ''' + """ + trace( + f"triage {self.identifier} complete {self.env.now:.3f}; " + f"waiting time was {self.wait_triage:.3f}" + ) + + def trauma_complete(self) -> None: + """ Patient stay in trauma is complete. - ''' - trace(f'stabilisation of patient {self.identifier} at ' - f'{self.env.now:.3f}') - - def treatment_complete(self): - ''' + """ + trace(f"stabilisation of patient {self.identifier} at " f"{self.env.now:.3f}") + + def treatment_complete(self) -> None: + """ Treatment complete event - ''' - trace(f'patient {self.identifier} treatment complete {self.env.now:.3f}; ' - f'waiting time was {self.wait_treat:.3f}') + """ + trace( + f"patient {self.identifier} treatment complete {self.env.now:.3f}; " + f"waiting time was {self.wait_treat:.3f}" + ) -class NonTraumaPathway(object): - ''' +class NonTraumaPathway: + """ Encapsulates the process a patient with minor injuries and illness. - - These patients are signed into the ED and triaged as having minor - complaints and streamed to registration and then examination. - - Post examination 40% are discharged while 60% proceed to treatment. + + These patients are signed into the ED and triaged as having minor + complaints and streamed to registration and then examination. + + Post examination 40% are discharged while 60% proceed to treatment. Following treatment they are discharged. - ''' - def __init__(self, identifier, env, args): - ''' + """ + + def __init__(self, identifier: int, env: simpy.Environment, args: Scenario) -> None: + """ Constructor method - + Params: ----- identifier: int a numeric identifier for the patient. - + env: SimPy.Environment the simulation environment - + args: Scenario Container class for the simulation parameters - - ''' + + """ self.identifier = identifier self.env = env self.args = args - + # triage resource self.triage = args.triage - + # metrics self.arrival = -np.inf self.wait_triage = -np.inf @@ -543,142 +585,154 @@ def __init__(self, identifier, env, args): self.wait_exam = -np.inf self.wait_treat = -np.inf self.total_time = -np.inf - + self.triage_duration = -np.inf self.reg_duration = -np.inf self.exam_duration = -np.inf self.treat_duration = -np.inf - - + def execute(self): - ''' + """ simulates the non-trauma/minor treatment process for a patient - + 1. request and wait for sign-in/triage 2. patient registration 3. examination 4.1 40% discharged 4.2 60% treatment then discharge - ''' + """ # record the time of arrival and entered the triage queue self.arrival = self.env.now - # request sign-in/triage + # request sign-in/triage with self.triage.request() as req: yield req - + # record the waiting time for triage self.wait_triage = self.env.now - self.arrival - trace(f'patient {self.identifier} triaged to minors ' - f'{self.env.now:.3f}') - + trace(f"patient {self.identifier} triaged to minors " f"{self.env.now:.3f}") + # sample triage duration. self.triage_duration = self.args.triage_dist.sample() yield self.env.timeout(self.triage_duration) - - trace(f'triage {self.identifier} complete {self.env.now:.3f}; ' - f'waiting time was {self.wait_triage:.3f}') - + + trace( + f"triage {self.identifier} complete {self.env.now:.3f}; " + f"waiting time was {self.wait_triage:.3f}" + ) + # record the time that entered the registration queue start_wait = self.env.now - - # request registration clert + + # request registration clert with self.args.registration.request() as req: yield req - + # record the waiting time for registration self.wait_reg = self.env.now - start_wait - trace(f'registration of patient {self.identifier} at ' - f'{self.env.now:.3f}') - + trace( + f"registration of patient {self.identifier} at " f"{self.env.now:.3f}" + ) + # sample registration duration. self.reg_duration = self.args.reg_dist.sample() yield self.env.timeout(self.reg_duration) - - trace(f'patient {self.identifier} registered at' - f'{self.env.now:.3f}; ' - f'waiting time was {self.wait_reg:.3f}') - + + trace( + f"patient {self.identifier} registered at" + f"{self.env.now:.3f}; " + f"waiting time was {self.wait_reg:.3f}" + ) + # record the time that entered the evaluation queue start_wait = self.env.now - + # request examination resource with self.args.exam.request() as req: yield req - + # record the waiting time for registration self.wait_exam = self.env.now - start_wait - trace(f'examination of patient {self.identifier} begins ' - f'{self.env.now:.3f}') - + trace( + f"examination of patient {self.identifier} begins " + f"{self.env.now:.3f}" + ) + # sample examination duration. self.exam_duration = self.args.exam_dist.sample() yield self.env.timeout(self.exam_duration) - - trace(f'patient {self.identifier} examination complete ' - f'at {self.env.now:.3f};' - f'waiting time was {self.wait_exam:.3f}') - + + trace( + f"patient {self.identifier} examination complete " + f"at {self.env.now:.3f};" + f"waiting time was {self.wait_exam:.3f}" + ) + # sample if patient requires treatment? self.require_treat = self.args.nt_p_treat_dist.sample() - + if self.require_treat: - + # record the time that entered the treatment queue start_wait = self.env.now - + # request treatment cubicle with self.args.cubicle_1.request() as req: yield req # record the waiting time for treatment self.wait_treat = self.env.now - start_wait - trace(f'treatment of patient {self.identifier} begins ' - f'{self.env.now:.3f}') + trace( + f"treatment of patient {self.identifier} begins " + f"{self.env.now:.3f}" + ) # sample treatment duration. self.treat_duration = self.args.nt_treat_dist.sample() yield self.env.timeout(self.treat_duration) - trace(f'patient {self.identifier} treatment complete ' - f'at {self.env.now:.3f};' - f'waiting time was {self.wait_treat:.3f}') - + trace( + f"patient {self.identifier} treatment complete " + f"at {self.env.now:.3f};" + f"waiting time was {self.wait_treat:.3f}" + ) + # total time in system - self.total_time = self.env.now - self.arrival + self.total_time = self.env.now - self.arrival class TreatmentCentreModel: - ''' + """ The treatment centre model - + Patients arrive at random to a treatment centre, are triaged and then processed in either a trauma or non-trauma pathway. - The main class that a user interacts with to run the model is + The main class that a user interacts with to run the model is `TreatmentCentreModel`. This implements a `.run()` method, contains a simple - algorithm for the non-stationary Poisson process for patients arrivals and - inits instances of `TraumaPathway` or `NonTraumaPathway` depending on the - arrival type. + algorithm for the non-stationary Poisson process for patients arrivals and + inits instances of `TraumaPathway` or `NonTraumaPathway` depending on the + arrival type. + + """ - ''' - def __init__(self, args): + def __init__(self, args: Scenario): self.env = simpy.Environment() self.args = args self.init_resources() - + self.patients = [] self.trauma_patients = [] self.non_trauma_patients = [] self.rc_period = None self.results = None - - def init_resources(self): - ''' + + def init_resources(self) -> None: + """ Init the number of resources and store in the arguments container object - + Resource list: 1. Sign-in/triage bays 2. registration clerks @@ -686,83 +740,77 @@ def init_resources(self): 4. trauma bays 5. non-trauma cubicles (1) 6. trauma cubicles (2) - - ''' + + """ # sign/in triage - self.args.triage = simpy.Resource(self.env, - capacity=self.args.n_triage) - + self.args.triage = simpy.Resource(self.env, capacity=self.args.n_triage) + # registration - self.args.registration = simpy.Resource(self.env, - capacity=self.args.n_reg) - + self.args.registration = simpy.Resource(self.env, capacity=self.args.n_reg) + # examination - self.args.exam = simpy.Resource(self.env, - capacity=self.args.n_exam) - + self.args.exam = simpy.Resource(self.env, capacity=self.args.n_exam) + # trauma - self.args.trauma = simpy.Resource(self.env, - capacity=self.args.n_trauma) - + self.args.trauma = simpy.Resource(self.env, capacity=self.args.n_trauma) + # non-trauma treatment - self.args.cubicle_1 = simpy.Resource(self.env, - capacity=self.args.n_cubicles_1) - + self.args.cubicle_1 = simpy.Resource(self.env, capacity=self.args.n_cubicles_1) + # trauma treatment - self.args.cubicle_2 = simpy.Resource(self.env, - capacity=self.args.n_cubicles_2) - - - - def run(self, results_collection_period=DEFAULT_RESULTS_COLLECTION_PERIOD): - ''' - Conduct a single run of the model in its current + self.args.cubicle_2 = simpy.Resource(self.env, capacity=self.args.n_cubicles_2) + + def run( + self, + results_collection_period: Optional[float] = DEFAULT_RESULTS_COLLECTION_PERIOD, + ) -> None: + """ + Conduct a single run of the model in its current configuration - - + + Parameters: ---------- results_collection_period, float, optional default = DEFAULT_RESULTS_COLLECTION_PERIOD - + warm_up, float, optional (default=0) - + length of initial transient period to truncate from results. - + Returns: -------- None - ''' + """ # setup the arrival generator process self.env.process(self.arrivals_generator()) - + # store rc period self.rc_period = results_collection_period - + # run self.env.run(until=results_collection_period) - - - def arrivals_generator(self): - ''' + + def arrivals_generator(self): + """ Simulate the arrival of patients to the model - + Patients either follow a TraumaPathway or NonTraumaPathway simpy process. - - Non stationary arrivals implemented via Thinning acceptance-rejection + + Non stationary arrivals implemented via Thinning acceptance-rejection algorithm. - ''' + """ for patient_count in itertools.count(): # this give us the index of dataframe to use t = int(self.env.now // 60) % self.args.arrivals.shape[0] - lambda_t = self.args.arrivals['arrival_rate'].iloc[t] + lambda_t = self.args.arrivals["arrival_rate"].iloc[t] # set to a large number so that at least 1 sample taken! - u = np.Inf - + u = np.inf + interarrival_time = 0.0 # reject samples if u >= lambda_t / lambda_max @@ -772,9 +820,9 @@ def arrivals_generator(self): # iat yield self.env.timeout(interarrival_time) - - trace(f'patient {patient_count} arrives at: {self.env.now:.3f}') - + + trace(f"patient {patient_count} arrives at: {self.env.now:.3f}") + # sample if the patient is trauma or non-trauma trauma = self.args.p_trauma_dist.sample() if trauma: @@ -783,226 +831,251 @@ def arrivals_generator(self): self.trauma_patients.append(new_patient) else: # create and store a non-trauma patient to update KPIs. - new_patient = NonTraumaPathway(patient_count, self.env, - self.args) + new_patient = NonTraumaPathway(patient_count, self.env, self.args) self.non_trauma_patients.append(new_patient) - + # start the pathway process for the patient self.env.process(new_patient.execute()) + class SimulationSummary: - ''' + """ End of run result processing logic of the simulation model - ''' - def __init__(self, model): - ''' + """ + + def __init__(self, model: TreatmentCentreModel) -> None: + """ Constructor - + Params: ------ model: TraumaCentreModel The model. - ''' + """ self.model = model self.args = model.args self.results = None - - def process_run_results(self): - ''' + + def process_run_results(self) -> None: + """ Calculates statistics at end of run. - ''' + """ self.results = {} - # list of all patients + # list of all patients patients = self.model.non_trauma_patients + self.model.trauma_patients - + # mean triage times (both types of patient) - mean_triage_wait = self.get_mean_metric('wait_triage', patients) - + mean_triage_wait = self.get_mean_metric("wait_triage", patients) + # triage utilisation (both types of patient) - triage_util = self.get_resource_util('triage_duration', - self.args.n_triage, - patients) - + triage_util = self.get_resource_util( + "triage_duration", self.args.n_triage, patients + ) + # mean waiting time for registration (non_trauma) - mean_reg_wait = self.get_mean_metric('wait_reg', - self.model.non_trauma_patients) - + mean_reg_wait = self.get_mean_metric("wait_reg", self.model.non_trauma_patients) + # registration utilisation (trauma) - reg_util = self.get_resource_util('reg_duration', - self.args.n_reg, - self.model.non_trauma_patients) - + reg_util = self.get_resource_util( + "reg_duration", self.args.n_reg, self.model.non_trauma_patients + ) + # mean waiting time for examination (non_trauma) - mean_wait_exam = self.get_mean_metric('wait_exam', - self.model.non_trauma_patients) - + mean_wait_exam = self.get_mean_metric( + "wait_exam", self.model.non_trauma_patients + ) + # examination utilisation (non-trauma) - exam_util = self.get_resource_util('exam_duration', - self.args.n_exam, - self.model.non_trauma_patients) - - - # mean waiting time for treatment (non-trauma) - mean_treat_wait = self.get_mean_metric('wait_treat', - self.model.non_trauma_patients) - + exam_util = self.get_resource_util( + "exam_duration", self.args.n_exam, self.model.non_trauma_patients + ) + + # mean waiting time for treatment (non-trauma) + mean_treat_wait = self.get_mean_metric( + "wait_treat", self.model.non_trauma_patients + ) + # treatment utilisation (non_trauma) - treat_util1 = self.get_resource_util('treat_duration', - self.args.n_cubicles_1, - self.model.non_trauma_patients) - + treat_util1 = self.get_resource_util( + "treat_duration", self.args.n_cubicles_1, self.model.non_trauma_patients + ) + # mean total time (non_trauma) - mean_total = self.get_mean_metric('total_time', - self.model.non_trauma_patients) - - # mean waiting time for trauma - mean_trauma_wait = self.get_mean_metric('wait_trauma', - self.model.trauma_patients) - + mean_total = self.get_mean_metric("total_time", self.model.non_trauma_patients) + + # mean waiting time for trauma + mean_trauma_wait = self.get_mean_metric( + "wait_trauma", self.model.trauma_patients + ) + # trauma utilisation (trauma) - trauma_util = self.get_resource_util('trauma_duration', - self.args.n_trauma, - self.model.trauma_patients) - - # mean waiting time for treatment (trauma) - mean_treat_wait2 = self.get_mean_metric('wait_treat', - self.model.trauma_patients) - + trauma_util = self.get_resource_util( + "trauma_duration", self.args.n_trauma, self.model.trauma_patients + ) + + # mean waiting time for treatment (trauma) + mean_treat_wait2 = self.get_mean_metric( + "wait_treat", self.model.trauma_patients + ) + # treatment utilisation (trauma) - treat_util2 = self.get_resource_util('treat_duration', - self.args.n_cubicles_2, - self.model.trauma_patients) + treat_util2 = self.get_resource_util( + "treat_duration", self.args.n_cubicles_2, self.model.trauma_patients + ) # mean total time (trauma) - mean_total2 = self.get_mean_metric('total_time', - self.model.trauma_patients) - - - self.results = {'00_arrivals':len(patients), - '01a_triage_wait': mean_triage_wait, - '01b_triage_util': triage_util, - '02a_registration_wait':mean_reg_wait, - '02b_registration_util': reg_util, - '03a_examination_wait':mean_wait_exam, - '03b_examination_util': exam_util, - '04a_treatment_wait(non_trauma)':mean_treat_wait, - '04b_treatment_util(non_trauma)':treat_util1, - '05_total_time(non-trauma)':mean_total, - '06a_trauma_wait':mean_trauma_wait, - '06b_trauma_util':trauma_util, - '07a_treatment_wait(trauma)':mean_treat_wait2, - '07b_treatment_util(trauma)':treat_util2, - '08_total_time(trauma)':mean_total2, - '09_throughput': self.get_throughput(patients)} - - - def get_mean_metric(self, metric, patients): - ''' + mean_total2 = self.get_mean_metric("total_time", self.model.trauma_patients) + + self.results = { + "00_arrivals": len(patients), + "01a_triage_wait": mean_triage_wait, + "01b_triage_util": triage_util, + "02a_registration_wait": mean_reg_wait, + "02b_registration_util": reg_util, + "03a_examination_wait": mean_wait_exam, + "03b_examination_util": exam_util, + "04a_treatment_wait(non_trauma)": mean_treat_wait, + "04b_treatment_util(non_trauma)": treat_util1, + "05_total_time(non-trauma)": mean_total, + "06a_trauma_wait": mean_trauma_wait, + "06b_trauma_util": trauma_util, + "07a_treatment_wait(trauma)": mean_treat_wait2, + "07b_treatment_util(trauma)": treat_util2, + "08_total_time(trauma)": mean_total2, + "09_throughput": self.get_throughput(patients), + } + + def get_mean_metric( + self, metric: str, patients: List[Union[TraumaPathway, NonTraumaPathway]] + ) -> float: + """ Calculate mean of the performance measure for the select cohort of patients, - - Only calculates metrics for patients where it has been + + Only calculates metrics for patients where it has been measured. - + Params: ------- metric: str The name of the metric e.g. 'wait_treat' - + patients: list - A list of patients - ''' - mean = np.array([getattr(p, metric) for p in patients - if getattr(p, metric) > -np.inf]).mean() + A list of TraumaPathway and NonTraumaPathway patients + + Returns: + ------- + float + """ + mean = np.array( + [getattr(p, metric) for p in patients if getattr(p, metric) > -np.inf] + ).mean() return mean - - - def get_resource_util(self, metric, n_resources, patients): - ''' + + def get_resource_util( + self, + metric: str, + n_resources: int, + patients: List[Union[TraumaPathway, NonTraumaPathway]], + ) -> float: + """ Calculate proportion of the results collection period where a resource was in use. - + Done by tracking the duration by patient. - - Only calculates metrics for patients where it has been + + Only calculates metrics for patients where it has been measured. - + Params: ------- metric: str The name of the metric e.g. 'treatment_duration' - + + n_resources: int + The number of resources available (e.g. number of trauma rooms) + patients: list - A list of patients - ''' - total = np.array([getattr(p, metric) for p in patients - if getattr(p, metric) > -np.inf]).sum() - + A list of TraumaPathway and NonTraumaPathway patients + + Returns: + ------- + float + """ + total = np.array( + [getattr(p, metric) for p in patients if getattr(p, metric) > -np.inf] + ).sum() + return total / (self.model.rc_period * n_resources) - - def get_throughput(self, patients): - ''' + + def get_throughput( + self, patients: List[Union[TraumaPathway, NonTraumaPathway]] + ) -> int: + """ Returns the total number of patients that have successfully been processed and discharged in the treatment centre (they have a total time record) - + Params: ------- patients: list list of all patient objects simulated. - + Returns: ------ - float - ''' + int + """ return len([p for p in patients if p.total_time > -np.inf]) - + def summary_frame(self): - ''' + """ Returns run results as a pandas.DataFrame - + Returns: ------- pd.DataFrame - ''' - #append to results df + """ + # append to results df if self.results is None: self.process_run_results() - df = pd.DataFrame({'1':self.results}) + df = pd.DataFrame({"1": self.results}) df = df.T - df.index.name = 'rep' + df.index.name = "rep" return df - # ## Executing a model -def single_run(scenario, rc_period=DEFAULT_RESULTS_COLLECTION_PERIOD, - random_no_set=DEFAULT_RNG_SET): - ''' + +def single_run( + scenario: Scenario, + rc_period: Optional[float] = DEFAULT_RESULTS_COLLECTION_PERIOD, + random_no_set: Optional[int] = DEFAULT_RNG_SET, +) -> pd.DataFrame: + """ Perform a single run of the model and return the results - + Parameters: ----------- - + scenario: Scenario object The scenario/parameters to run - + rc_period: int The length of the simulation run that collects results - + random_no_set: int or None, optional (default=DEFAULT_RNG_SET) - Controls the set of random seeds used by the stochastic parts of the + Controls the set of random seeds used by the stochastic parts of the model. Set to different ints to get different results. Set to None for a random set of seeds. - + Returns: -------- pandas.DataFrame: results from single run. - ''' - + """ + # set random number set - this controls sampling for the run. scenario.set_random_no_set(random_no_set) @@ -1011,135 +1084,144 @@ def single_run(scenario, rc_period=DEFAULT_RESULTS_COLLECTION_PERIOD, # run the model model.run(results_collection_period=rc_period) - + # run results summary = SimulationSummary(model) summary_df = summary.summary_frame() - + return summary_df -def multiple_replications(scenario, rc_period=DEFAULT_RESULTS_COLLECTION_PERIOD, - n_reps=5): - ''' +def multiple_replications( + scenario: Scenario, + rc_period: Optional[int] = DEFAULT_RESULTS_COLLECTION_PERIOD, + n_reps: Optional[int] = 5, +) -> pd.DataFrame: + """ Perform multiple replications of the model. - + Params: ------ scenario: Scenario Parameters/arguments to configure the model - + rc_period: float, optional (default=DEFAULT_RESULTS_COLLECTION_PERIOD) - results collection period. + results collection period. the number of minutes to run the model to collect results n_reps: int, optional (default=DEFAULT_N_REPS) Number of independent replications to run. - + Returns: -------- pandas.DataFrame - ''' + """ + + results = [ + single_run(scenario, rc_period, random_no_set=rep) for rep in range(n_reps) + ] - results = [single_run(scenario, rc_period, random_no_set=rep) - for rep in range(n_reps)] - - #format and return results in a dataframe + # format and return results in a dataframe df_results = pd.concat(results) - df_results.index = np.arange(1, len(df_results)+1) - df_results.index.name = 'rep' + df_results.index = np.arange(1, len(df_results) + 1) + df_results.index.name = "rep" return df_results - - # ##Β Scenario Analysis -def get_scenarios(): - ''' + +def get_scenarios() -> Dict[str, Scenario]: + """ Creates a dictionary object containing objects of type `Scenario` to run. - + Returns: -------- dict Contains the scenarios for the model - ''' + """ scenarios = {} - scenarios['base'] = Scenario() - + scenarios["base"] = Scenario() + # extra triage capacity - scenarios['triage+1'] = Scenario(n_triage=DEFAULT_N_TRIAGE+1) - + scenarios["triage+1"] = Scenario(n_triage=DEFAULT_N_TRIAGE + 1) + # extra examination capacity - scenarios['exam+1'] = Scenario(n_exam=DEFAULT_N_EXAM+1) - + scenarios["exam+1"] = Scenario(n_exam=DEFAULT_N_EXAM + 1) + # extra non-trauma treatment capacity - scenarios['treat+1'] = Scenario(n_cubicles_1=DEFAULT_N_CUBICLES_1+1) - + scenarios["treat+1"] = Scenario(n_cubicles_1=DEFAULT_N_CUBICLES_1 + 1) + # swap over 1 exam room for extra treat cubicle - scenarios['swap_exam_treat'] = Scenario(n_triage=DEFAULT_N_TRIAGE+1, - n_exam=DEFAULT_N_EXAM-1) - + scenarios["swap_exam_treat"] = Scenario( + n_triage=DEFAULT_N_TRIAGE + 1, n_exam=DEFAULT_N_EXAM - 1 + ) + # scenario + 3 min short mean exam times. - scenarios['short_exam'] = Scenario(n_triage=DEFAULT_N_TRIAGE+1, - n_exam=DEFAULT_N_EXAM-1, - exam_mean=12.0) + scenarios["short_exam"] = Scenario( + n_triage=DEFAULT_N_TRIAGE + 1, n_exam=DEFAULT_N_EXAM - 1, exam_mean=12.0 + ) return scenarios - -def run_scenario_analysis(scenarios, rc_period, n_reps): - ''' +def run_scenario_analysis( + scenarios: Dict[str, Scenario], rc_period: float, n_reps: int +) -> Dict[str, pd.DataFrame]: + """ Run each of the scenarios for a specified results collection period and replications. - + Params: ------ scenarios: dict dictionary of Scenario objects - + rc_period: float model run length - + n_rep: int Number of replications - - ''' - print('Scenario Analysis') - print(f'No. Scenario: {len(scenarios)}') - print(f'Replications: {n_reps}') + + Returns: + ------- + Dict + + """ + print("Scenario Analysis") + print(f"No. Scenario: {len(scenarios)}") + print(f"Replications: {n_reps}") scenario_results = {} for sc_name, scenario in scenarios.items(): - - print(f'Running {sc_name}', end=' => ') - replications = multiple_replications(scenario, rc_period=rc_period, - n_reps=n_reps) - print('done.\n') - - #save the results + + print(f"Running {sc_name}", end=" => ") + replications = multiple_replications( + scenario, rc_period=rc_period, n_reps=n_reps + ) + print("done.\n") + + # save the results scenario_results[sc_name] = replications - - print('Scenario analysis complete.') - return scenario_results + print("Scenario analysis complete.") + return scenario_results -def scenario_summary_frame(scenario_results): - ''' +def scenario_summary_frame(scenario_results: Dict[str, pd.DataFrame]) -> pd.DataFrame: + """ Mean results for each performance measure by scenario - + Parameters: ---------- scenario_results: dict - dictionary of replications. + dictionary of replications. Key identifies the performance measure - + Returns: ------- pd.DataFrame - ''' + """ columns = [] summary = pd.DataFrame() for sc_name, replications in scenario_results.items(): @@ -1148,7 +1230,3 @@ def scenario_summary_frame(scenario_results): summary.columns = columns return summary - - - -