Merge pull request #17 from Boehringer-Ingelheim/update-analysis-vign…

…ette

Update for Analysis Vignette

.Rbuildignore

@@ -11,3 +11,4 @@ docs
 ^cran-comments\.md$
 ^cran_submission_script\.R$
 ^CRAN-SUBMISSION$
+^vignettes/*_files$

.github/workflows/linux.yml

@@ -24,6 +24,7 @@ jobs:
       R_KEEP_PKG_SOURCE: yes
 
     steps:
+      - uses: quarto-dev/quarto-actions/setup@v2
       - uses: actions/checkout@v2
 
       - uses: r-lib/actions/setup-pandoc@v2

.github/workflows/macos.yml

@@ -26,6 +26,7 @@ jobs:
       R_KEEP_PKG_SOURCE: yes
 
     steps:
+      - uses: quarto-dev/quarto-actions/setup@v2
       - uses: actions/checkout@v2
 
       - uses: r-lib/actions/setup-pandoc@v2

.github/workflows/pkgdown.yml

@@ -20,14 +20,15 @@ jobs:
     env:
       GITHUB_PAT: ${{ secrets.GITHUB_TOKEN }}
     steps:
+      - uses: quarto-dev/quarto-actions/setup@v2
       - uses: actions/checkout@v2
 
       - uses: r-lib/actions/setup-pandoc@v2
 
       - uses: r-lib/actions/setup-r@v2
         with:
           use-public-rspm: true
-      
+
       - name: Install JAGS
         run: |
           sudo apt-get update -y

.github/workflows/windows.yml

@@ -25,6 +25,7 @@ jobs:
       R_KEEP_PKG_SOURCE: yes
 
     steps:
+      - uses: quarto-dev/quarto-actions/setup@v2
       - uses: actions/checkout@v2
 
       - uses: r-lib/actions/setup-pandoc@v2
@@ -40,6 +41,10 @@ jobs:
           extra-packages: any::rcmdcheck
           needs: check
 
+      - name: Install BayesianMCPMod
+        shell: bash
+        run: R CMD INSTALL --preclean .
+
       - uses: r-lib/actions/check-r-package@v2
         with:
           upload-snapshots: true

DESCRIPTION

@@ -6,8 +6,8 @@ Authors@R: c(
     person("Boehringer Ingelheim Pharma GmbH & Co. KG", role = c("cph", "fnd")),
     person("Stephan", "Wojciekowski", , "[email protected]", role = c("aut", "cre")),
     person("Lars", "Andersen", , "[email protected]", role = "aut"),
-    person("Steven", "Brooks", , "steven.brooks@boehringer-ingelheim.com", role = "ctb"),
-    person("Sebastian", "Bossert", , "sebastian.bossert@boehringer-ingelheim.com", role = "aut")
+    person("Sebastian", "Bossert", , "sebastian.bossert@boehringer-ingelheim.com", role = "aut"),
+    person("Steven", "Brooks", , "steven.brooks@boehringer-ingelheim.com", role = "ctb")
   )
 Description: Bayesian MCPMod (Fleischer et al. (2022)
     <doi:10.1002/pst.2193>) is an innovative method that improves the
@@ -38,17 +38,19 @@ Imports:
     RBesT,
     stats
 Suggests:
+    reactable,
+    tibble,
+    quarto,
     clinDR,
     dplyr,
     knitr,
     rmarkdown,
     spelling,
     testthat (>= 3.0.0)
-VignetteBuilder: 
-    knitr
+VignetteBuilder: quarto
 Config/testthat/edition: 3
 Encoding: UTF-8
 Language: en-US
 LazyData: true
 Roxygen: list(markdown = TRUE)
-RoxygenNote: 7.2.3
+RoxygenNote: 7.3.1

NAMESPACE

@@ -9,7 +9,6 @@ S3method(print,postList)
 S3method(summary,postList)
 export(assessDesign)
 export(getBootstrapQuantiles)
-export(getBootstrapSamples)
 export(getContr)
 export(getCritProb)
 export(getESS)

R/BMCPMod.R

@@ -361,7 +361,7 @@ getCritProb <- function (
 #'                       crit_prob_adj  = critVal,
 #'                       simple         = FALSE)
 #' 
-#' @return Bayesian MCP test result as well as modelling result.
+#' @return Bayesian MCP test result as well as modeling result.
 #' 
 #' @export
 performBayesianMCPMod <- function (
@@ -399,7 +399,7 @@ performBayesianMCPMod <- function (
     stop ("Argument 'contr' must be of type 'optContr'")
 
   }
-  
+
   b_mcp <- performBayesianMCP(
     posterior_list = posterior_list,
     contr          = contr,
@@ -408,8 +408,8 @@ performBayesianMCPMod <- function (
   fits_list <- lapply(seq_along(posterior_list), function (i) {
 
     if (b_mcp[i, 1]) {
-      
-      sign_models <- b_mcp[i, -c(1, 2)] > attr(b_mcp, "critProbAdj")
+
+      sign_models <- b_mcp[i, -c(1, 2)] > attr(b_mcp, "crit_prob_adj")
 
       model_fits  <- getModelFits(
         models      = model_shapes,

R/bootstrapping.R

@@ -1,20 +1,22 @@
-#' @title getBootstrapSamples
+#' @title getBootstrapQuantiles
 #'
-#' @description A function for the calculation of bootstrapped model predictions. 
-#' Samples from the posterior distribution are drawn (via the RBesT function rmix()) and for every sample the simplified fitting step (see getModelFits() function) and a prediction is performed. 
-#' These fits are then used to identify the specified quantiles. 
+#' @description A function for the calculation of bootstrapped model predictions.
+#' Samples from the posterior distribution are drawn (via the RBesT function rmix()) and for every sample the simplified fitting step (see getModelFits() function) and a prediction is performed.
+#' These fits are then used to identify the specified quantiles.
 #' This approach can be considered as the Bayesian equivalent of the frequentist bootstrap approach described in O'Quigley et al. (2017).
 #' Instead of drawing n bootstrap samples from the sampling distribution of the trial dose-response estimates, here the samples are directly taken from the posterior distribution.
 #' @references O'Quigley J, Iasonos A, Bornkamp B. 2017. Handbook of Methods for Designing, Monitoring, and Analyzing Dose-Finding Trials (1st ed.). Chapman and Hall/CRC. doi:10.1201/9781315151984
-#' @param model_fits An object of class modelFits, i.e. information about fitted models & corresponding model coefficients as well as the posterior distribution that was the basis for the model fitting 
+#'
+#' @param model_fits An object of class modelFits, i.e. information about fitted models & corresponding model coefficients as well as the posterior distribution that was the basis for the model fitting
+#' @param quantiles A vector of quantiles that should be evaluated
 #' @param n_samples Number of samples that should be drawn as basis for the bootstrapped quantiles
 #' @param doses A vector of doses for which a prediction should be performed
 #' @param avg_fit Boolean variable, defining whether an average fit (based on generalized AIC weights) should be performed in addition to the individual models. Default TRUE.
 #'
 #' @examples
 #' posterior_list <- list(Ctrl = RBesT::mixnorm(comp1 = c(w = 1, m = 0, s = 1), sigma = 2),
 #'                        DG_1 = RBesT::mixnorm(comp1 = c(w = 1, m = 3, s = 1.2), sigma = 2),
-#'                        DG_2 = RBesT::mixnorm(comp1 = c(w = 1, m = 4, s = 1.5), sigma = 2) ,  
+#'                        DG_2 = RBesT::mixnorm(comp1 = c(w = 1, m = 4, s = 1.5), sigma = 2) ,
 #'                        DG_3 = RBesT::mixnorm(comp1 = c(w = 1, m = 6, s = 1.2), sigma = 2) ,
 #'                        DG_4 = RBesT::mixnorm(comp1 = c(w = 1, m = 6.5, s = 1.1), sigma = 2))
 #' models         <- c("exponential", "linear")
@@ -23,40 +25,43 @@
 #'                                posterior   = posterior_list,
 #'                                dose_levels = dose_levels,
 #'                                simple      = TRUE)
-#'                                
-#' getBootstrapSamples(model_fits = fit,
+#'
+#' getBootstrapQuantiles(model_fits = fit,
+#'                     quantiles  = c(0.025, 0.5, 0.975),
 #'                     n_samples  = 10, # speeding up example run time
 #'                     doses      = c(0, 6, 8))
-#' 
-#' @return  A data frame with columns for model, dose, and bootstrapped samples 
-#' 
+#'
+#' @return  A data frame with columns for model, dose, and bootstrapped samples
+#'
 #' @export
-getBootstrapSamples <- function (
-    
+getBootstrapQuantiles <- function (
+
   model_fits,
+  quantiles,
   n_samples = 1e3,
   doses     = NULL,
   avg_fit   = TRUE
-  
+
 ) {
-  
+
   mu_hat_samples <- sapply(attr(model_fits, "posterior"),
                            RBesT::rmix, n = n_samples)
   sd_hat         <- summary.postList(attr(model_fits, "posterior"))[, 2]
-  
+
   dose_levels    <- model_fits[[1]]$dose_levels
   model_names    <- names(model_fits)
-
+  quantile_probs <- sort(unique(quantiles))
+
   if (is.null(doses)) {
-    
+
     doses <- seq(min(dose_levels), max(dose_levels), length.out = 100L)
-    
+
   }
-  
+
   preds <- apply(mu_hat_samples, 1, function (mu_hat) {
-    
+
     preds_mu_hat <- sapply(model_names, function (model) {
-      
+
       fit <- DoseFinding::fitMod(
         dose  = model_fits[[1]]$dose_levels,
         resp  = mu_hat,
@@ -65,101 +70,50 @@ getBootstrapSamples <- function (
         type  = "general",
         bnds  = DoseFinding::defBnds(
           mD = max(model_fits[[1]]$dose_levels))[[model]])
-      
+
       preds <- stats::predict(fit, doseSeq = doses, predType = "ls-means")
       attr(preds, "gAIC") <- DoseFinding::gAIC(fit)
-      
+
       return (preds)
-      
+
     }, simplify = FALSE)
-    
+
     preds_mu_mat <- do.call(rbind, preds_mu_hat)
-    
+
     if (avg_fit) {
-      
+
       avg_fit_indx <- which.min(sapply(preds_mu_hat, attr, "gAIC"))
       preds_mu_mat <- rbind(preds_mu_mat, avgFit = preds_mu_mat[avg_fit_indx, ])
-      
+
     }
-    
+
     return (preds_mu_mat)
-    
+
   })
-  
+
   if (avg_fit) {
-    
+
     model_names <- c(model_names, "avgFit")
-    
+
   }
-
-  bs_samples <- as.data.frame(preds[order(rep(seq_along(model_names), length(doses))), ])
-  colnames(bs_samples) <- paste0("preds_", seq_len(n_samples))
-
-  bs_samples_data <- cbind(
-    model = rep(
-      x    = factor(model_names,
-                    levels = c("linear", "emax", "exponential",
-                               "sigEmax", "logistic", "quadratic",
-                               "avgFit")),
-      each = length(doses)),
-    dose  = doses,
-    bs_samples)
-
-  # tidyr::pivot_longer(bs_samples_data, cols = contains("preds"))
-
-  return (bs_samples_data)
-
-}
 
-#' @title getBootstrapQuantiles
-#'
-#' @description Calculates quantiles from bootstrapped dose predictions.
-#' Can be used to derive credible intervals to assess the uncertainty for the model fit.
-#' @param bs_samples An object of class bootstrappedSample as created by getBootstrapSamples
-#' @param quantiles A vector of quantiles that should be evaluated 
-#'
-#' @examples
-#' posterior_list <- list(Ctrl = RBesT::mixnorm(comp1 = c(w = 1, m = 0, s = 1), sigma = 2),
-#'                        DG_1 = RBesT::mixnorm(comp1 = c(w = 1, m = 3, s = 1.2), sigma = 2),
-#'                        DG_2 = RBesT::mixnorm(comp1 = c(w = 1, m = 4, s = 1.5), sigma = 2) ,  
-#'                        DG_3 = RBesT::mixnorm(comp1 = c(w = 1, m = 6, s = 1.2), sigma = 2) ,
-#'                        DG_4 = RBesT::mixnorm(comp1 = c(w = 1, m = 6.5, s = 1.1), sigma = 2))
-#' models         <- c("exponential", "linear")
-#' dose_levels    <- c(0, 1, 2, 4, 8)
-#' fit            <- getModelFits(models      = models,
-#'                                posterior   = posterior_list,
-#'                                dose_levels = dose_levels,
-#'                                simple      = TRUE)
-#'                                
-#' bs_samples     <- getBootstrapSamples(model_fits = fit,
-#'                                       n_samples  = 10, # speeding up example run time
-#'                                       doses      = c(0, 6, 8))
-#'                            
-#' getBootstrapQuantiles(bs_samples = bs_samples,
-#'                       quantiles  = c(0.025, 0.5, 0.975))
-#' @return  A data frame with entries doses, models, and quantiles 
-#' 
-#' @export
-getBootstrapQuantiles <- function (
+  sort_indx <- order(rep(seq_along(model_names), length(doses)))
+  quant_mat <- t(apply(X      = preds[sort_indx, ],
+                       MARGIN = 1,
+                       FUN    = stats::quantile,
+                       probs  = quantile_probs))
 
-  bs_samples,
-  quantiles
 
-) {
-
-  quantile_probs <- sort(unique(quantiles))
-
-  bs_quantiles <- t(apply(X      = bs_samples[, -c(1, 2)],
-                          MARGIN = 1,
-                          FUN    = stats::quantile,
-                          probs  = quantile_probs))
-
-  bs_quantiles_data <- cbind(
-    bs_samples[, c(1, 2)],
-    as.data.frame(bs_quantiles))
-
-  return (bs_quantiles_data)
+  cr_bounds_data <- cbind(
+    doses  = doses,
+    models = rep(
+      x    = factor(model_names, levels = model_names),
+      each = length(doses)),
+    as.data.frame(quant_mat))
+
+  return (cr_bounds_data)
 
 }
 
 
+