Dev

.Rbuildignore

@@ -10,4 +10,5 @@
 ^\.lintr$
 ^outputs$
 ^rmarkdown$
-^CITATION\.cff$
+^CITATION\.cff$
+^run_rmarkdown\.sh$

.github/workflows/R-CMD-check.yaml

@@ -2,8 +2,7 @@
 # Need help debugging build failures? Start at https://github.com/r-lib/actions#where-to-find-help
 on:
   push:
-    branches: [main, master]
-  pull_request:
+    branches: [main]
   workflow_dispatch:
 
 name: R-CMD-check.yaml

.github/workflows/lint.yaml

@@ -4,7 +4,6 @@
 on:
   push:
     branches: [main]
-  pull_request:
   workflow_dispatch:
 
 name: lint

.gitignore

@@ -1,4 +1,49 @@
-.Rproj.user
+# History files
 .Rhistory
+.Rapp.history
+
+# Session Data files
 .RData
+.RDataTmp
+
+# User-specific files
 .Ruserdata
+
+# Example code in package build process
+*-Ex.R
+
+# Output files from R CMD build
+/*.tar.gz
+
+# Output files from R CMD check
+/*.Rcheck/
+
+# RStudio files
+.Rproj.user/
+
+# produced vignettes
+vignettes/*.html
+vignettes/*.pdf
+
+# OAuth2 token, see https://github.com/hadley/httr/releases/tag/v0.3
+.httr-oauth
+
+# knitr and R markdown default cache directories
+*_cache/
+/cache/
+
+# Temporary files created by R markdown
+*.utf8.md
+*.knit.md
+
+# R Environment Variables
+.Renviron
+
+# pkgdown site
+docs/
+
+# translation temp files
+po/*~
+
+# RStudio Connect folder
+rsconnect/

DESCRIPTION

@@ -23,15 +23,15 @@ Imports:
     purrr,
     rlang,
     tidyr,
-    R6,
     tidyselect,
     future,
-    future.apply
+    future.apply,
+    ggplot2
 Suggests:
     testthat (>= 3.0.0),
+    patrick,
     lintr,
     devtools,
-    ggplot2,
     xtable,
     data.table
 Config/testthat/edition: 3

NAMESPACE

@@ -1,29 +1,42 @@
 # Generated by roxygen2: do not edit by hand
 
-export(Defaults)
-export(defaults)
+export(confidence_interval_method)
+export(get_run_results)
 export(model)
-export(process_replications)
-export(trial)
+export(parameters)
+export(run_scenarios)
+export(runner)
 export(valid_inputs)
-importFrom(R6,R6Class)
+importFrom(dplyr,filter)
 importFrom(dplyr,full_join)
 importFrom(dplyr,group_by)
 importFrom(dplyr,lead)
 importFrom(dplyr,mutate)
 importFrom(dplyr,n_distinct)
+importFrom(dplyr,pull)
+importFrom(dplyr,select)
+importFrom(dplyr,slice_head)
 importFrom(dplyr,summarise)
 importFrom(future,multisession)
 importFrom(future,plan)
 importFrom(future,sequential)
 importFrom(future.apply,future_lapply)
+importFrom(ggplot2,aes)
+importFrom(ggplot2,geom_line)
+importFrom(ggplot2,geom_ribbon)
+importFrom(ggplot2,geom_vline)
+importFrom(ggplot2,ggplot)
+importFrom(ggplot2,ggsave)
+importFrom(ggplot2,labs)
+importFrom(ggplot2,theme_minimal)
 importFrom(magrittr,"%>%")
 importFrom(purrr,reduce)
 importFrom(rlang,.data)
 importFrom(simmer,add_generator)
 importFrom(simmer,add_resource)
 importFrom(simmer,get_mon_arrivals)
 importFrom(simmer,get_mon_resources)
+importFrom(simmer,now)
 importFrom(simmer,release)
 importFrom(simmer,run)
 importFrom(simmer,seize)
@@ -32,6 +45,8 @@ importFrom(simmer,timeout)
 importFrom(simmer,trajectory)
 importFrom(simmer,wrap)
 importFrom(stats,rexp)
+importFrom(stats,sd)
+importFrom(stats,t.test)
 importFrom(tidyr,drop_na)
 importFrom(tidyr,pivot_wider)
 importFrom(tidyselect,any_of)

R/choose_replications.R

@@ -0,0 +1,113 @@
+#' Use the confidence interval method to select the number of replications.
+#'
+#' @param replications Number of times to run the model.
+#' @param desired_precision Desired mean deviation from confidence interval.
+#' @param metric Name of performance metric to assess.
+#' @param yaxis_title Label for y axis.
+#' @param path Path inc. filename to save figure to.
+#' @param min_rep A suggested minimum number of replications (default=NULL).
+#'
+#' @importFrom stats sd t.test
+#' @importFrom dplyr filter slice_head select pull
+#' @importFrom ggplot2 ggplot aes geom_line geom_ribbon geom_vline labs
+#' theme_minimal ggsave
+#' @importFrom rlang .data
+#'
+#' @return Dataframe with results from each replication.
+#' @export
+
+confidence_interval_method <- function(replications, desired_precision, metric,
+                                       yaxis_title, path, min_rep = NULL) {
+  # Run model for specified number of replications
+  param <- parameters(number_of_runs = replications)
+  raw_results <- runner(param)
+  results <- get_run_results(raw_results)
+
+  # If mean of metric is less than 1, multiply by 100
+  if (mean(results[[metric]]) < 1L) {
+    results[[paste0("adj_", metric)]] <- results[[metric]] * 100L
+    metric <- paste0("adj_", metric)
+  }
+
+  # Initialise list to store the results
+  cumulative_list <- list()
+
+  # For each row in the dataframe, filter to rows up to the i-th replication
+  # then perform calculations
+  for (i in 1L:replications) {
+
+    # Filter rows up to the i-th replication
+    subset <- results[[metric]][1L:i]
+
+    # Calculate mean
+    mean <- mean(subset)
+
+    # Some calculations require more than 1 observation else will error...
+    if (i == 1L) {
+      # When only one observation, set to NA
+      std_dev <- NA
+      ci_lower <- NA
+      ci_upper <- NA
+      deviation <- NA
+    } else {
+      # Else, calculate standard deviation, 95% confidence interval, and
+      # percentage deviation
+      std_dev <- sd(subset)
+      ci <- t.test(subset)[["conf.int"]]
+      ci_lower <- ci[[1L]]
+      ci_upper <- ci[[2L]]
+      deviation <- ((ci_upper - mean) / mean) * 100L
+    }
+
+    # Append to the cumulative list
+    cumulative_list[[i]] <- data.frame(
+      replications = i,
+      cumulative_mean = mean,
+      cumulative_std = std_dev,
+      ci_lower = ci_lower,
+      ci_upper = ci_upper,
+      perc_deviation = deviation
+    )
+  }
+
+  # Combine the list into a single data frame
+  cumulative <- do.call(rbind, cumulative_list)
+
+  # Get the minimum number of replications where deviation is less than target
+  compare <- cumulative %>%
+    filter(.data[["perc_deviation"]] <= desired_precision * 100L)
+  if (nrow(compare) > 0L) {
+    # Get minimum number
+    n_reps <- compare %>%
+      slice_head() %>%
+      dplyr::select(replications) %>%
+      pull()
+    print(paste0("Reached desired precision (", desired_precision, ") in ",
+                 n_reps, " replications."))
+  } else {
+    warning("Running ", replications, " replications did not reach ",
+            "desired precision (", desired_precision, ").")
+  }
+
+  # Plot the cumulative mean and confidence interval
+  p <- ggplot(cumulative, aes(x = .data[["replications"]],
+                              y = .data[["cumulative_mean"]])) +
+    geom_line() +
+    geom_ribbon(aes(ymin = ci_lower, ymax = ci_upper), alpha = 0.2)
+
+  # If specified, plot the minimum suggested number of replications
+  if (!is.null(min_rep)) {
+    p <- p +
+      geom_vline(xintercept = min_rep, linetype = "dashed", color = "red")
+  }
+
+  # Modify labels and style
+  p <- p +
+    labs(x = "Replications", y = yaxis_title) +
+    theme_minimal()
+
+  # Save the plot
+  ggsave(filename = path, width = 6.5, height = 4L, bg = "white")
+
+  return(cumulative)
+}

R/process_replications.R → R/get_run_results.R

@@ -1,4 +1,4 @@
-#' Process results from each replication.
+#' Get results from each replication.
 #'
 #' For each replication (there can be one or many), calculate the: (1) number
 #' of arrivals, (2) mean wait time for each resource, (3) mean activity time
@@ -22,14 +22,14 @@
 #' @importFrom dplyr group_by summarise n_distinct mutate lead full_join
 #' @importFrom purrr reduce
 #' @importFrom rlang .data
-#' @importFrom simmer get_mon_resources get_mon_arrivals
+#' @importFrom simmer get_mon_resources get_mon_arrivals now
 #' @importFrom tidyr pivot_wider drop_na
 #' @importFrom tidyselect any_of
 #'
 #' @return Tibble with results from each replication.
 #' @export
 
-process_replications <- function(results) {
+get_run_results <- function(results) {
 
   # Remove patients who were still waiting and had not completed
   results[["arrivals"]] <- results[["arrivals"]] %>%