From ce75a932b06189c304b01f55f16bd8a5d219b307 Mon Sep 17 00:00:00 2001 From: Stu Field Date: Tue, 7 Jul 2026 17:36:36 -0600 Subject: [PATCH] Major bug fix to AUC of Summary method - major bugfix! - trapezoidal area calculation of AUC within the S3 summary() method was incorrect - resulted in clearly incorrect AUCs reported and added to plotting routine - fixed in this commit - updated affected unit tests --- R/utils.R | 16 +++---- tests/testthat/_snaps/stability-selection.md | 48 ++++++++++---------- 2 files changed, 31 insertions(+), 33 deletions(-) diff --git a/R/utils.R b/R/utils.R index e6fe82c..7b6df4f 100644 --- a/R/utils.R +++ b/R/utils.R @@ -54,25 +54,23 @@ get_cores <- function() { restore_rng_kind <- getFromNamespace("restore_rng_kind", "withr") #' Internal for calculating the AUC for all stability paths +#' using the trapezoidal rule. Must use `abs()` because decreasing. #' #' @param x A `stabpath_matrix` entry of a `stab_path` object. #' #' @param values A vector of values to iterate over, -#' typically `lambda_norm`. +#' typically `lambda_norm` (`lambda / max(lambda)`). #' #' @importFrom utils head tail #' @noRd calc_path_auc <- function(x, values) { stopifnot( - "Must pass `x` as a `stability path matrix`." = is_stabpath_matrix(x) - ) - sum_diffs <- head(values, -1L) + tail(values, -1L) - # flip matrix to perform row-wise operations column-wise - vapply( - data.frame(t(x)), - function(.x) as.double(diff(.x) %*% sum_diffs) / 2, - 0.1 + "`x` must be a `stability path matrix`." = is_stabpath_matrix(x) ) + auc_trapz <- function(x, y) { + sum(abs(diff(x)) * (head(y, -1L) + tail(y, -1L)) / 2) + } + apply(x, 1L, auc_trapz, x = values) } log_rfu <- function(x) { diff --git a/tests/testthat/_snaps/stability-selection.md b/tests/testthat/_snaps/stability-selection.md index 1da6e39..a53acb1 100644 --- a/tests/testthat/_snaps/stability-selection.md +++ b/tests/testthat/_snaps/stability-selection.md @@ -91,26 +91,26 @@ # A tibble: 20 x 4 feature MaxSelectProb AUC FDRbound - 1 feat_t 0.93 0.404939 NA - 2 feat_d 0.91 0.360777 NA - 3 feat_q 0.905 0.308470 NA - 4 feat_a 0.9 0.307907 NA - 5 feat_g 0.895 0.273354 NA - 6 feat_j 0.89 0.416035 NA - 7 feat_s 0.89 0.322614 NA - 8 feat_r 0.885 0.342585 NA - 9 feat_m 0.88 0.408485 NA - 10 feat_n 0.87 0.270069 NA - 11 feat_b 0.865 0.290262 NA - 12 feat_e 0.865 0.277968 NA - 13 feat_c 0.85 0.335890 NA - 14 feat_f 0.85 0.257492 NA - 15 feat_h 0.845 0.291972 NA - 16 feat_l 0.845 0.259685 NA - 17 feat_o 0.825 0.250749 NA - 18 feat_p 0.815 0.223145 NA - 19 feat_k 0.81 0.257189 NA - 20 feat_i 0.8 0.270370 NA + 1 feat_t 0.93 0.436532 NA + 2 feat_d 0.91 0.408626 NA + 3 feat_q 0.905 0.256633 NA + 4 feat_a 0.9 0.286384 NA + 5 feat_g 0.895 0.222145 NA + 6 feat_j 0.89 0.440141 NA + 7 feat_s 0.89 0.291720 NA + 8 feat_r 0.885 0.342005 NA + 9 feat_m 0.88 0.533218 NA + 10 feat_n 0.87 0.220431 NA + 11 feat_b 0.865 0.245937 NA + 12 feat_e 0.865 0.253644 NA + 13 feat_c 0.85 0.287507 NA + 14 feat_f 0.85 0.224109 NA + 15 feat_h 0.845 0.243903 NA + 16 feat_l 0.845 0.221616 NA + 17 feat_o 0.825 0.198936 NA + 18 feat_p 0.815 0.186961 NA + 19 feat_k 0.81 0.216318 NA + 20 feat_i 0.8 0.235128 NA --- @@ -120,10 +120,10 @@ # A tibble: 4 x 4 feature MaxSelectProb AUC FDRbound - 1 feat_t 0.93 0.404939 0.003125 - 2 feat_d 0.91 0.360777 0.00625 - 3 feat_q 0.905 0.308470 0.009375 - 4 feat_a 0.9 0.307907 0.0125 + 1 feat_t 0.93 0.436532 0.003125 + 2 feat_d 0.91 0.408626 0.00625 + 3 feat_q 0.905 0.256633 0.009375 + 4 feat_a 0.9 0.286384 0.0125 # `stab_sel` S3 print method returns expected known output