Skip to content

Commit 834ef50

Browse files
authored
Merge pull request #31 from InsightRX/stat-accuracy
Add accuracy metric to results stats summary
2 parents 1f51072 + 070a7bb commit 834ef50

14 files changed

Lines changed: 371 additions & 11 deletions

NAMESPACE

Lines changed: 7 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -8,6 +8,7 @@ S3method(print,mipdeval_results)
88
S3method(print,mipdeval_results_bayesian_impact)
99
S3method(print,mipdeval_results_shrinkage)
1010
S3method(print,mipdeval_results_stats_summ)
11+
export(accuracy)
1112
export(add_grouping_column)
1213
export(calculate_bayesian_impact)
1314
export(calculate_shrinkage)
@@ -17,13 +18,19 @@ export(compare_psn_proseval_results)
1718
export(group_by_dose)
1819
export(group_by_time)
1920
export(install_default_literature_model)
21+
export(is_accurate)
22+
export(is_accurate_abs)
23+
export(is_accurate_rel)
2024
export(new_ode_model)
2125
export(parse_psn_proseval_results)
2226
export(reldiff_psn_execute_results)
2327
export(reldiff_psn_proseval_results)
2428
export(run_eval)
29+
export(stats_summ_options)
2530
export(vpc_options)
2631
importFrom(PKPDsim,install_default_literature_model)
2732
importFrom(PKPDsim,new_ode_model)
2833
importFrom(rlang,.data)
34+
importFrom(rlang,caller_arg)
35+
importFrom(rlang,caller_env)
2936
importFrom(utils,read.csv)

NEWS.md

Lines changed: 2 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -1,5 +1,7 @@
11
# mipdeval (development version)
22

3+
* Added `accuracy()` statistic for calculating the proportion of predicted drug concentrations that fall within specified absolute and relative error margins. This statistic will only be calculated in `run_eval()` when the absolute and relative error margins have been specified in the `.stats_summ_options` argument (see `stats_summ_options()`). (#30)
4+
35
# mipdeval 0.1.0
46

57
* Initial version.

R/calculate_stats.R

Lines changed: 65 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -3,13 +3,19 @@
33
#' @param res output object (`mipdeval_results`) from `run_eval()`, or
44
#' `data.frame` with raw results.
55
#' @param rounding number of decimals to round to.
6+
#' @param acc_error_abs,acc_error_rel For calculating [accuracy()]: Positive number
7+
#' providing an absolute or relative error margin. The cutoff is exclusive of
8+
#' the error margin. When `NULL` (the default), accuracy will not be
9+
#' calculated and will return `NA` instead.
610
#'
711
#' @returns tibble
812
#'
913
#' @export
1014
calculate_stats <- function(
1115
res,
12-
rounding = 3
16+
rounding = 3,
17+
acc_error_abs = NULL,
18+
acc_error_rel = NULL
1319
) {
1420
if(inherits(res, "mipdeval_results")) {
1521
res <- res$results
@@ -23,10 +29,67 @@ calculate_stats <- function(
2329
rmse = rmse(.data$dv, .data$value),
2430
nrmse = nrmse(.data$dv, .data$value),
2531
mpe = mpe(.data$dv, .data$value),
26-
mape = mape(.data$dv, .data$value)
32+
mape = mape(.data$dv, .data$value),
33+
accuracy = dplyr::if_else(
34+
!is.null(acc_error_abs) & !is.null(acc_error_rel),
35+
accuracy(.data$dv, .data$value, acc_error_abs, acc_error_rel),
36+
NA
37+
)
2738
) |>
2839
dplyr::mutate(dplyr::across(dplyr::everything(), round, rounding)) |>
2940
dplyr::as_tibble()
3041
class(out) <- c("mipdeval_results_stats_summ", class(out))
3142
out
3243
}
44+
45+
#' Options for summary statistics
46+
#'
47+
#' @param ... These dots are reserved for future extensibility and must be empty.
48+
#' @inheritParams calculate_stats
49+
#'
50+
#' @returns A list.
51+
#' @export
52+
stats_summ_options <- function(
53+
...,
54+
rounding = 3,
55+
acc_error_abs = NULL,
56+
acc_error_rel = NULL
57+
) {
58+
rlang::check_dots_empty()
59+
check_required_accuracy(acc_error_abs, acc_error_rel)
60+
out <- list(
61+
rounding = vctrs::vec_assert(
62+
rounding, ptype = numeric(), size = 1L, arg = "rounding"
63+
),
64+
acc_error_abs = vec_assert_or_null(
65+
acc_error_abs, ptype = numeric(), size = 1L, arg = "acc_error_abs"
66+
),
67+
acc_error_rel = vec_assert_or_null(
68+
acc_error_rel, ptype = numeric(), size = 1L, arg = "acc_error_rel"
69+
)
70+
)
71+
structure(out, class = "mipdeval_stats_summ_options")
72+
}
73+
74+
check_required_accuracy <- function(
75+
acc_error_abs = NULL,
76+
acc_error_rel = NULL,
77+
call = caller_env()
78+
) {
79+
acc_error_abs_null <- is.null(acc_error_abs) & !is.null(acc_error_rel)
80+
acc_error_rel_null <- !is.null(acc_error_abs) & is.null(acc_error_rel)
81+
if (!acc_error_abs_null & !acc_error_rel_null) return()
82+
cli::cli_abort(
83+
message = c(
84+
paste0(
85+
"{.arg acc_error_abs} and {.arg acc_error_rel} must both be specified ",
86+
"to calculate accuracy."
87+
),
88+
"x" = dplyr::case_when(
89+
acc_error_abs_null ~ "{.arg acc_error_abs} is NULL.",
90+
acc_error_rel_null ~ "{.arg acc_error_rel} is NULL."
91+
)
92+
),
93+
call = call
94+
)
95+
}

R/mipdeval-package.R

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -2,7 +2,7 @@
22
"_PACKAGE"
33

44
## usethis namespace: start
5-
#' @importFrom rlang .data
5+
#' @importFrom rlang .data caller_arg caller_env
66
#' @importFrom utils read.csv
77
## usethis namespace: end
88
NULL

R/misc.R

Lines changed: 77 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -26,6 +26,25 @@ is_timevarying <- function(.data, .cols) {
2626
)
2727
}
2828

29+
#' Assert an argument has known prototype and/or size or is NULL
30+
#'
31+
#' @inheritParams vctrs::vec_assert
32+
#'
33+
#' @returns Either throws an error or returns `x`, invisibly.
34+
vec_assert_or_null <- function(
35+
x,
36+
ptype = NULL,
37+
size = NULL,
38+
arg = caller_arg(x),
39+
call = caller_env()
40+
) {
41+
if (!is.null(x)) {
42+
vctrs::vec_assert(x = x, ptype = ptype, size = size, arg = arg, call = call)
43+
} else {
44+
NULL
45+
}
46+
}
47+
2948
#' Root-mean-squared error
3049
#'
3150
#' @param obs observations vector
@@ -68,6 +87,63 @@ mpe <- function (obs, pred) {
6887
sum((obs - pred)/obs)/length(obs)
6988
}
7089

90+
#' Accuracy
91+
#'
92+
#' Accuracy provides a measure of clinical suitability, defined by whether model
93+
#' predicted drug concentrations fall within an absolute OR relative error
94+
#' margin of the measured concentrations.
95+
#'
96+
#' @param obs Observations vector.
97+
#' @param pred Predictions vector.
98+
#' @param error_abs,error_rel Positive number providing an absolute or relative
99+
#' error margin. The cutoff is exclusive of the error margin. Defaults to `0`,
100+
#' meaning no predictions fall within the error margin.
101+
#'
102+
#' @returns For `is_accurate()`, `is_accurate_abs()`, and `is_accurate_rel()`: A
103+
#' logical vector indicating whether or not each predicted drug concentration
104+
#' was considered accurate according to the specified absolute or relative
105+
#' error margin(s).
106+
#'
107+
#' For `accuracy()`: A single value between 0 and 1 indicating the proportion
108+
#' of predicted drug concentrations that fell within the specified absolute
109+
#' and relative error margins.
110+
#'
111+
#' @examples
112+
#' # Does the predicted drug concentration fall within 0.5 mg/L error margin?
113+
#' is_accurate_abs(6, 5, error_abs = 0.5)
114+
#'
115+
#' # Does the predicted drug concentration fall within 25% error margin?
116+
#' is_accurate_rel(6, 5, error_rel = 0.25)
117+
#'
118+
#' # Does the predicted drug concentration fall within 0.5 mg/L OR 25%?
119+
#' is_accurate(6, 5, error_abs = 0.5, error_rel = 0.25)
120+
#'
121+
#' # What proportion of predicted drug concentrations fell within 0.5 mg/L OR 25%?
122+
#' accuracy(rnorm(10, 6), rnorm(10, 5), error_abs = 0.5, error_rel = 0.25)
123+
#'
124+
#' @export
125+
accuracy <- function(obs, pred, error_abs = 0, error_rel = 0) {
126+
mean(is_accurate(obs, pred, error_abs, error_rel))
127+
}
128+
129+
#' @rdname accuracy
130+
#' @export
131+
is_accurate <- function(obs, pred, error_abs = 0, error_rel = 0) {
132+
is_accurate_abs(obs, pred, error_abs) | is_accurate_rel(obs, pred, error_rel)
133+
}
134+
135+
#' @rdname accuracy
136+
#' @export
137+
is_accurate_abs <- function(obs, pred, error_abs = 0) {
138+
abs(pred - obs) < error_abs
139+
}
140+
141+
#' @rdname accuracy
142+
#' @export
143+
is_accurate_rel <- function(obs, pred, error_rel = 0) {
144+
(pred/obs > 1 - error_rel) & (pred/obs < 1 + error_rel)
145+
}
146+
71147
#' Weighted sum-of-squares of residuals
72148
#'
73149
#' @inheritParams rmse
@@ -83,3 +159,4 @@ ss <- function(obs, pred, w = NULL) {
83159
if(sum(w) == 0) return(NA)
84160
sum(w * (obs - pred)^2)
85161
}
162+

R/run_eval.R

Lines changed: 12 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -28,6 +28,8 @@
2828
#' approach has been called "model predictive control (MPC)"
2929
#' (www.page-meeting.org/?abstract=9076) and may be more predictive than
3030
#' "regular" MAP in some scenarios. Default is `FALSE`.
31+
#' @param .stats_summ_options Options for summary statistics. This must be the
32+
#' result from a call to [stats_summ_options()].
3133
#' @param .vpc_options Options for VPC simulations. This must be the result from
3234
#' a call to [vpc_options()].
3335
#' @param threads number of threads to divide computations on. Default is 1,
@@ -59,11 +61,15 @@ run_eval <- function(
5961
weight_prior = 1,
6062
censor_covariates = TRUE,
6163
incremental = FALSE,
64+
.stats_summ_options = stats_summ_options(),
6265
.vpc_options = vpc_options(),
6366
threads = 1,
6467
progress = TRUE,
6568
verbose = TRUE
6669
) {
70+
# Evaluate options at the start so any errors are triggered immediately.
71+
.stats_summ_options <- .stats_summ_options
72+
.vpc_options <- .vpc_options
6773

6874
if(progress) { # configure progressbars
6975
progressr::handlers(global = TRUE)
@@ -169,7 +175,12 @@ run_eval <- function(
169175
# res is NULL when vpc_options(..., vpc_only = TRUE).
170176
if (!is.null(res)) {
171177
if(verbose) cli::cli_progress_step("Calculating forecasting statistics")
172-
out$stats_summ <- calculate_stats(out)
178+
out$stats_summ <- calculate_stats(
179+
out,
180+
rounding = .stats_summ_options$rounding,
181+
acc_error_abs = .stats_summ_options$acc_error_abs,
182+
acc_error_rel = .stats_summ_options$acc_error_rel
183+
)
173184
out$shrinkage <- calculate_shrinkage(out)
174185
out$bayesian_impact <- calculate_bayesian_impact(out)
175186
cli::cli_progress_done()

man/accuracy.Rd

Lines changed: 55 additions & 0 deletions
Some generated files are not rendered by default. Learn more about customizing how changed files appear on GitHub.

man/calculate_stats.Rd

Lines changed: 6 additions & 1 deletion
Some generated files are not rendered by default. Learn more about customizing how changed files appear on GitHub.

man/run_eval.Rd

Lines changed: 4 additions & 0 deletions
Some generated files are not rendered by default. Learn more about customizing how changed files appear on GitHub.

man/stats_summ_options.Rd

Lines changed: 29 additions & 0 deletions
Some generated files are not rendered by default. Learn more about customizing how changed files appear on GitHub.

0 commit comments

Comments
 (0)