Dev

DESCRIPTION

@@ -20,6 +20,5 @@ RoxygenNote: 7.3.2
 Imports:
     simmer
 Suggests:
-    testthat (>= 3.0.0),
-    lintr
+    devtools
 Config/testthat/edition: 3

NAMESPACE

@@ -0,0 +1,29 @@
+# Generated by roxygen2: do not edit by hand
+
+export(add_patient_generator)
+export(create_asu_arrivals)
+export(create_asu_los)
+export(create_asu_routing)
+export(create_asu_trajectory)
+export(create_parameters)
+export(create_rehab_arrivals)
+export(create_rehab_los)
+export(create_rehab_routing)
+export(create_rehab_trajectory)
+export(model)
+export(sample_routing)
+export(transform_to_lnorm)
+importFrom(simmer,add_generator)
+importFrom(simmer,branch)
+importFrom(simmer,get_attribute)
+importFrom(simmer,get_mon_arrivals)
+importFrom(simmer,get_mon_resources)
+importFrom(simmer,log_)
+importFrom(simmer,set_attribute)
+importFrom(simmer,simmer)
+importFrom(simmer,timeout)
+importFrom(simmer,trajectory)
+importFrom(simmer,wrap)
+importFrom(stats,rexp)
+importFrom(stats,rlnorm)
+importFrom(utils,capture.output)

R/add_patient_generator.R

@@ -0,0 +1,35 @@
+#' Add patient generator to Simmer environment.
+#'
+#' Creates a patient generator using an exponential inter-arrival distribution.
+#' The generator name is automatically constructed as \{unit\}_\{patient_type\}.
+#'
+#' @param env Simmer environment object. The simulation environment where
+#' generators will be added.
+#' @param trajectory Simmer trajectory object. Defines patient journey logic
+#' through the healthcare system.
+#' @param unit Character string specifying the care unit. Must be either "asu"
+#' (Acute Stroke Unit) or "rehab" (Rehabilitation Unit). Used to access correct
+#' parameter set and name the generator.
+#' @param patient_type Character string specifying patient category. Must be
+#' one of: "stroke", "tia", "neuro", or "other". Determines which arrival rate
+#' parameter is used.
+#' @param param Nested list containing simulation parameters. Must have
+#' structure \code{param$<unit>_arrivals$<patient_type>} containing numeric
+#' arrival intervals (e.g., \code{param$asu_arrivals$stroke = 10}).
+#'
+#' @importFrom simmer add_generator
+#' @importFrom stats rexp
+#'
+#' @return The modified Simmer environment with the new patient generator added.
+#' @export
+
+add_patient_generator <- function(env, trajectory, unit, patient_type, param) {
+  add_generator(
+    .env = env,
+    name_prefix = paste0(unit, "_", patient_type),
+    trajectory = trajectory,
+    distribution = function() {
+      rexp(1L, 1L / param[[paste0(unit, "_arrivals")]][[patient_type]])
+    }
+  )
+}

R/create_asu_trajectory.R

@@ -0,0 +1,94 @@
+#' Create acute stroke unit (ASU) patient trajectory.
+#'
+#' Represents patient stay in the ASU - samples their (1) destination after
+#' the ASU, and (2) length of stay (LOS) on the ASU.
+#'
+#' @param env Simmer environment object. The simulation environment where
+#' generators will be added.
+#' @param patient_type Character string specifying patient category. Must be
+#' one of: "stroke", "tia", "neuro", or "other". Determines which arrival rate
+#' parameter is used.
+#' @param param Nested list containing simulation parameters. Must have
+#' structure \code{param$asu_routing$<patient_type>} containing the probability
+#' of routing to each destination (e.g.
+#' \code{param$asu_routing$stroke$rehab = 0.24}).
+#'
+#' @importFrom simmer branch get_attribute log_ set_attribute timeout trajectory
+#' @importFrom stats rlnorm
+#'
+#' @return Simmer trajectory object. Defines patient journey logic through the
+#' healthcare system.
+#' @export
+
+create_asu_trajectory <- function(env, patient_type, param) {
+
+  # Set up simmer trajectory object...
+  trajectory(paste0("ASU_", patient_type, "_path")) |>
+
+    log_("🚶 Arrived at ASU") |>
+
+    # Sample destination after ASU (as destination influences length of stay)
+    set_attribute("post_asu_destination", function() {
+      sample_routing(prob_list = param[["asu_routing"]][[patient_type]])
+    }) |>
+
+    log_(function() {
+      # Retrieve attribute, and use to get post-ASU destination as a string
+      dest_index <- get_attribute(env, "post_asu_destination")
+      dest_names <- names(param[["asu_routing"]][[patient_type]])
+      dest <- dest_names[dest_index]
+      # Create log message
+      paste0("🎯 Planned ASU -> ", dest_index, " (", dest, ")")
+    }) |>
+
+    set_attribute("asu_los", function() {
+      # Retrieve attribute, and use to get post-ASU destination as a string
+      dest_index <- get_attribute(env, "post_asu_destination")
+      dest_names <- names(param[["asu_routing"]][[patient_type]])
+      dest <- dest_names[dest_index]
+
+      # Determine which LOS distribution to use
+      if (patient_type == "stroke") {
+        los_params <- switch(
+          dest,
+          esd = param[["asu_los_lnorm"]][["stroke_esd"]],
+          rehab = param[["asu_los_lnorm"]][["stroke_noesd"]],
+          other = param[["asu_los_lnorm"]][["stroke_mortality"]],
+          stop("Stroke post-asu destination '", dest, "' invalid",
+               call. = FALSE)
+        )
+      } else {
+        los_params <- param[["asu_los_lnorm"]][[patient_type]]
+      }
+
+      # Sample LOS from lognormal
+      rlnorm(
+        n = 1L,
+        meanlog = los_params[["meanlog"]],
+        sdlog = los_params[["sdlog"]]
+      )
+    }) |>
+
+    log_(function() {
+      paste0("⏳ ASU length of stay: ",
+             round(get_attribute(env, "asu_los"), 3L))
+    }) |>
+
+    timeout(function() get_attribute(env, "asu_los")) |>
+
+    log_("🏁 ASU stay completed") |>
+
+    # If that patient's destination is rehab, then start on that trajectory
+    branch(
+      option = function() {
+        # Retrieve attribute, and use to get post-ASU destination as a string
+        dest_index <- get_attribute(env, "post_asu_destination")
+        dest_names <- names(param[["asu_routing"]][[patient_type]])
+        dest <- dest_names[dest_index]
+        # Return 1 for rehab and 0 otherwise
+        if (dest == "rehab") 1L else 0L
+      },
+      continue = FALSE,  # Do not continue main trajectory after branch
+      create_rehab_trajectory(env, patient_type, param)
+    )
+}

R/create_rehab_trajectory.R

@@ -0,0 +1,79 @@
+#' Create rehab patient trajectory.
+#'
+#' Represents patient stay on the rehabilitation unit - samples their (1)
+#' destination after rehab, and (2) length of stay (LOS) on the unit.
+#'
+#' @param env Simmer environment object. The simulation environment where
+#' generators will be added.
+#' @param patient_type Character string specifying patient category. Must be
+#' one of: "stroke", "tia", "neuro", or "other". Determines which arrival rate
+#' parameter is used.
+#' @param param Nested list containing simulation parameters. Must have
+#' structure \code{param$rehab_routing$<patient_type>} containing the
+#' probability of routing to each destination (e.g.
+#' \code{param$rehab_routing$stroke$esd = 0.40}).
+#'
+#' @importFrom simmer get_attribute log_ set_attribute timeout trajectory
+#' @importFrom stats rlnorm
+#'
+#' @return Simmer trajectory object. Defines patient journey logic through the
+#' healthcare system.
+#' @export
+
+create_rehab_trajectory <- function(env, patient_type, param) {
+
+  # Set up simmer trajectory object...
+  trajectory(paste0("rehab_", patient_type, "_path")) |>
+
+    log_("🚶 Arrived at rehab") |>
+
+    # Sample destination after rehab (as destination influences length of stay)
+    set_attribute("post_rehab_destination", function() {
+      sample_routing(prob_list = param[["rehab_routing"]][[patient_type]])
+    }) |>
+
+    log_(function() {
+      # Retrieve attribute, and use to get post-rehab destination as a string
+      dest_index <- get_attribute(env, "post_rehab_destination")
+      dest_names <- names(param[["rehab_routing"]][[patient_type]])
+      dest <- dest_names[dest_index]
+      # Create log message
+      paste0("🎯 Planned rehab -> ", dest_index, " (", dest, ")")
+    }) |>
+
+    set_attribute("rehab_los", function() {
+      # Retrieve attribute, and use to get post-rehab destination as a string
+      dest_index <- get_attribute(env, "post_rehab_destination")
+      dest_names <- names(param[["rehab_routing"]][[patient_type]])
+      dest <- dest_names[dest_index]
+
+      # Determine which LOS distribution to use
+      if (patient_type == "stroke") {
+        los_params <- switch(
+          dest,
+          esd = param[["rehab_los_lnorm"]][["stroke_esd"]],
+          other = param[["rehab_los_lnorm"]][["stroke_noesd"]],
+          stop("Stroke post-rehab destination '", dest, "' invalid",
+               call. = FALSE)
+        )
+      } else {
+        los_params <- param[["rehab_los_lnorm"]][[patient_type]]
+      }
+
+      # Sample LOS from lognormal
+      rlnorm(
+        n = 1L,
+        meanlog = los_params[["meanlog"]],
+        sdlog = los_params[["sdlog"]]
+      )
+    }) |>
+
+    log_(function() {
+      paste0("⏳ Rehab length of stay: ",
+             round(get_attribute(env, "rehab_los"), 3L))
+    }) |>
+
+    timeout(function() get_attribute(env, "rehab_los")) |>
+
+    log_("🏁 Rehab stay completed")
+}

R/model.R

@@ -0,0 +1,86 @@
+#' Run simulation.
+#'
+#' @param run_number Integer representing index of current simulation run.
+#' @param param Named list of model parameters.
+#' @param set_seed Whether to set seed within the model function (which we
+#' may not wish to do if being set elsewhere - such as done in \code{runner()}).
+#' Default is TRUE.
+#'
+#' @importFrom simmer get_mon_arrivals get_mon_resources simmer wrap
+#' @importFrom utils capture.output
+#'
+#' @return TBC
+#' @export
+
+model <- function(run_number, param, set_seed = TRUE) {
+
+  # Set random seed based on run number
+  if (set_seed) {
+    set.seed(run_number)
+  }
+
+  # Determine whether to get verbose activity logs
+  verbose <- any(c(param[["log_to_console"]], param[["log_to_file"]]))
+
+  # Transform LOS parameters to lognormal scale
+  param[["asu_los_lnorm"]] <- transform_to_lnorm(param[["asu_los"]])
+  param[["rehab_los_lnorm"]] <- transform_to_lnorm(param[["rehab_los"]])
+
+  # Create simmer environment - set verbose to FALSE as using custom logs
+  env <- simmer("simulation", verbose = FALSE)
+
+  # Add ASU and rehab direct admission patient generators
+  for (unit in c("asu", "rehab")) {
+    for (patient_type in names(param[[paste0(unit, "_arrivals")]])) {
+
+      # Create patient trajectory
+      traj <- if (unit == "asu") {
+        create_asu_trajectory(env, patient_type, param)
+      } else {
+        create_rehab_trajectory(env, patient_type, param)
+      }
+
+      # Add patient generator using the created trajectory
+      sim_log <- capture.output(
+        env <- add_patient_generator( # nolint
+          env = env,
+          trajectory = traj,
+          unit = unit,
+          patient_type = patient_type,
+          param = param
+        )
+      )
+    }
+  }
+
+  # Run the model
+  sim_log <- capture.output(
+    env <- env |> # nolint
+      simmer::run(20L) |>
+      wrap()
+  )
+
+  # Save and/or display the log
+  if (isTRUE(verbose)) {
+    # Create full log message by adding parameters
+    param_string <- paste(names(param), param, sep = "=", collapse = ";\n ")
+    full_log <- append(c("Parameters:", param_string, "Log:"), sim_log)
+    # Print to console
+    if (isTRUE(param[["log_to_console"]])) {
+      print(full_log)
+    }
+    # Save to file
+    if (isTRUE(param[["log_to_file"]])) {
+      writeLines(full_log, param[["file_path"]])
+    }
+  }
+
+  # Extract the monitored arrivals and resources information from the simmer
+  # environment object
+  result <- list(
+    arrivals = get_mon_arrivals(env, per_resource = TRUE, ongoing = TRUE),
+    resources = get_mon_resources(env)
+  )
+
+  return(result)
+}