Bioequivalence and Bioavailability Forum

Dear all,

have you ever wondered how to find a suitable adjusted α if the desired combination of target power, expected T/R-ratio (and even the acceptance range) is not given in any of the many publications? In such a case we have to validate the framework within the desired range of n₁/CV-combinations in order to demonstrate that the overall type I error is preserved.
Here is my recipe – inspired by Fuglsang’s¹ ideas:

1. Explore the desired grid of n₁/CV-combinations.
   Use Pocock’s α 0.0294 as a starter; run 10⁵ simulations for speed reasons.
   Find the location of maximum inflation. For variants of “Method C” the maximum inflation is generally seen at the minimum n₁ and CV ~10–25%; for variants of “Method B” it can be anywhere…
   
2. At this location simulate a vector of adjusted alphas (10⁶ sim’s). Lower limit 0.025 (Bonferroni), upper limit α_max/2+10% (a feasible guess in most cases).
   Perform a linear regression of type I errors vs. α_adj (THX to Anders!); the intersection of the regression line with 0.05 is the “best” α_adj.
   
3. Validate the grid with the chosen α_adj. I suggest a narrower grid of n₁/CV-combinations (at least around the planned n₁ and expected CV). 10⁶ sim’s per grid-point are mandatory.

In R by means of package Power2Stage:

library(Power2Stage)
##############
# First step #
##############
n1.min  <- 12
n1.max  <- 60
n1.step <- 12
CV.min  <- 0.1
CV.max  <- 0.8
CV.step <- 0.05
method  <- "B"
GMR     <- 0.95 # assumed T/R-ratio
theta0  <- 1.25 # upper BE-limit for simulating alpha
target  <- 0.80 # target power (generally 0.8–0.9)
usePE   <- FALSE
alpha   <- c(0.0294, 0.0294)
           # Pocock's alpha.
           # Change only if you know what you are doing!
Nmax    <- Inf
           # No upper sample size limit (futility).
           # Change only if you know what you are doing!
           # Generally alpha is not affected, but have a
           # close look at power!
# Potvin's defaults: n1 12-60 (step 12), CV 0.1-1.0 (step 0.1), GMR 0.95, target 0.8
n1 <- seq(n1.min, n1.max, by=n1.step)
CV <- seq(CV.min, CV.max, by=CV.step)
nsims <- 1e5
pBE <- NULL
for(j in 1:length(CV)){
  for(k in 1:length(n1)){
    pBE <- c(pBE, power.2stage(method=method, alpha=alpha,
            n1=n1[k], GMR=GMR, CV=CV[j], targetpower=target,
            pmethod="nct", usePE=usePE, Nmax=Nmax, theta0=theta0,
            nsims=nsims, print=F, details=F)$pBE)
  }
}
grid <- matrix(data=pBE, nrow=length(CV), ncol=length(n1), byrow=T,
  dimnames=list(CV, n1))
pos <- as.numeric(which(grid == max(grid), arr.ind=TRUE))
CV.pos <- as.numeric(rownames(grid)[pos[1]])
n1.pos <- as.numeric(colnames(grid)[pos[2]])
alpha.max <- grid[pos[1], pos[2]]
grid
cat(paste0("\nFor adj. \u03b1 of ", alpha[1], " a maximum inflation of ",
  alpha.max, " is detected at CV=", CV.pos*100, "% and n1=", n1.pos, ".\n"))
LL <- 0.025           # Bonferroni (most conservative)
UL <- alpha.max/2*1.1 # guess!
if(UL <= alpha[1]) UL <- alpha[1]*1.1 # modification for small inflation
###############
# Second step #
###############
nsims <- 1e6
x <- seq(LL, UL, length.out=10)
pBE <- NULL
for(j in 1:length(x)){
    pBE <- c(pBE, power.2stage(method=method, alpha=rep(x[j], 2),
            n1=n1.pos, GMR=GMR, CV=CV.pos, targetpower=target,
            pmethod="nct", usePE=usePE, Nmax=Nmax, theta0=theta0,
            nsims=nsims, print=F, details=F)$pBE)
}
mod <- lm(pBE ~ x)
alpha.ad <- rep(round(as.numeric((0.05-coef(mod)[1])/coef(mod)[2]), 4), 2)
sig <- binom.test(0.05*nsims, nsims, alternative="less")$conf.int[2]
plot(x, pBE, xlab="adj. alpha", ylab="emp. alpha",
  main=paste0("“Method ", method, "”, Inflation at CV = ",
  CV.pos*100, "% and n1 = ", n1.pos))
  abline(h=0.05, col="red")         # target alpha
  abline(h=sig, lty="dotted", col="red") # significance limit
  abline(mod)                       # regression line
  abline(v=alpha[1], lty="dotted")  # first guess (e.g., 0.0294)
  abline(v=alpha.ad[1], col="blue") # estimated adjusted alpha
cat("estimated adjusted \u03b1:", alpha.ad[1], "\n")
##############
# Third step #
##############
pBE <- NULL
for(j in 1:length(CV)){
  for(k in 1:length(n1)){
    pBE <- c(pBE, power.2stage(method=method, alpha=alpha.ad,
            n1=n1[k], GMR=GMR, CV=CV[j], targetpower=target,
            pmethod="nct", usePE=usePE, Nmax=Nmax, theta0=theta0,
            nsims=nsims, print=F, details=F)$pBE)
  }
}
grid <- matrix(data=pBE, nrow=length(CV), ncol=length(n1), byrow=T,
  dimnames=list(CV, n1))
pos <- as.numeric(which(grid == max(grid), arr.ind=TRUE))
CV.pos <- as.numeric(rownames(grid)[pos[1]])
n1.pos <- as.numeric(colnames(grid)[pos[2]])
alpha.max <- grid[pos[1], pos[2]]
cat(paste0("\nValidation of chosen adjusted \u03b1 (", alpha.ad[1], "):\n",
  "\u201cMethod ", method, "\u201d, assumed T/R ", GMR*100, "%, ",
  "target power ", target*100, "%\n"))
grid
cat(paste0("\nWith an adj. \u03b1 of ", alpha.ad[1], " a maximum inflation of ",
  alpha.max, " is seen at CV=", CV.pos*100, "% and n1=", n1.pos, ".\n"))
if(alpha.max <= sig){
  cat(paste0("\nSignificance limit for ", nsims, " simulations: ",
  round(sig, 5), "; no significant inflation seen.",
  "\nThe chosen adjusted \u03b1 is justified within the assessed range of ",
  "stage 1 sample sizes (", n1.min, "\u002d", n1.max,") and CVs (", CV.min*100,
  "\u002d", CV.max*100, "%).\n"))
}else{
  cat(paste0("\nSignificance limit for ", nsims, " simulations: ",
  round(sig, 5), "; significant inflation seen.\n"))
  if(alpha.max < 0.052){
    cat(paste0("However, the maximum inflation of ", alpha.max,
    " does not exceed Potvin\u2019s \u201cnegligible inflation\u201d of 0.052. ",
    "Therefore, the chosen adjusted \u03b1 is justified within the assessed range of ",
    "stage 1 sample sizes (", n1.min, "\u002d", n1.max,") and CVs (", CV.min*100,
    "\u002d", CV.max*100, "%).\n"))
  }else{
    cat(paste("The maximum inflation of", alpha.max,
    "exceeds Potvin\u2019s \u201cnegligible inflation\u201d of 0.052.",
    "Therefore, a smaller adjusted \u03b1 should be tried!\n"))
  }
}

I got:
          12      24      36      48      60
0.1  0.03995 0.02985 0.02976 0.02978 0.02988
0.15 0.05234 0.04226 0.03234 0.02987 0.02988
0.2  0.05408 0.05068 0.04561 0.03912 0.03279
0.25 0.05011 0.05369 0.05047 0.04814 0.04406
0.3  0.04262 0.05398 0.05317 0.05117 0.04947
0.35 0.03668 0.05042 0.05272 0.05364 0.05152
0.4  0.03228 0.04438 0.05266 0.05389 0.05275
0.45 0.03112 0.03727 0.04858 0.05410 0.05367
0.5  0.02959 0.03287 0.04243 0.05190 0.05319
0.55 0.02940 0.03146 0.03642 0.04705 0.05209
0.6  0.02901 0.03020 0.03251 0.04075 0.04876
0.65 0.02888 0.03017 0.03036 0.03586 0.04347
0.7  0.02849 0.03006 0.02983 0.03249 0.03820
0.75 0.02839 0.02964 0.02867 0.03095 0.03348
0.8  0.02847 0.03000 0.02953 0.03005 0.03117

For adj. α of 0.0294 a maximum inflation of 0.0541 is detected at CV=45% and n1=48.



estimated adjusted α: 0.0274

           12       24       36       48       60
0.1  0.037296 0.027386 0.027369 0.027309 0.027285
0.15 0.049151 0.039538 0.030016 0.027398 0.027285
0.2  0.051255 0.047733 0.042950 0.036368 0.030491
0.25 0.047249 0.049676 0.047458 0.044682 0.041669
0.3  0.040369 0.050304 0.049422 0.047641 0.045841
0.35 0.034418 0.047104 0.050399 0.049285 0.048116
0.4  0.031017 0.040444 0.049617 0.049995 0.049238
0.45 0.029327 0.034318 0.045297 0.049996 0.050201
0.5  0.028601 0.030571 0.039118 0.047798 0.050117
0.55 0.028188 0.029017 0.033528 0.042823 0.048769
0.6  0.027837 0.028173 0.030263 0.036796 0.045379
0.65 0.027771 0.027985 0.028479 0.032260 0.039910
0.7  0.027626 0.027742 0.027986 0.029644 0.034807
0.75 0.027706 0.027803 0.027488 0.028406 0.031237
0.8  0.027586 0.027720 0.027595 0.027846 0.028965

With an adj. α of 0.0274 a maximum inflation of 0.051255 is seen at CV=20% and n1=12.

Significance limit for 1e+06 simulations: 0.05036; significant inflation seen.
However, the maximum inflation of 0.051255 does not exceed Potvin’s “negligible inflation” of 0.052. Therefore, the chosen adjusted α is justified within the assessed range of stage 1 sample sizes (12-60) and CVs (10-80%).

Get a decent cup of coffee – it takes a while (on my machine 11 min for Step 1 and 10 min for Step 2) Depending on the chosen grid expect up to a couple of days^a for Step 3. The simple grid of 75·10⁶ sim’s took six hours to complete.
For “Method B”, T/R-ratio 0.9, 90% power I got 0.0274. Slightly larger than the 0.0269 Fuglsang reported² for “Method C”. Not surprising, since B is always more conservative than C. In other words a slightly larger α is expected to lead to a similar inflation.
BTW, for “Method C” I got 0.0268 (10 & 17 min). IMHO, a nice agreement^b (different software: C vs. R, different seeds of the pseudo-random generator, different power-methods: shifted t vs. noncentral t).

Don’t forget the third step – regulators want to see only that (EMA: “appropriate steps must be taken to preserve the overall type I error of the experiment” and “the choice of how much alpha to spend at the interim analysis is at the company’s discretion”). If you want to introduce a futility criterion (e.g., an upper total sample size or even fiddle with usePE=TRUE), simulating power is crucial in order to avoid a nasty surprise.

---

1. Fuglsang A. Controlling type I errors for two-stage bioequivalence study designs. Clin Res Regul Aff. 2011;28(4):100–5. doi:10.3109/10601333.2011.631547
   
2. Fuglsang A. Sequential Bioequivalence Trial Designs with Increased Power and Controlled Type I Error Rates. AAPS J. 2013;15(3):659–61. doi:10.1208/s12248-013-9475-5

---

1. Hint: Split the grid and run more than one R-session simultaneously. Two sessions on my machine suck up only 25% of the CPU’s resources.
   
   
   
   
2. A comparison of “Method C”, 40·10⁶ sim’s. HB = my homebrew (0.0268), AF = Anders Fuglsang (0.0269):
              12             24             36             48             60     
CV     HB     AF      HB     AF      HB     AF      HB     AF      HB     AF 
0.1  0.0488 0.0484  0.0496 0.0498  0.0499 0.0498  0.0498 0.0500  0.0497 0.0499
0.2  0.0500 0.0501  0.0468 0.0472  0.0444 0.0445  0.0456 0.0459  0.0487 0.0486
0.3  0.0392 0.0394  0.0492 0.0496  0.0485 0.0486  0.0468 0.0474  0.0451 0.0452
0.4  0.0304 0.0307  0.0394 0.0399  0.0482 0.0486  0.0491 0.0494  0.0484 0.0485
0.5  0.0279 0.0282  0.0301 0.0296  0.0379 0.0382  0.0466 0.0469  0.0489 0.0492
0.6  0.0272 0.0273  0.0276 0.0275  0.0293 0.0296  0.0357 0.0362  0.0441 0.0442
0.7  0.0271 0.0271  0.0271 0.0273  0.0273 0.0272  0.0287 0.0290  0.0339 0.0339
0.8  0.0269 0.0267  0.0269 0.0272  0.0268 0.0272  0.0271 0.0273  0.0283 0.0284
   31 (78%) smaller, 2 (5%) identical, 7 (18%) larger.