## -----------------------------------------------------------------------------
library(ssPilot)


## -----------------------------------------------------------------------------
standard_sample_size(
  sd = 1,
  effect = 0.5,
  power = 0.90,
  alpha = 0.05
)


## -----------------------------------------------------------------------------
pilot_n <- 12
k <- 2 * pilot_n - 2
k


## -----------------------------------------------------------------------------
sd <- 1
variance <- sd^2
conf_level <- 0.80

variance_ucl <-
  k * variance /
  qchisq(1 - conf_level, df = k)

variance_ucl


## -----------------------------------------------------------------------------
sd_ucl <- sqrt(variance_ucl)

sd_ucl


## -----------------------------------------------------------------------------
effect <- 0.50
power <- 0.90
alpha <- 0.05
allocation <- 1

main_n <- standard_sample_size(
  sd = sd_ucl,
  effect = effect,
  power = power,
  alpha = alpha,
  allocation = allocation
)

main_n


## -----------------------------------------------------------------------------
pilot_n <- 12
sd <- 1
conf_level <- 0.80
effect <- 0.50
power <- 0.90
alpha <- 0.05
allocation <- 1

# Degrees of freedom
k <- 2 * pilot_n - 2

# Upper confidence limit for the variance
variance_ucl <-
  k * sd^2 /
  qchisq(1 - conf_level, df = k)

# Upper confidence limit for the standard deviation
sd_ucl <- sqrt(variance_ucl)

# Main-trial sample size
main_n <- standard_sample_size(
  sd = sd_ucl,
  effect = effect,
  power = power,
  alpha = alpha,
  allocation = allocation
)

main_n


## -----------------------------------------------------------------------------
ucl_sample_size(
  pilot_n = 12,
  sd = 1,
  effect = 0.50,
  power = 0.90,
  alpha = 0.05,
  conf_level = 0.80
)


## -----------------------------------------------------------------------------
pilot_sizes <- 10:50

pilot_sizes


## -----------------------------------------------------------------------------
sd <- 1
conf_level <- 0.80

variance_ucl <- sapply(pilot_sizes, function(m) {
  k <- 2 * m - 2
  k * sd^2 /
    qchisq(1 - conf_level, df = k)
})

sd_ucl <- sqrt(variance_ucl)
sd_ucl


## -----------------------------------------------------------------------------
effect <- 0.50
power <- 0.90
alpha <- 0.05

main_sizes <- sapply(sd_ucl, function(s) {
  standard_sample_size(
    sd = s,
    effect = effect,
    power = power,
    alpha = alpha,
    allocation = 1
  )
})

main_sizes


## -----------------------------------------------------------------------------
total_sizes <- pilot_sizes + main_sizes

total_sizes


## -----------------------------------------------------------------------------
optimal_index <- which.min(total_sizes)
optimal_pilot_n <- pilot_sizes[optimal_index]
optimal_main_n <- main_sizes[optimal_index]
optimal_total_n <- total_sizes[optimal_index]

optimal_pilot_n
optimal_main_n
optimal_total_n


## -----------------------------------------------------------------------------
optimization_results <- data.frame(
  pilot_n_per_arm = pilot_sizes,
  main_n_per_arm = main_sizes,
  total_n_per_arm = total_sizes
)

optimization_results


## -----------------------------------------------------------------------------
sd <- 1
effect <- 0.50
power <- 0.90
alpha <- 0.05
conf_level <- 0.80

standard_result <- standard_sample_size(
  sd = sd,
  effect = effect,
  power = power,
  alpha = alpha
)

ucl_result <- ucl_sample_size(
  pilot_n = 12, # as suggested by Julious (2005)
  sd = sd,
  effect = effect,
  power = power,
  alpha = alpha,
  conf_level = conf_level
)

optimized_ucl_result <- optimized_ucl_sample_size(
  sd = sd,
  effect = effect,
  power = power,
  alpha = alpha,
  conf_level = conf_level
)




## -----------------------------------------------------------------------------
comparison <- data.frame(
  method = c(
    "Standard",
    "UCL",
    "Optimized UCL"
  ),
  main_n_per_arm = c(
    standard_result,
    ucl_result$main_n_per_arm,
    optimized_ucl_result$main_n_per_arm
  )
)

comparison