From e33e8489539d3edaaa41a04e70096b46c6c0dd9c Mon Sep 17 00:00:00 2001 From: Aaron Caldwell Date: Fri, 17 Apr 2026 07:19:30 -0500 Subject: [PATCH 1/2] fix for potential errors --- .claude/settings.local.json | 14 ++++++- R/cor_test.R | 21 ++++++++-- tests/testthat/test-corr.R | 76 +++++++++++++++++++++++++++++++++++++ 3 files changed, 107 insertions(+), 4 deletions(-) diff --git a/.claude/settings.local.json b/.claude/settings.local.json index c7c82df..94bc906 100644 --- a/.claude/settings.local.json +++ b/.claude/settings.local.json @@ -128,7 +128,19 @@ "Bash(ls:*)", "Bash(\"C:/Program Files/R/R-4.5.1/bin/Rscript.exe\" -e \"git diff HEAD -- tests/testthat/test-smd_calc_htest.R\")", "Bash(\"C:/Program Files/R/R-4.5.1/bin/Rscript.exe\" junk/diag_bca.R)", - "Bash(git -C /Users/aaroncaldwell/Documents/GitHub/TOSTER check-ignore vignettes/robustTOST.html)" + "Bash(git -C /Users/aaroncaldwell/Documents/GitHub/TOSTER check-ignore vignettes/robustTOST.html)", + "Bash(\"/c/Program Files/R/R-4.5.3/bin/Rscript.exe\" -e \"devtools::load_all\\(\\); testthat::test_file\\('tests/testthat/test-corr.R', reporter = 'summary'\\)\")", + "Bash('/c/Program Files/R/R-4.5.3/bin/Rscript.exe' -e ':*)", + "Bash(\"/c/Program Files/R/R-4.5.3/bin/Rscript.exe\" -e \"cat\\('hello\\\\n'\\)\")", + "Bash(\"/c/Program Files/R/R-4.5.3/bin/Rscript.exe\" -e \"devtools::load_all\\(\\); cat\\('loaded\\\\n'\\)\")", + "Bash(\"/c/Program Files/R/R-4.5.3/bin/Rscript.exe\" -e \"devtools::load_all\\(\\); x <- c\\(1,2,3,NA,5\\); y <- c\\(2,4,NA,8,10\\); cat\\(cor\\(x,y,use='complete.obs'\\),'\\\\n'\\)\")", + "Bash(\"/c/Program Files/R/R-4.5.3/bin/Rscript.exe\" -e \"devtools::load_all\\(\\); x <- c\\(1,2,3,NA,5,6,7,8,9,10\\); y <- c\\(2,4,NA,8,10,12,14,16,18,20\\); res <- z_cor_test\\(x,y,method='pearson'\\); print\\(res\\)\")", + "Bash(\"/c/Program Files/R/R-4.5.3/bin/Rscript.exe\" -e \"devtools::load_all\\(\\); testthat::test_file\\('tests/testthat/test-corr.R', reporter = 'progress'\\)\")", + "Bash('/c/Program Files/R/R-4.5.3/bin/Rscript.exe' -e ' *)", + "Bash(\"/c/Program Files/R/R-4.5.3/bin/Rscript.exe\" -e \"cat\\('atanh\\(1\\) =', atanh\\(1\\), '\\\\n'\\); cat\\('atanh\\(-1\\) =', atanh\\(-1\\), '\\\\n'\\)\")", + "Bash(\"/c/Program Files/R/R-4.5.3/bin/Rscript.exe\" -e \"x <- 1:10; y <- 2*\\(1:10\\); cat\\('cor=', cor\\(x,y\\), '\\\\n'\\); jk <- sapply\\(1:10, function\\(i\\) atanh\\(cor\\(x[-i], y[-i]\\)\\)\\); cat\\('all Inf:', all\\(is.infinite\\(jk\\)\\), '\\\\n'\\); cat\\('mean-jk[1]:', mean\\(jk\\)-jk[1], '\\\\n'\\)\")", + "Bash(\"/c/Program Files/R/R-4.5.3/bin/Rscript.exe\" -e \"x <- 1:10; y <- 2*\\(1:10\\); cat\\('pearson:', cor\\(x,y,method='pearson'\\), '\\\\n'\\); cat\\('spearman:', cor\\(x,y,method='spearman'\\), '\\\\n'\\); cat\\('kendall:', cor\\(x,y,method='kendall'\\), '\\\\n'\\); cat\\('identical to 1:', cor\\(x,y,method='spearman'\\)==1, cor\\(x,y,method='kendall'\\)==1, '\\\\n'\\)\")", + "Bash('/c/Program Files/R/R-4.5.3/bin/Rscript.exe' -e 'devtools::load_all\\(\\); x <- 1:10; y <- 2*\\(1:10\\); res <- z_cor_test\\(x, y, method='\\\\''spearman'\\\\'', se_method='\\\\''jackknife'\\\\''\\); cat\\('\\\\''z.se='\\\\'', res$stderr['\\\\''z.se'\\\\''], '\\\\''\\\\nestimate='\\\\'', res$estimate, '\\\\''\\\\np='\\\\'', res$p.value, '\\\\''\\\\n'\\\\''\\)')" ] } } diff --git a/R/cor_test.R b/R/cor_test.R index 854f8bb..a60b896 100644 --- a/R/cor_test.R +++ b/R/cor_test.R @@ -148,15 +148,17 @@ z_cor_test = function(x, intmult = c(NA,1) } } - r_xy = cor(x,y, - method = method) + DNAME <- paste(deparse(substitute(x)), "and", deparse(substitute(y))) df = data.frame(x=x, y=y) df = na.omit(df) + x = df$x + y = df$y n_obs = nrow(df) + r_xy = cor(x,y, + method = method) z_xy = rho_to_z(r_xy) - DNAME <- paste(deparse(substitute(x)), "and", deparse(substitute(y))) NVAL = null znull = rho_to_z(null) # capture raw method string before overwrite -------- @@ -189,10 +191,23 @@ z_cor_test = function(x, # jackknife SE -------- if (se_method == "jackknife") { + if (!is.finite(r_xy) || abs(r_xy) >= 1) { + stop("Jackknife SE cannot be computed when |r| = 1 ", + "(Fisher z-transform is infinite). ", + "Use se_method = 'analytic' instead.", + call. = FALSE) + } jk_cors <- vapply(seq_len(n_obs), function(i) { rho_to_z(cor(x[-i], y[-i], method = method_raw)) }, numeric(1)) z.se <- sqrt(((n_obs - 1) / n_obs) * sum((jk_cors - mean(jk_cors))^2)) + if (!is.finite(z.se) || z.se <= 0) { + stop("Jackknife SE cannot be computed: result is non-finite or zero. ", + "This commonly occurs when |r| is numerically at or near 1 ", + "(e.g. strictly monotonic data under Spearman/Kendall). ", + "Use se_method = 'analytic' instead.", + call. = FALSE) + } } # append SE label to method name -------- diff --git a/tests/testthat/test-corr.R b/tests/testthat/test-corr.R index cf9a376..906072c 100644 --- a/tests/testthat/test-corr.R +++ b/tests/testthat/test-corr.R @@ -914,3 +914,79 @@ test_that("boot_cor_test: stud results similar to z_cor_test for large n Pearson # p-values should be in the same ballpark expect_equal(boot_res$p.value, z_res$p.value, tolerance = 0.1) }) + +# z_cor_test missing data handling ----- + +test_that("z_cor_test handles missing data correctly", { + + set.seed(54321) + x <- c(1.1, 2.3, 3.0, NA, 5.2, 6.1, 7.8, 8.0, 9.4, 10.7) + y <- c(2.2, 3.8, NA, 7.9, 10.5, 11.6, 14.3, 15.1, 17.9, 21.0) + + # Should not error with NAs present + for (m in c("pearson", "spearman", "kendall")) { + res <- z_cor_test(x, y, method = m) + expect_true(is.finite(res$p.value), + label = paste("finite p with NAs for", m)) + expect_true(is.finite(unname(res$estimate)), + label = paste("finite estimate with NAs for", m)) + # N should reflect complete cases only (8 of 10) + expect_equal(unname(res$parameter), 8, + label = paste("N reflects complete cases for", m)) + } + + # Results should match manually removing NAs + complete <- complete.cases(x, y) + x_clean <- x[complete] + y_clean <- y[complete] + + res_na <- z_cor_test(x, y, method = "pearson") + res_clean <- z_cor_test(x_clean, y_clean, method = "pearson") + + expect_equal(res_na$p.value, res_clean$p.value) + expect_equal(unname(res_na$estimate), unname(res_clean$estimate)) + expect_equal(unname(res_na$statistic), unname(res_clean$statistic)) + + # Jackknife SE should also work with NAs + res_jk <- z_cor_test(x, y, method = "pearson", se_method = "jackknife") + expect_true(is.finite(res_jk$p.value), + label = "jackknife with NAs gives finite p") + expect_equal(unname(res_jk$parameter), 8) + + # Equivalence test with NAs + res_eq <- z_cor_test(x, y, method = "pearson", + alternative = "equivalence", null = 0.99) + expect_true(is.finite(res_eq$p.value), + label = "equivalence with NAs gives finite p") +}) + +test_that("z_cor_test errors when |r| = 1 with jackknife SE", { + # Perfect positive correlation: atanh(1) = Inf -> jackknife SE becomes NaN + x <- 1:10 + y <- 2 * (1:10) + + # Analytic path still works (returns Inf z / p ~ 0) + expect_silent(z_cor_test(x, y, method = "pearson")) + + # Jackknife should error with a clear message + expect_error( + z_cor_test(x, y, method = "pearson", se_method = "jackknife"), + regexp = "Jackknife SE cannot be computed" + ) + + # Perfect negative correlation also triggers the guard + expect_error( + z_cor_test(x, -y, method = "pearson", se_method = "jackknife"), + regexp = "Jackknife SE cannot be computed" + ) + + # Spearman/Kendall on strictly monotonic data also yield |r| = 1 + expect_error( + z_cor_test(x, y, method = "spearman", se_method = "jackknife"), + regexp = "Jackknife SE cannot be computed" + ) + expect_error( + z_cor_test(x, y, method = "kendall", se_method = "jackknife"), + regexp = "Jackknife SE cannot be computed" + ) +}) From 0ad55df9dfa72e1241805919bfee8681afeddf13 Mon Sep 17 00:00:00 2001 From: Aaron Caldwell Date: Fri, 17 Apr 2026 10:37:25 -0500 Subject: [PATCH 2/2] fix min n requirement --- R/cor_test.R | 11 ++++++++++ tests/testthat/test-corr.R | 42 ++++++++++++++++++++++++++++++++++++++ 2 files changed, 53 insertions(+) diff --git a/R/cor_test.R b/R/cor_test.R index a60b896..4d13601 100644 --- a/R/cor_test.R +++ b/R/cor_test.R @@ -155,6 +155,17 @@ z_cor_test = function(x, x = df$x y = df$y n_obs = nrow(df) + # minimum sample size for SE formulas -------- + # Pearson/Spearman SE use (n - 3); Kendall uses (n - 4). + min_n <- if (method == "kendall") 5 else 4 + if (n_obs < min_n) { + stop("Not enough complete cases to compute the ", + method, " correlation standard error: ", + "found ", n_obs, " complete observation(s), ", + "need at least ", min_n, ". ", + "Check for missing values in x/y.", + call. = FALSE) + } r_xy = cor(x,y, method = method) diff --git a/tests/testthat/test-corr.R b/tests/testthat/test-corr.R index 906072c..7fad4dd 100644 --- a/tests/testthat/test-corr.R +++ b/tests/testthat/test-corr.R @@ -990,3 +990,45 @@ test_that("z_cor_test errors when |r| = 1 with jackknife SE", { regexp = "Jackknife SE cannot be computed" ) }) + +test_that("z_cor_test errors when too few complete cases for SE formula", { + # Pearson/Spearman require n >= 4; Kendall requires n >= 5. + + # --- Raw sample too small --- + x3 <- c(1.1, 2.3, 4.7) + y3 <- c(0.5, 1.9, 3.2) + expect_error( + z_cor_test(x3, y3, method = "pearson"), + regexp = "Not enough complete cases" + ) + expect_error( + z_cor_test(x3, y3, method = "spearman"), + regexp = "Not enough complete cases" + ) + + # n = 4 is below Kendall's threshold but OK for Pearson/Spearman + x4 <- c(1.1, 2.3, 4.7, 6.8) + y4 <- c(0.5, 1.9, 3.2, 5.1) + expect_error( + z_cor_test(x4, y4, method = "kendall"), + regexp = "Not enough complete cases" + ) + expect_silent(z_cor_test(x4, y4, method = "pearson")) + expect_silent(z_cor_test(x4, y4, method = "spearman")) + + # --- NAs reducing complete cases below threshold --- + # 10 obs but only 3 complete -> fails Pearson/Spearman + xn <- c(1.1, 2.3, 4.7, NA, NA, NA, NA, NA, NA, NA) + yn <- c(0.5, 1.9, 3.2, NA, NA, NA, NA, NA, NA, NA) + expect_error( + z_cor_test(xn, yn, method = "pearson"), + regexp = "Not enough complete cases" + ) + + # Error message mentions method and observed N + err <- tryCatch(z_cor_test(x3, y3, method = "kendall"), + error = function(e) conditionMessage(e)) + expect_match(err, "kendall") + expect_match(err, "3 complete observation") + expect_match(err, "at least 5") +})