Untitled

# Define the function to calculate the power of the test
power_test <- function(n, p, c) {
  pbinom(c - 1, size = n, prob = p, lower.tail = FALSE)
}

# Define the quantile function for the normal approximation
qnorm_approximation <- function(alpha, n, p) {
  qnorm(1 - alpha) * sqrt(n * p * (1 - p)) + n * p
}

# Initial parameters
p1 <- 1/2
p2 <- 2/3
alpha <- 0.10
beta <- 1 - 0.95

# Approximate the sample size n and critical value c
n_approx <- (qnorm_approximation(alpha, 1, p1) + qnorm_approximation(beta, 1, p2))^2 / (2 * (p2 - p1))^2
n <- ceiling(n_approx) # Sample size should be an integer
c_approx <- (n + qnorm_approximation(alpha, n, p1))/2
c <- ceiling(c_approx) # Critical value should be an integer

# Calculate the actual power at p=1/2 and p=2/3 using the binomial distribution
power_at_p1 <- power_test(n, p1, c)
power_at_p2 <- power_test(n, p2, c)

# Print the results
cat("Sample size n:", n, "\n")
cat("Critical value c:", c, "\n")
cat("Power at p=1/2:", power_at_p1, "\n")
cat("Power at p=2/3:", power_at_p2, "\n")

# If the power at p=1/2 is greater than desired, increase c by 1 and recalculate
if (power_at_p1 > alpha) {
  c <- c + 1
  power_at_p1 <- power_test(n, p1, c)
  power_at_p2 <- power_test(n, p2, c)
  cat("Adjusted critical value c:", c, "\n")
  cat("Adjusted power at p=1/2:", power_at_p1, "\n")
  cat("Adjusted power at p=2/3:", power_at_p2, "\n")
}