R/CABCD.R

In covadap: Implement Covariate-Adaptive Randomization

# Covariate-Adjusted Biased Coin Design
CABCD <- function(data,
                  a = 3,
                  print.results = TRUE) {
  if (missing(data)) {
    stop("data must be provided")
  }
  if (!is.numeric(a) | a <= 0) {
    stop("a must be a positive number")
  }

  dr = data_preproc(data)
  n = dr$n
  if(n < 3){
    stop("at least 3 observations are needed")
  }
  delta = numeric(n)
  A.strata = N.strata = D.strata = numeric(dr$n.strata)
  strata = dr$strata
  obs.strata = dr$obs.strata
  n_cov = dr$n_cov
  n_levels = dr$n_levels
  for (i in 1:n) {
    sel = strata[i]
    Dca = D.strata[sel]
    if (Dca == 0) {
      delta[i] = rbinom(1, 1, .5)
    } else if (Dca > 0) {
      delta[i] =  rbinom(1, 1, 1 / (Dca ^ a + 1))
    } else{
      delta[i] =  rbinom(1, 1, (abs(Dca) ^ a) / (abs(Dca) ^ a + 1))
    }
    N.strata[sel] = N.strata[sel] + 1
    A.strata[sel] = A.strata[sel] + delta[i]
    D.strata = 2 * A.strata - N.strata
  }
  num_strata = dr$num_strata
  whit.cov = list()
  for (nc in 1:n_cov) {
    imb.temp = numeric()
    for (j in 1:n_levels[nc]) {
      N.tot = sum(num_strata[, nc] == j)
      A.tot = sum(delta[num_strata[, nc] == j])
      imb.temp[j] = 2*A.tot - N.tot
    }
    whit.cov[[nc]] = imb.temp
    names(whit.cov[[nc]]) = paste(dr$var_levels[[nc]][-1])
  }
  names(whit.cov) = colnames(dr$dm)
  Imb.measures = Imb.m(delta, dr$Fn)

  res = list(
    summary.info = list(
      Design = "CABCD",
      Sample_size = n,
      n_cov = n_cov,
      n_levels = paste(n_levels, collapse = " "),
      var_names = paste(colnames(obs.strata), collapse = " "),
      cov_levels_names = dr$var_levels,
      parameter_a = a
    ),
    Assignments = factor(delta, labels = c("B", "A"), levels = 0:1),
    Imbalances.summary = Imb.measures,
    Strata.measures = cbind(obs.strata, N.strata, A.strata, D.strata),
    Imbalances = list(
      Imb.measures = Imb.measures[-3],
      Overall.imb = Imb.measures[[3]],
      Within.strata = cbind(obs.strata, D.strata),
      Within.cov = whit.cov
    ),
    data = data,
    observed.strata = obs.strata
  )
  class(res) = "covadap"
  if (print.results) {
    print_covadap(res)
  }
  invisible(res)
}