From cf41e28ee07b5bd2aae7cae21adf87cecdf9d2fb Mon Sep 17 00:00:00 2001 From: fweber144 Date: Mon, 5 May 2025 12:25:26 +0200 Subject: [PATCH 1/4] fix `sqrt()` of negative value in MSE calculation --- R/summary_funs.R | 10 +++++++++- 1 file changed, 9 insertions(+), 1 deletion(-) diff --git a/R/summary_funs.R b/R/summary_funs.R index 4655d79f1..baa776ab0 100644 --- a/R/summary_funs.R +++ b/R/summary_funs.R @@ -377,7 +377,15 @@ get_stat <- function(summaries, summaries_baseline = NULL, ((mu_baseline - y)^2 - mse_b)) / (n_full - 1) } if (stat != "rmse") { - value_se <- sqrt(value_se^2 - 2 * cov_mse_e_b + var_mse_b) + value_se_sq <- value_se^2 - 2 * cov_mse_e_b + var_mse_b + if (!is.na(value_se_sq) && sign(value_se_sq) == -1) { + if (abs(value_se_sq) < sqrt(.Machine$double.eps)) { + value_se_sq <- 0 + } else { + stop("Negative (and numerically non-zero) `value_se_sq`.") + } + } + value_se <- sqrt(value_se_sq) } } if (stat == "mse") { From f98216cf8ce1244fb4454af248d7cca04aaf2ea7 Mon Sep 17 00:00:00 2001 From: fweber144 Date: Mon, 5 May 2025 12:28:47 +0200 Subject: [PATCH 2/4] fix possible `sqrt()` of negative value in RMSE calculation --- R/summary_funs.R | 16 ++++++++++++---- 1 file changed, 12 insertions(+), 4 deletions(-) diff --git a/R/summary_funs.R b/R/summary_funs.R index baa776ab0..9c41629db 100644 --- a/R/summary_funs.R +++ b/R/summary_funs.R @@ -395,12 +395,20 @@ get_stat <- function(summaries, summaries_baseline = NULL, value <- sqrt(mse_e) - ifelse(is.null(summaries_baseline), 0, sqrt(mse_b)) # the first-order Taylor approximation of the variance if (is.null(summaries_baseline)) { - value_se <- sqrt(value_se^2 / mse_e / 4) + value_se_sq <- value_se^2 / mse_e / 4 } else { - value_se <- sqrt((value_se^2 / mse_e - - 2 * cov_mse_e_b / sqrt(mse_e * mse_b) + - var_mse_b / mse_b) / 4) + value_se_sq <- (value_se^2 / mse_e - + 2 * cov_mse_e_b / sqrt(mse_e * mse_b) + + var_mse_b / mse_b) / 4 } + if (!is.na(value_se_sq) && sign(value_se_sq) == -1) { + if (abs(value_se_sq) < sqrt(.Machine$double.eps)) { + value_se_sq <- 0 + } else { + stop("Negative (and numerically non-zero) `value_se_sq`.") + } + } + value_se <- sqrt(value_se_sq) } else if (stat == "R2") { y_mean_w <- mean(wobs * y) # simple transformation of mse From 1ce889ad130ea5990ab7fa749f891fdfb8614dd6 Mon Sep 17 00:00:00 2001 From: fweber144 Date: Mon, 5 May 2025 12:34:26 +0200 Subject: [PATCH 3/4] add function `sqrt_cut0()` to avoid redundancies in the code --- R/misc.R | 12 ++++++++++++ R/summary_funs.R | 40 ++++++++-------------------------------- 2 files changed, 20 insertions(+), 32 deletions(-) diff --git a/R/misc.R b/R/misc.R index 0a45af5a0..6712f6bee 100644 --- a/R/misc.R +++ b/R/misc.R @@ -743,3 +743,15 @@ use_progressr <- function() { interactive() && identical(foreach::getDoParName(), "doFuture")) } + +sqrt_cut0 <- function(x) { + if (!is.na(x) && sign(x) == -1) { + if (abs(x) < sqrt(.Machine$double.eps)) { + x <- 0 + } else { + stop("Negative (and numerically non-zero) value used as input to ", + "sqrt_cut0().") + } + } + return(sqrt(x)) +} diff --git a/R/summary_funs.R b/R/summary_funs.R index 9c41629db..938d5207b 100644 --- a/R/summary_funs.R +++ b/R/summary_funs.R @@ -377,15 +377,7 @@ get_stat <- function(summaries, summaries_baseline = NULL, ((mu_baseline - y)^2 - mse_b)) / (n_full - 1) } if (stat != "rmse") { - value_se_sq <- value_se^2 - 2 * cov_mse_e_b + var_mse_b - if (!is.na(value_se_sq) && sign(value_se_sq) == -1) { - if (abs(value_se_sq) < sqrt(.Machine$double.eps)) { - value_se_sq <- 0 - } else { - stop("Negative (and numerically non-zero) `value_se_sq`.") - } - } - value_se <- sqrt(value_se_sq) + value_se <- sqrt_cut0(value_se^2 - 2 * cov_mse_e_b + var_mse_b) } } if (stat == "mse") { @@ -395,20 +387,12 @@ get_stat <- function(summaries, summaries_baseline = NULL, value <- sqrt(mse_e) - ifelse(is.null(summaries_baseline), 0, sqrt(mse_b)) # the first-order Taylor approximation of the variance if (is.null(summaries_baseline)) { - value_se_sq <- value_se^2 / mse_e / 4 + value_se <- sqrt_cut0(value_se^2 / mse_e / 4) } else { - value_se_sq <- (value_se^2 / mse_e - - 2 * cov_mse_e_b / sqrt(mse_e * mse_b) + - var_mse_b / mse_b) / 4 + value_se <- sqrt_cut0((value_se^2 / mse_e - + 2 * cov_mse_e_b / sqrt(mse_e * mse_b) + + var_mse_b / mse_b) / 4) } - if (!is.na(value_se_sq) && sign(value_se_sq) == -1) { - if (abs(value_se_sq) < sqrt(.Machine$double.eps)) { - value_se_sq <- 0 - } else { - stop("Negative (and numerically non-zero) `value_se_sq`.") - } - } - value_se <- sqrt(value_se_sq) } else if (stat == "R2") { y_mean_w <- mean(wobs * y) # simple transformation of mse @@ -454,17 +438,9 @@ get_stat <- function(summaries, summaries_baseline = NULL, # delta=TRUE mse_e <- mse_e - mse_b } - value_se_sq <- (value_se^2 - - 2 * mse_e / mse_y * cov_mse_e_y + - (mse_e / mse_y)^2 * var_mse_y) / mse_y^2 - if (!is.na(value_se_sq) && sign(value_se_sq) == -1) { - if (abs(value_se_sq) < sqrt(.Machine$double.eps)) { - value_se_sq <- 0 - } else { - stop("Negative (and numerically non-zero) `value_se_sq`.") - } - } - value_se <- sqrt(value_se_sq) + value_se <- sqrt_cut0((value_se^2 - + 2 * mse_e / mse_y * cov_mse_e_y + + (mse_e / mse_y)^2 * var_mse_y) / mse_y^2) } } else if (stat %in% c("acc", "pctcorr", "auc")) { y <- y_wobs_test$y From 8ca624cc2a18e85555fa84bb297e2b9e3687ea84 Mon Sep 17 00:00:00 2001 From: fweber144 Date: Mon, 5 May 2025 21:25:00 +0200 Subject: [PATCH 4/4] use `sqrt()` instead of `sqrt_cut0()` where possible --- R/summary_funs.R | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/R/summary_funs.R b/R/summary_funs.R index 938d5207b..e6f45ec36 100644 --- a/R/summary_funs.R +++ b/R/summary_funs.R @@ -387,7 +387,7 @@ get_stat <- function(summaries, summaries_baseline = NULL, value <- sqrt(mse_e) - ifelse(is.null(summaries_baseline), 0, sqrt(mse_b)) # the first-order Taylor approximation of the variance if (is.null(summaries_baseline)) { - value_se <- sqrt_cut0(value_se^2 / mse_e / 4) + value_se <- sqrt(value_se^2 / mse_e / 4) } else { value_se <- sqrt_cut0((value_se^2 / mse_e - 2 * cov_mse_e_b / sqrt(mse_e * mse_b) +