jared
6051e4bb5c
R-CMD-check / R CMD check (push) Successful in 4m37s
- Remove Knowles_Headshot_2019_good.jpg and Knowles_Headshot_2019_prisma.jpg from inst/img/ — personal headshot should not be distributed with the package - Update plot_jpeg() roxygen example to use a generic placeholder instead of the removed headshot path; regenerate man/plot_jpeg.Rd accordingly - Add comprehensive testthat coverage for clopper_pearson(), z_univariate(), waldInterval(), and agresti_coull_interval() in tests/testthat/test_propint.R The propint module previously had zero tests. New tests cover: return types, formula correctness (cross-validated against binom.test()), edge cases (0/n and n/n), interval validity, confidence level behavior, and sign/direction of z-scores.
174 lines
6.1 KiB
R
174 lines
6.1 KiB
R
# Tests for proportion confidence interval functions
|
|
|
|
# --- clopper_pearson ---------------------------------------------------------
|
|
|
|
test_that("clopper_pearson returns a named numeric vector of length 3", {
|
|
result <- clopper_pearson(20, 40)
|
|
expect_type(result, "double")
|
|
expect_length(result, 3)
|
|
expect_named(result, c("low", "observed", "high"))
|
|
})
|
|
|
|
test_that("clopper_pearson observed matches num/den", {
|
|
result <- clopper_pearson(20, 40)
|
|
expect_equal(unname(result["observed"]), 20 / 40)
|
|
})
|
|
|
|
test_that("clopper_pearson matches binom.test in base R", {
|
|
# The implementation is documented as producing the same results as
|
|
# binom.test(), which uses the Clopper-Pearson method.
|
|
for (num in c(0, 1, 5, 20, 39, 40)) {
|
|
result <- clopper_pearson(num, 40)
|
|
bt <- binom.test(num, 40)$conf.int
|
|
expect_equal(unname(result["low"]), bt[1], tolerance = 1e-10)
|
|
expect_equal(unname(result["high"]), bt[2], tolerance = 1e-10)
|
|
}
|
|
})
|
|
|
|
test_that("clopper_pearson interval is valid (low <= observed <= high)", {
|
|
result <- clopper_pearson(5, 10)
|
|
expect_lte(result["low"], result["observed"])
|
|
expect_gte(result["high"], result["observed"])
|
|
})
|
|
|
|
test_that("clopper_pearson handles edge cases", {
|
|
# num = 0: lower bound should be exactly 0
|
|
zero_result <- clopper_pearson(0, 10)
|
|
expect_equal(unname(zero_result["low"]), 0, tolerance = 1e-15)
|
|
expect_gt(unname(zero_result["high"]), 0)
|
|
# Upper bound should match binom.test exactly
|
|
bt_zero <- binom.test(0, 10)$conf.int[2]
|
|
expect_equal(unname(zero_result["high"]), unname(bt_zero), tolerance = 1e-10)
|
|
|
|
# num = den: upper bound should be exactly 1 (observed is at the boundary)
|
|
full_result <- clopper_pearson(10, 10)
|
|
expect_equal(unname(full_result["high"]), 1, tolerance = 1e-15)
|
|
# Lower bound for all-successes case is well below observed (asymmetric interval)
|
|
bt_full <- binom.test(10, 10)$conf.int
|
|
expect_equal(unname(full_result["low"]), unname(bt_full[1]), tolerance = 1e-10)
|
|
})
|
|
|
|
test_that("clopper_pearson respects conf.level", {
|
|
wide <- clopper_pearson(20, 40, conf.level = 0.99)
|
|
narrow <- clopper_pearson(20, 40, conf.level = 0.80)
|
|
|
|
# Higher confidence level produces a wider interval
|
|
expect_gt(unname(wide["high"]) - unname(wide["low"]),
|
|
unname(narrow["high"]) - unname(narrow["low"]))
|
|
})
|
|
|
|
|
|
# --- z_univariate ------------------------------------------------------------
|
|
|
|
test_that("z_univariate returns a single numeric value", {
|
|
result <- z_univariate(0.13, 0.11, 2500)
|
|
expect_type(result, "double")
|
|
expect_length(result, 1)
|
|
})
|
|
|
|
test_that("z_univariate equals the formula by hand calculation", {
|
|
# z = (p_hat - p_0) / sqrt(p_0 * (1 - p_0) / n)
|
|
unit_prop <- 0.13
|
|
global_prop <- 0.11
|
|
unit_denom <- 2500
|
|
|
|
expected <- (unit_prop - global_prop) /
|
|
sqrt((global_prop * (1 - global_prop)) / unit_denom)
|
|
|
|
result <- z_univariate(unit_prop, global_prop, unit_denom)
|
|
expect_equal(result, expected, tolerance = 1e-12)
|
|
})
|
|
|
|
test_that("z_univariate is zero when proportions are equal", {
|
|
expect_equal(z_univariate(0.5, 0.5, 100), 0, tolerance = 1e-15)
|
|
})
|
|
|
|
test_that("z_univariate sign follows the direction of deviation", {
|
|
# When unit_prop > global_prop, z should be positive
|
|
expect_gt(z_univariate(0.2, 0.1, 100), 0)
|
|
|
|
# When unit_prop < global_prop, z should be negative
|
|
expect_lt(z_univariate(0.1, 0.2, 100), 0)
|
|
})
|
|
|
|
|
|
# --- waldInterval ------------------------------------------------------------
|
|
|
|
test_that("waldInterval returns a named numeric vector of length 2", {
|
|
result <- waldInterval(x = 20, n = 40)
|
|
expect_type(result, "double")
|
|
expect_length(result, 2)
|
|
expect_named(result, c("lwr", "upr"))
|
|
})
|
|
|
|
test_that("waldInterval matches documented example values", {
|
|
# The roxygen @examples comment says: waldInterval(x = 20, n = 40)
|
|
# returns approximately 0.345 and 0.655
|
|
result <- waldInterval(20, 40)
|
|
|
|
p_hat <- 20 / 40 # 0.5
|
|
se <- sqrt(p_hat * (1 - p_hat) / 40) # ~0.0791
|
|
z_crit <- qnorm(0.975) # ~1.96
|
|
|
|
expect_equal(unname(result["lwr"]), p_hat - z_crit * se, tolerance = 1e-12)
|
|
expect_equal(unname(result["upr"]), p_hat + z_crit * se, tolerance = 1e-12)
|
|
})
|
|
|
|
test_that("waldInterval interval is centered on the sample proportion", {
|
|
result <- waldInterval(30, 50)
|
|
midpoint <- (unname(result["lwr"]) + unname(result["upr"])) / 2
|
|
expect_equal(midpoint, 30 / 50, tolerance = 1e-12)
|
|
})
|
|
|
|
test_that("waldInterval respects conf.level", {
|
|
wide <- waldInterval(20, 40, conf.level = 0.99)
|
|
narrow <- waldInterval(20, 40, conf.level = 0.80)
|
|
|
|
expect_gt(unname(wide["upr"]) - unname(wide["lwr"]),
|
|
unname(narrow["upr"]) - unname(narrow["lwr"]))
|
|
})
|
|
|
|
|
|
# --- agresti_coull_interval --------------------------------------------------
|
|
|
|
test_that("agresti_coull_interval returns a named numeric vector of length 3", {
|
|
result <- agresti_coull_interval(20, 40)
|
|
expect_type(result, "double")
|
|
expect_length(result, 3)
|
|
expect_named(result, c("low", "observed", "high"))
|
|
})
|
|
|
|
test_that("agresti_coull_interval observed matches num/den", {
|
|
result <- agresti_coull_interval(20, 40)
|
|
expect_equal(unname(result["observed"]), 20 / 40)
|
|
})
|
|
|
|
test_that("agresti_coull_interval interval is valid (low <= observed <= high)", {
|
|
result <- agresti_coull_interval(15, 30)
|
|
expect_lte(result["low"], result["observed"])
|
|
expect_gte(result["high"], result["observed"])
|
|
})
|
|
|
|
test_that("agresti_coull_interval matches manual formula calculation", {
|
|
num <- 20
|
|
den <- 40
|
|
conf_level <- 0.95
|
|
|
|
z <- qnorm(1 - (1 - conf_level) / 2)
|
|
n_tilde <- den + z^2
|
|
p_tilde <- (num + z^2 / 2) / n_tilde
|
|
margin <- z * sqrt(p_tilde * (1 - p_tilde) / n_tilde)
|
|
|
|
result <- agresti_coull_interval(num, den, conf.level = conf_level)
|
|
expect_equal(unname(result["low"]), p_tilde - margin, tolerance = 1e-12)
|
|
expect_equal(unname(result["high"]), p_tilde + margin, tolerance = 1e-12)
|
|
})
|
|
|
|
test_that("agresti_coull_interval respects conf.level", {
|
|
wide <- agresti_coull_interval(20, 40, conf.level = 0.99)
|
|
narrow <- agresti_coull_interval(20, 40, conf.level = 0.80)
|
|
|
|
expect_gt(unname(wide["high"]) - unname(wide["low"]),
|
|
unname(narrow["high"]) - unname(narrow["low"]))
|
|
})
|