From e7c74250594f3d1acc9bc9f298ed35fea0426360 Mon Sep 17 00:00:00 2001 From: Dominic Tong Date: Tue, 2 Dec 2025 01:56:52 +0000 Subject: [PATCH 1/7] initial commit with dose_cmts --- DESCRIPTION | 2 +- R/parse_input_data.R | 9 +++++++-- R/run_eval.R | 2 ++ 3 files changed, 10 insertions(+), 3 deletions(-) diff --git a/DESCRIPTION b/DESCRIPTION index 5e0ae21..ca32725 100644 --- a/DESCRIPTION +++ b/DESCRIPTION @@ -42,6 +42,6 @@ Remotes: Config/testthat/edition: 3 Encoding: UTF-8 Roxygen: list(markdown = TRUE) -RoxygenNote: 7.3.2 +RoxygenNote: 7.3.3 VignetteBuilder: knitr LazyData: true diff --git a/R/parse_input_data.R b/R/parse_input_data.R index 3264cf9..0a5c4d3 100644 --- a/R/parse_input_data.R +++ b/R/parse_input_data.R @@ -9,6 +9,7 @@ parse_input_data <- function( data, covariates, ids, + dose_cmts, group = NULL, verbose = TRUE ) { @@ -31,6 +32,7 @@ parse_input_data <- function( new_data, parse_nm_data, covariates = covariates, + dose_cmts = dose_cmts, group = group ) @@ -39,10 +41,13 @@ parse_input_data <- function( out } -parse_nm_data <- function(data, covariates, group = NULL) { +parse_nm_data <- function(data, covariates, dose_cmts, group = NULL) { out <- list() out$covariates <- make_covariates_object(data, covariates) - out$regimen <- PKPDsim::nm_to_regimen(data) + out$regimen <- PKPDsim::nm_to_regimen( + data, + dose_cmts = dose_cmts + ) obs <- data |> dplyr::filter(.data$EVID == 0) if(!is.null(group)) { diff --git a/R/run_eval.R b/R/run_eval.R index 23d738f..61e9348 100644 --- a/R/run_eval.R +++ b/R/run_eval.R @@ -49,6 +49,7 @@ run_eval <- function( model, data, + dose_cmts, ids = NULL, parameters = NULL, fixed = NULL, @@ -102,6 +103,7 @@ run_eval <- function( data = input_data, covariates = covariates, ids = ids, + dose_cmts = dose_cmts, group = group, verbose = verbose ) From 82a607c8158478ddb93dc59a968cf515011499f3 Mon Sep 17 00:00:00 2001 From: Dominic Tong Date: Wed, 3 Dec 2025 00:45:11 +0000 Subject: [PATCH 2/7] rm dose_cmts --- R/parse_input_data.R | 9 ++------- R/run_eval.R | 2 -- 2 files changed, 2 insertions(+), 9 deletions(-) diff --git a/R/parse_input_data.R b/R/parse_input_data.R index 0a5c4d3..3264cf9 100644 --- a/R/parse_input_data.R +++ b/R/parse_input_data.R @@ -9,7 +9,6 @@ parse_input_data <- function( data, covariates, ids, - dose_cmts, group = NULL, verbose = TRUE ) { @@ -32,7 +31,6 @@ parse_input_data <- function( new_data, parse_nm_data, covariates = covariates, - dose_cmts = dose_cmts, group = group ) @@ -41,13 +39,10 @@ parse_input_data <- function( out } -parse_nm_data <- function(data, covariates, dose_cmts, group = NULL) { +parse_nm_data <- function(data, covariates, group = NULL) { out <- list() out$covariates <- make_covariates_object(data, covariates) - out$regimen <- PKPDsim::nm_to_regimen( - data, - dose_cmts = dose_cmts - ) + out$regimen <- PKPDsim::nm_to_regimen(data) obs <- data |> dplyr::filter(.data$EVID == 0) if(!is.null(group)) { diff --git a/R/run_eval.R b/R/run_eval.R index 61e9348..23d738f 100644 --- a/R/run_eval.R +++ b/R/run_eval.R @@ -49,7 +49,6 @@ run_eval <- function( model, data, - dose_cmts, ids = NULL, parameters = NULL, fixed = NULL, @@ -103,7 +102,6 @@ run_eval <- function( data = input_data, covariates = covariates, ids = ids, - dose_cmts = dose_cmts, group = group, verbose = verbose ) From cb2583ba39983678b5e8d46a267be42ab6fdbc09 Mon Sep 17 00:00:00 2001 From: Ron Keizer Date: Wed, 3 Dec 2025 18:09:06 +0000 Subject: [PATCH 3/7] skip in case of fit error --- R/calculate_stats.R | 1 + R/run_eval_core.R | 47 ++++++++++++++++++++++++++++++--------------- 2 files changed, 33 insertions(+), 15 deletions(-) diff --git a/R/calculate_stats.R b/R/calculate_stats.R index fd94c7f..776488b 100644 --- a/R/calculate_stats.R +++ b/R/calculate_stats.R @@ -25,6 +25,7 @@ calculate_stats <- function( cols = c("pred", "map_ipred", "iter_ipred"), names_to = "type" ) |> dplyr::group_by(.data$type, .data$apriori) |> + dplyr::filter(!is.na(value)) |> dplyr::summarise( rmse = rmse(.data$dv, .data$value), nrmse = nrmse(.data$dv, .data$value), diff --git a/R/run_eval_core.R b/R/run_eval_core.R index 08bbdf9..be62a22 100644 --- a/R/run_eval_core.R +++ b/R/run_eval_core.R @@ -67,22 +67,39 @@ run_eval_core <- function( iov_bins = mod_obj$bins, verbose = FALSE ) + if(inherits(fit, "error")) { + ## create NA records for this fit + pred_data <- tibble::tibble( + id = obs_data$id, + t = obs_data$t, + dv = NA, + ipred = NA, + pred = NA, + ofv = NA, + ss_w = NA, + `_iteration` = iterations[i], + `_grouper` = obs_data$`_grouper` + ) + par_dummy <- as.data.frame(mod_upd$parameters) + par_dummy[, 1:ncol(par_dummy)] <- NA + fit_pars <- dplyr::mutate(as.data.frame(par_dummy), id = obs_data$id[1]) + } else { + ## Data frame with predictive data + pred_data <- tibble::tibble( + id = obs_data$id, + t = obs_data$t, + dv = fit$dv, + ipred = fit$ipred, + pred = fit$pred, + ofv = fit$fit$value, + ss_w = ss(fit$dv, fit$ipred, weights), + `_iteration` = iterations[i], + `_grouper` = obs_data$`_grouper` + ) + ## Add parameter estimates + fit_pars <- dplyr::mutate(as.data.frame(fit$parameters), id = obs_data$id[1]) + } - ## Data frame with predictive data - pred_data <- tibble::tibble( - id = obs_data$id, - t = obs_data$t, - dv = fit$dv, - ipred = fit$ipred, - pred = fit$pred, - ofv = fit$fit$value, - ss_w = ss(fit$dv, fit$ipred, weights), - `_iteration` = iterations[i], - `_grouper` = obs_data$`_grouper` - ) - - ## Add parameter estimates - fit_pars <- dplyr::mutate(as.data.frame(fit$parameters), id = obs_data$id[1]) comb <- dplyr::bind_rows( comb, dplyr::left_join(pred_data, fit_pars, by = "id") ) From 20f5d2aff584a1a284450991c9849d03bbdd23d2 Mon Sep 17 00:00:00 2001 From: Ron Keizer Date: Wed, 3 Dec 2025 18:22:32 +0000 Subject: [PATCH 4/7] add warning when errors --- R/calculate_stats.R | 7 ++++++- 1 file changed, 6 insertions(+), 1 deletion(-) diff --git a/R/calculate_stats.R b/R/calculate_stats.R index 776488b..45939c6 100644 --- a/R/calculate_stats.R +++ b/R/calculate_stats.R @@ -20,12 +20,17 @@ calculate_stats <- function( if(inherits(res, "mipdeval_results")) { res <- res$results } + ## Check for errors during fits / predictions + errors <- res |> + dplyr::filter(is.na(pred) | (is.na(map_ipred) & !apriori) | is.na(iter_ipred)) + if(nrow(errors) > 0) { + cli::cli_warn("Errors where encountered in {nrow(errors)} out of {nrow(res)} evaluated predictions. The problems occurred in patient(s) {unique(errors$id)}.") + } out <- res |> tidyr::pivot_longer( cols = c("pred", "map_ipred", "iter_ipred"), names_to = "type" ) |> dplyr::group_by(.data$type, .data$apriori) |> - dplyr::filter(!is.na(value)) |> dplyr::summarise( rmse = rmse(.data$dv, .data$value), nrmse = nrmse(.data$dv, .data$value), From 02be7896b9e0f2673cfc326f073c021d4b794529 Mon Sep 17 00:00:00 2001 From: dominic-irx <59584945+dominic-irx@users.noreply.github.com> Date: Wed, 3 Dec 2025 10:45:20 -0800 Subject: [PATCH 5/7] fix typo --- R/calculate_stats.R | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/R/calculate_stats.R b/R/calculate_stats.R index 45939c6..a01c462 100644 --- a/R/calculate_stats.R +++ b/R/calculate_stats.R @@ -24,7 +24,7 @@ calculate_stats <- function( errors <- res |> dplyr::filter(is.na(pred) | (is.na(map_ipred) & !apriori) | is.na(iter_ipred)) if(nrow(errors) > 0) { - cli::cli_warn("Errors where encountered in {nrow(errors)} out of {nrow(res)} evaluated predictions. The problems occurred in patient(s) {unique(errors$id)}.") + cli::cli_warn("Errors were encountered in {nrow(errors)} out of {nrow(res)} evaluated predictions. The problems occurred in patient(s) {unique(errors$id)}.") } out <- res |> tidyr::pivot_longer( From a0a4f9a39e2cc88ddb3a3e26dd9125edb3fe42a0 Mon Sep 17 00:00:00 2001 From: Ron Keizer Date: Wed, 3 Dec 2025 20:35:04 +0000 Subject: [PATCH 6/7] add arg to pass MAP options to get_map_estimates() --- R/run_eval.R | 6 ++++++ R/run_eval_core.R | 37 ++++++++++++++++++++++--------------- man/run_eval.Rd | 6 ++++++ man/run_eval_core.Rd | 8 +++++++- 4 files changed, 41 insertions(+), 16 deletions(-) diff --git a/R/run_eval.R b/R/run_eval.R index 23d738f..c190463 100644 --- a/R/run_eval.R +++ b/R/run_eval.R @@ -32,6 +32,10 @@ #' result from a call to [stats_summ_options()]. #' @param .vpc_options Options for VPC simulations. This must be the result from #' a call to [vpc_options()]. +#' @param .fit_options Options for controlling MAP Bayesian fit. E.g. +#' `.map_options = list(control = list(reltol = 1e-04))` will perform a slightly +#' coarser fit (but more stable) which can be useful in some case (default +#' `reltol = 1e-05`). #' @param threads number of threads to divide computations on. Default is 1, #' i.e. no parallel execution #' @param ruv residual error variability magnitude, specified as list. @@ -63,6 +67,7 @@ run_eval <- function( incremental = FALSE, .stats_summ_options = stats_summ_options(), .vpc_options = vpc_options(), + .fit_options = list(), threads = 1, progress = TRUE, verbose = TRUE @@ -127,6 +132,7 @@ run_eval <- function( weight_prior = weight_prior, incremental = incremental, progress_function = p, + .fit_options = .fit_options, .threads = threads, .skip = .vpc_options$vpc_only ) diff --git a/R/run_eval_core.R b/R/run_eval_core.R index be62a22..42eec80 100644 --- a/R/run_eval_core.R +++ b/R/run_eval_core.R @@ -14,7 +14,8 @@ run_eval_core <- function( weight_prior = 1, censor_covariates = TRUE, incremental = FALSE, - progress_function = function() {} + progress_function = function() {}, + .fit_options = NULL ) { progress_function() @@ -52,20 +53,26 @@ run_eval_core <- function( mod_upd$parameters <- fit$parameters # take params from previous fit mod_upd$omega <- fit$vcov } - fit <- PKPDmap::get_map_estimates( - model = mod_obj$model, - parameters = mod_upd$parameters, - omega = mod_upd$omega, - error = mod_obj$ruv, - fixed = mod_obj$fixed, - as_eta = mod_obj$kappa, - data = data$observations, - covariates = cov_data, - regimen = data$regimen, - weight_prior = weight_prior, - weights = weights, - iov_bins = mod_obj$bins, - verbose = FALSE + fit <- do.call( + PKPDmap::get_map_estimates, + c( + list( + model = mod_obj$model, + parameters = mod_upd$parameters, + omega = mod_upd$omega, + error = mod_obj$ruv, + fixed = mod_obj$fixed, + as_eta = mod_obj$kappa, + data = data$observations, + covariates = cov_data, + regimen = data$regimen, + weight_prior = weight_prior, + weights = weights, + iov_bins = mod_obj$bins, + verbose = FALSE + ), + .fit_options + ) ) if(inherits(fit, "error")) { ## create NA records for this fit diff --git a/man/run_eval.Rd b/man/run_eval.Rd index 1c49a8d..73a71a9 100644 --- a/man/run_eval.Rd +++ b/man/run_eval.Rd @@ -21,6 +21,7 @@ run_eval( incremental = FALSE, .stats_summ_options = stats_summ_options(), .vpc_options = vpc_options(), + .fit_options = list(), threads = 1, progress = TRUE, verbose = TRUE @@ -88,6 +89,11 @@ result from a call to \code{\link[=stats_summ_options]{stats_summ_options()}}.} \item{.vpc_options}{Options for VPC simulations. This must be the result from a call to \code{\link[=vpc_options]{vpc_options()}}.} +\item{.fit_options}{Options for controlling MAP Bayesian fit. E.g. +\code{.map_options = list(control = list(reltol = 1e-04))} will perform a slightly +coarser fit (but more stable) which can be useful in some case (default +\code{reltol = 1e-05}).} + \item{threads}{number of threads to divide computations on. Default is 1, i.e. no parallel execution} diff --git a/man/run_eval_core.Rd b/man/run_eval_core.Rd index 23de731..2314866 100644 --- a/man/run_eval_core.Rd +++ b/man/run_eval_core.Rd @@ -13,7 +13,8 @@ run_eval_core( censor_covariates = TRUE, incremental = FALSE, progress_function = function() { - } + }, + .fit_options = NULL ) } \arguments{ @@ -46,6 +47,11 @@ approach has been called "model predictive control (MPC)" "regular" MAP in some scenarios. Default is \code{FALSE}.} \item{progress_function}{function to increment progress bar} + +\item{.fit_options}{Options for controlling MAP Bayesian fit. E.g. +\code{.map_options = list(control = list(reltol = 1e-04))} will perform a slightly +coarser fit (but more stable) which can be useful in some case (default +\code{reltol = 1e-05}).} } \value{ a \code{data.frame} with individual predictions From 428fffdc266ebd63afe5ff3e964e32857ba49776 Mon Sep 17 00:00:00 2001 From: Michael McCarthy <51542091+mccarthy-m-g@users.noreply.github.com> Date: Thu, 22 Jan 2026 16:56:32 -0800 Subject: [PATCH 7/7] use `fit_options()` for `.fit_options` arg input --- NAMESPACE | 1 + R/run_eval.R | 31 ++++++++++++++++++++++++++----- man/fit_options.Rd | 20 ++++++++++++++++++++ man/run_eval.Rd | 8 +++----- man/run_eval_core.Rd | 6 ++---- 5 files changed, 52 insertions(+), 14 deletions(-) create mode 100644 man/fit_options.Rd diff --git a/NAMESPACE b/NAMESPACE index 916d326..683b32d 100644 --- a/NAMESPACE +++ b/NAMESPACE @@ -15,6 +15,7 @@ export(calculate_shrinkage) export(calculate_stats) export(compare_psn_execute_results) export(compare_psn_proseval_results) +export(fit_options) export(group_by_dose) export(group_by_time) export(install_default_literature_model) diff --git a/R/run_eval.R b/R/run_eval.R index c190463..c04c3cf 100644 --- a/R/run_eval.R +++ b/R/run_eval.R @@ -32,10 +32,8 @@ #' result from a call to [stats_summ_options()]. #' @param .vpc_options Options for VPC simulations. This must be the result from #' a call to [vpc_options()]. -#' @param .fit_options Options for controlling MAP Bayesian fit. E.g. -#' `.map_options = list(control = list(reltol = 1e-04))` will perform a slightly -#' coarser fit (but more stable) which can be useful in some case (default -#' `reltol = 1e-05`). +#' @param .fit_options Options for controlling MAP Bayesian fit. This must be +#' the result from a call to [fit_options()]. #' @param threads number of threads to divide computations on. Default is 1, #' i.e. no parallel execution #' @param ruv residual error variability magnitude, specified as list. @@ -67,7 +65,7 @@ run_eval <- function( incremental = FALSE, .stats_summ_options = stats_summ_options(), .vpc_options = vpc_options(), - .fit_options = list(), + .fit_options = fit_options(), threads = 1, progress = TRUE, verbose = TRUE @@ -195,3 +193,26 @@ run_eval <- function( ## 5. Return results out } + +#' Options for controlling MAP Bayesian fit +#' +#' @param ... These dots are reserved for future extensibility and must be empty. +#' @param reltol Relative convergence tolerance. `reltol = 1e-04` will perform a +#' slightly coarser but more stable fit, which can be useful in some case. +#' +#' @returns A list. +#' @export +fit_options <- function( + ..., + reltol = 1e-05 +) { + rlang::check_dots_empty() + out <- list( + control = list( + reltol = vctrs::vec_assert( + reltol, ptype = numeric(), size = 1L, arg = "reltol" + ) + ) + ) + structure(out, class = "mipdeval_fit_options") +} diff --git a/man/fit_options.Rd b/man/fit_options.Rd new file mode 100644 index 0000000..47f320c --- /dev/null +++ b/man/fit_options.Rd @@ -0,0 +1,20 @@ +% Generated by roxygen2: do not edit by hand +% Please edit documentation in R/run_eval.R +\name{fit_options} +\alias{fit_options} +\title{Options for controlling MAP Bayesian fit} +\usage{ +fit_options(..., reltol = 1e-05) +} +\arguments{ +\item{...}{These dots are reserved for future extensibility and must be empty.} + +\item{reltol}{Relative convergence tolerance. \code{reltol = 1e-04} will perform a +slightly coarser but more stable fit, which can be useful in some case.} +} +\value{ +A list. +} +\description{ +Options for controlling MAP Bayesian fit +} diff --git a/man/run_eval.Rd b/man/run_eval.Rd index 73a71a9..9c0c94c 100644 --- a/man/run_eval.Rd +++ b/man/run_eval.Rd @@ -21,7 +21,7 @@ run_eval( incremental = FALSE, .stats_summ_options = stats_summ_options(), .vpc_options = vpc_options(), - .fit_options = list(), + .fit_options = fit_options(), threads = 1, progress = TRUE, verbose = TRUE @@ -89,10 +89,8 @@ result from a call to \code{\link[=stats_summ_options]{stats_summ_options()}}.} \item{.vpc_options}{Options for VPC simulations. This must be the result from a call to \code{\link[=vpc_options]{vpc_options()}}.} -\item{.fit_options}{Options for controlling MAP Bayesian fit. E.g. -\code{.map_options = list(control = list(reltol = 1e-04))} will perform a slightly -coarser fit (but more stable) which can be useful in some case (default -\code{reltol = 1e-05}).} +\item{.fit_options}{Options for controlling MAP Bayesian fit. This must be +the result from a call to \code{\link[=fit_options]{fit_options()}}.} \item{threads}{number of threads to divide computations on. Default is 1, i.e. no parallel execution} diff --git a/man/run_eval_core.Rd b/man/run_eval_core.Rd index 2314866..a9c6020 100644 --- a/man/run_eval_core.Rd +++ b/man/run_eval_core.Rd @@ -48,10 +48,8 @@ approach has been called "model predictive control (MPC)" \item{progress_function}{function to increment progress bar} -\item{.fit_options}{Options for controlling MAP Bayesian fit. E.g. -\code{.map_options = list(control = list(reltol = 1e-04))} will perform a slightly -coarser fit (but more stable) which can be useful in some case (default -\code{reltol = 1e-05}).} +\item{.fit_options}{Options for controlling MAP Bayesian fit. This must be +the result from a call to \code{\link[=fit_options]{fit_options()}}.} } \value{ a \code{data.frame} with individual predictions