context("Binomial data distribution")
test_that("Constructor works", {
expect_true(
Binomial(.2)@two_armed)
expect_true(
!Binomial(.4, two_armed = FALSE)@two_armed)
expect_error(
Binomial(-.1))
expect_error(
Binomial(1.1))
}) # end 'constructor works'
test_that("pdf is defined correctly", {
dist <- Binomial(.1, TRUE)
x <- seq(-3, 3, by = .1)
n <- 22
theta <- .35
r_A <- theta + dist@rate_control
sigma_A <- sqrt(sum(c(r_A * (1 - r_A), dist@rate_control * (1 - dist@rate_control))))
sigma_0 <- sqrt(2 * mean(c(r_A, dist@rate_control)) * (1 - mean(c(r_A, dist@rate_control))))
expect_equal(
probability_density_function(dist, x, n, theta),
stats::dnorm(x, mean = sqrt(n) * theta / sigma_0, sd = sigma_A / sigma_0),
tolerance = 1e-6, scale = 1)
})
test_that("cdf is defined correctly", {
dist <- Binomial(.3, TRUE)
x <- seq(-3, 3, by = .1)
n <- 120
theta <- .02
r_A <- theta + dist@rate_control
sigma_A <- sqrt(sum(c(r_A * (1 - r_A), dist@rate_control * (1 - dist@rate_control))))
sigma_0 <- sqrt(2 * mean(c(r_A, dist@rate_control)) * (1 - mean(c(r_A, dist@rate_control))))
expect_equal(
cumulative_distribution_function(dist, x, n, theta),
stats::pnorm(x, mean = sqrt(n) * theta / sigma_0, sd = sigma_A / sigma_0),
tolerance = 1e-6, scale = 1)
})
test_that("quantile is defined correctly", {
dist <- Binomial(.6, FALSE)
probs <- seq(.01, 1, by = .01)
n <- 99
theta <- .125
r_A <- theta + dist@rate_control
sigma_A <- sqrt(r_A * (1 - r_A))
sigma_0 <- sqrt(2 * r_A * (1 - r_A))
expect_equal(
quantile(dist, probs, n, theta),
stats::qnorm(probs, mean = sqrt(n) * theta / sigma_0, sd = sigma_A / sigma_0),
tolerance = 1e-6, scale = 1)
expect_warning(
quantile(dist, -1, n, theta)
)
})
test_that("vectorization works", {
dist <- Binomial(.3)
prior <- PointMassPrior(c(.3, .4), c(.4, .6))
expect_equal(
predictive_pdf(dist, prior, c(0, 1), 10),
sapply(c(0, 1), function(x) predictive_pdf(dist, prior, x, 10))
)
})
# Let us check if the statistical tests are working correctly
r_c <- .4
r_e <- .6
r_diff <- r_e - r_c
dist <- Binomial(r_c)
alpha <- 0.025
beta <- 0.2
n <- ceiling(rpact::getSampleSizeRates(alpha = alpha, beta = beta,
groups = 2, normalApproximation = TRUE,
pi1 = r_e, pi2 = r_c)$nFixed1)
c <- qnorm(1 - alpha)
test_that("Approximation works in cdf", {
iters <- 1e6
p_c <- stats::rbinom(iters, n, r_c) / n
p_e <- stats::rbinom(iters, n, r_e) / n
p_e_2 <- stats::rbinom(iters, n, r_c) / n
sigma_A <- sqrt(p_e * (1 - p_e) + p_c * (1 - p_c))
sigma_0 <- sqrt(p_e_2 * (1 - p_e_2) + p_c * (1 - p_c))
expect_equal(
mean(1 - cumulative_distribution_function(dist, sqrt(n) * (p_e - p_c) / sigma_A, n, 0) <= alpha),
0.807, # reference: M. Kieser (2018), Fallzahlberechnung in der medizinischen Forschung, Tab. 5.2
tolerance = 1e-2, scale = 1
)
expect_equal(
mean(1 - cumulative_distribution_function(dist, sqrt(n) * (p_e_2 - p_c) / sigma_0, n, 0) <= alpha),
alpha,
tolerance = 1 / sqrt(iters), scale = 1
)
})
test_that("Approximation works in evaluate", {
design <- OneStageDesign(n, c)
H_0 <- PointMassPrior(0, 1)
H_1 <- PointMassPrior(r_diff, 1)
expect_equal(
evaluate(Power(dist, H_1), design),
1 - beta,
tolerance = 1e-3, scale = 1
)
expect_equal(
evaluate(Power(dist, H_0), design),
alpha,
tolerance = 1e-3, scale = 1
)
})
test_that("Approximation works in simulate", {
iters <- 1e6
test_statistic_H0 <- simulate(object = dist, nsim = iters, n = n, theta = 0, seed = 42)
test_statistic_H1 <- simulate(object = dist, nsim = iters, n = n, theta = r_diff, seed = 42)
expect_equal(
mean(test_statistic_H1 > c),
1 - beta,
tolerance = 1 / sqrt(iters), scale = 1
)
expect_equal(
mean(test_statistic_H0 > c),
alpha,
tolerance = 1 / sqrt(iters), scale = 1
)
})
test_that("show method", {
expect_equal(
capture.output(show(Binomial(.2))),
"Binomial<two-armed>, response rate in control group: 0.2 "
)
expect_equal(
capture.output(show(Binomial(.8, two_armed = FALSE))),
"Binomial<single-armed>, response rate in control group: 0.8 "
)
})
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.