New commit for additional statistics

R/Design-methods.R

@@ -417,6 +417,19 @@ setMethod("simulate",
           data = thisData
         )
 
+        # Get the MTD estimate from the samples.
+        sample_mtd <- dose(
+          mean(object@nextBest@target),
+          model = object@model,
+          samples = thisSamples
+        )
+
+        # Create list with median MTD and CV.
+        additional_stats <- list(
+          median_mtd = median(sample_mtd),
+          cv_mtd = mad(sample_mtd) / median(sample_mtd)
+        )
+
         ## return the results
         thisResult <-
           list(
@@ -431,7 +444,8 @@ setMethod("simulate",
                 stopit,
                 "message"
               ),
-            report_results = stopit_results
+            report_results = stopit_results,
+            additional_stats = additional_stats
           )
         return(thisResult)
       }
@@ -470,13 +484,17 @@ setMethod("simulate",
       stopResults <- lapply(resultList, "[[", "report_results")
       stop_matrix <- as.matrix(do.call(rbind, stopResults))
 
+      # Result list of additional stats such as median MTD and mean CV MTD values.
+      additional_stats <- lapply(resultList, "[[", "additional_stats")
+
       ## return the results in the Simulations class object
       ret <- Simulations(
         data = dataList,
         doses = recommendedDoses,
         fit = fitList,
         stop_report = stop_matrix,
         stop_reasons = stopReasons,
+        additional_stats = additional_stats,
         seed = RNGstate
       )
 
@@ -1075,6 +1093,9 @@ setMethod("simulate",
       ## for dual simulations as it would fail in summary otherwise (for dual simulations reporting is not implemented)
       stop_report <- matrix(TRUE, nrow = nsim)
 
+      ## For dual simulations summary of additional statistics like median MTD and mean CV MTD
+      additional_stats <- vector(mode = "list", length = nsim)
+
       ## return the results in the DualSimulations class object
       ret <- DualSimulations(
         data = dataList,
@@ -1085,6 +1106,7 @@ setMethod("simulate",
         fit_biomarker = fitBiomarkerList,
         stop_report = stop_report,
         stop_reasons = stopReasons,
+        additional_stats = additional_stats,
         seed = RNGstate
       )
 

R/Simulations-class.R

@@ -107,6 +107,7 @@ GeneralSimulations <- function(data,
 #' @slot fit (`list`)\cr final fits
 #' @slot stop_reasons (`list`)\cr stopping reasons for each simulation run
 #' @slot stop_report matrix of stopping rule outcomes
+#' @slot additional_stats list of median MTD and mean CV MTD values
 #' @aliases Simulations
 #' @export
 .Simulations <-
@@ -115,7 +116,8 @@ GeneralSimulations <- function(data,
     slots = c(
       fit = "list",
       stop_report = "matrix",
-      stop_reasons = "list"
+      stop_reasons = "list",
+      additional_stats = "list"
     ),
     prototype = prototype(
       fit =
@@ -125,7 +127,9 @@ GeneralSimulations <- function(data,
         ),
       stop_report = matrix(TRUE, nrow = 2),
       stop_reasons =
-        list("A", "A")
+        list("A", "A"),
+      additional_stats =
+        list(1, 0.4)
     ),
     contains = "GeneralSimulations",
     validity =
@@ -153,6 +157,12 @@ GeneralSimulations <- function(data,
           "stop_report must be a matrix of mode logical in which the number of rows equals the number of simulations
       and which must not contain any missing values"
         )
+
+        o$check(
+          identical(length(object@additional_stats), nSims),
+          "additional_stats must have same length as data"
+        )
+
         o$result()
       }
   )
@@ -164,19 +174,22 @@ GeneralSimulations <- function(data,
 #' @param fit (`list`)\cr see slot definition.
 #' @param stop_reasons (`list`)\cr see slot definition.
 #' @param stop_report see [`Simulations`]
+#' @param additional_stats (`list`)\cr see slot definition.
 #' @param \dots additional parameters from [`GeneralSimulations`]
 #'
 #' @example examples/Simulations-class-Simulations.R
 #' @export
 Simulations <- function(fit,
                         stop_reasons,
                         stop_report,
+                        additional_stats,
                         ...) {
   start <- GeneralSimulations(...)
   .Simulations(start,
     fit = fit,
     stop_report = stop_report,
-    stop_reasons = stop_reasons
+    stop_reasons = stop_reasons,
+    additional_stats = additional_stats
   )
 }
 
@@ -331,6 +344,7 @@ DualSimulations <- function(rho_est,
     representation(
       stop_report = "matrix",
       fitAtDoseMostSelected = "numeric",
+      additional_stats = "list",
       meanFit = "list"
     ),
     contains = "GeneralSimulationsSummary"

R/Simulations-methods.R

@@ -634,6 +634,7 @@ setMethod("summary",
       ret <- .SimulationsSummary(
         start,
         stop_report = object@stop_report,
+        additional_stats = object@additional_stats,
         fitAtDoseMostSelected = fitAtDoseMostSelected,
         meanFit = meanFit
       )
@@ -985,6 +986,17 @@ setMethod("show",
         "Fitted toxicity rate at dose most often selected"
       )
 
+      # Report results of additional statistics such as median MTD and mean CV MTD
+
+      cat(
+        "Median MTD : ",
+        round(unlist(additional_stats)[1], 2), "\n"
+      )
+
+      cat(
+        "Mean CV (MTD) [%] : ",
+        round(unlist(additional_stats)[2] * 100, 2), "\n"
+      )
 
       # Report individual stopping rules with non-<NA> labels.
 

tests/testthat/test-Design-methods.R

@@ -1,5 +1,39 @@
 # simulate ----
 
+test_that("Test if simulate generate the expected output.", {
+  data <- h_get_data(placebo = FALSE)
+  model <- h_get_logistic_normal()
+  increments <- h_increments_relative()
+  next_best <- h_next_best_ncrm()
+  size <- CohortSizeConst(size = 3)
+
+  # Extreme truth function, which has constant probability 1 in dose grid range.
+  truth <- probFunction(model, alpha0 = 175, alpha1 = 5)
+  stop_rule <- StoppingMinPatients(nPatients = 5)
+  design <- Design(
+    model = model,
+    stopping = stop_rule,
+    increments = increments,
+    nextBest = next_best,
+    cohort_size = size,
+    data = data,
+    startingDose = 25
+  )
+
+  my_options <- McmcOptions(burnin = 100, step = 2, samples = 5, rng_kind = "Mersenne-Twister", rng_seed = 3)
+
+  sim <- simulate(
+    design,
+    nsim = 1,
+    truth = truth,
+    seed = 819,
+    mcmcOptions = my_options
+  )
+
+  expect_snapshot(sim)
+})
+
+
 ## NextBestInfTheory ----
 
 test_that("NextBestInfTheory produces consistent results for empty data", {