10 Group sequential design boundary with MaxCombo test

In this chapter, we calculate boundaries based on the error spending approach following Gordon Lan and DeMets (1983).

We continue use the same example scenario in the last chapter.

enrollRates <- tibble::tibble(Stratum = "All", duration = 12, rate = 500 / 12)

failRates <- tibble::tibble(
  Stratum = "All",
  duration = c(4, 100),
  failRate = log(2) / 15, # Median survival is 15 months
  hr = c(1, .6),
  dropoutRate = 0.001
)
# Randomization Ratio is 1:1
ratio <- 1

# Type I error (one-sided)
alpha <- 0.025

# Power (1 - beta)
beta <- 0.2
power <- 1 - beta

# Interim analysis time
analysisTimes <- c(12, 24, 36)
  • Example 1 in the previous chapter
fh_test <- rbind(
  data.frame(
    rho = 0, gamma = 0, tau = -1,
    test = 1,
    Analysis = 1:3,
    analysisTimes = c(12, 24, 36)
  ),
  data.frame(
    rho = c(0, 0.5), gamma = 0.5, tau = -1,
    test = 2:3,
    Analysis = 3, analysisTimes = 36
  )
)
fh_test
#>   rho gamma tau test Analysis analysisTimes
#> 1 0.0   0.0  -1    1        1            12
#> 2 0.0   0.0  -1    1        2            24
#> 3 0.0   0.0  -1    1        3            36
#> 4 0.0   0.5  -1    2        3            36
#> 5 0.5   0.5  -1    3        3            36

10.1 Sample size calculation based on spending function

  • Implementation in gsdmvn
gs_design_combo(enrollRates,
  failRates,
  fh_test,
  alpha = 0.025,
  beta = 0.2,
  ratio = 1,
  binding = FALSE, # test.type = 4 non-binding futility bound
  upper = gs_spending_combo,
  upar = list(sf = gsDesign::sfLDOF, total_spend = 0.025), # alpha spending
  lower = gs_spending_combo,
  lpar = list(sf = gsDesign::sfLDOF, total_spend = 0.2), # beta spending
) %>%
  mutate(`Probability_Null (%)` = Probability_Null * 100) %>%
  select(-Probability_Null) %>%
  mutate_if(is.numeric, round, digits = 2)
#>   Analysis Bound Time      N Events     Z Probability Probability_Null (%)
#> 1        1 Upper   12 301.31  64.72  6.18        0.00                 0.00
#> 3        2 Upper   24 301.31 148.41  2.80        0.22                 0.26
#> 5        3 Upper   36 301.31 199.64  2.10        0.80                 2.49
#> 2        1 Lower   12 301.31  64.72 -2.72        0.00                   NA
#> 4        2 Lower   24 301.31 148.41  0.65        0.08                   NA
#> 6        3 Lower   36 301.31 199.64  2.10        0.20                   NA
  • Simulation results based on 10,000 replications.
#>     n      d analysisTimes lower upper
#> 7 302  64.78            12  0.00  0.00
#> 8 302 148.75            24  0.08  0.22
#> 9 302 200.18            36  0.20  0.80
  • Compared with \(FH(0, 0)\) using boundary based on test.type = 4
x <- gsdmvn::gs_design_wlr(
  enrollRates = enrollRates, failRates = failRates,
  ratio = ratio, alpha = alpha, beta = beta,
  weight = function(x, arm0, arm1) {
    gsdmvn::wlr_weight_fh(x, arm0, arm1, rho = 0, gamma = 0)
  },
  upper = gs_spending_bound,
  lower = gs_spending_bound,
  upar = list(sf = gsDesign::sfLDOF, total_spend = alpha),
  lpar = list(sf = gsDesign::sfLDOF, total_spend = beta),
  analysisTimes = analysisTimes
)$bounds %>%
  mutate_if(is.numeric, round, digits = 2)
x
#> # A tibble: 6 × 11
#>   Analysis Bound  Time     N Events     Z Probability   AHR theta  info info0
#>      <dbl> <chr> <dbl> <dbl>  <dbl> <dbl>       <dbl> <dbl> <dbl> <dbl> <dbl>
#> 1        1 Upper    12  377.    81   3.79        0     0.84  0.17  20.2  20.3
#> 2        2 Upper    24  377.   186.  2.36        0.46  0.72  0.33  46.3  46.8
#> 3        3 Upper    36  377.   250.  2.01        0.8   0.68  0.38  61.8  63.3
#> 4        1 Lower    12  377.    81  -1.18        0.03  0.84  0.17  20.2  20.3
#> 5        2 Lower    24  377.   186.  1.15        0.14  0.72  0.33  46.3  46.8
#> 6        3 Lower    36  377.   250.  2.01        0.2   0.68  0.38  61.8  63.3
  • Compared with \(FH(0.5, 0.5)\) using boundary based on test.type = 4
x <- gsdmvn::gs_design_wlr(
  enrollRates = enrollRates, failRates = failRates,
  ratio = ratio, alpha = alpha, beta = beta,
  weight = function(x, arm0, arm1) {
    gsdmvn::wlr_weight_fh(x, arm0, arm1, rho = 0.5, gamma = 0.5)
  },
  upper = gs_spending_bound,
  lower = gs_spending_bound,
  upar = list(sf = gsDesign::sfLDOF, total_spend = alpha),
  lpar = list(sf = gsDesign::sfLDOF, total_spend = beta),
  analysisTimes = analysisTimes
)$bounds %>%
  mutate_if(is.numeric, round, digits = 2)
x
#> # A tibble: 6 × 11
#>   Analysis Bound  Time     N Events     Z Probability   AHR theta  info info0
#>      <dbl> <chr> <dbl> <dbl>  <dbl> <dbl>       <dbl> <dbl> <dbl> <dbl> <dbl>
#> 1        1 Upper    12  304.   65.3  5.05        0     0.79  0.68  1.76  1.77
#> 2        2 Upper    24  304.  150.   2.51        0.42  0.68  0.93  6.17  6.28
#> 3        3 Upper    36  304.  201.   1.99        0.8   0.65  0.97  9.16  9.44
#> 4        1 Lower    12  304.   65.3 -1.8         0     0.79  0.68  1.76  1.77
#> 5        2 Lower    24  304.  150.   1.12        0.12  0.68  0.93  6.17  6.28
#> 6        3 Lower    36  304.  201.   1.99        0.2   0.65  0.97  9.16  9.44
  • Compared with \(FH(0, 0.5)\) using boundary based on test.type = 4
x <- gsdmvn::gs_design_wlr(
  enrollRates = enrollRates, failRates = failRates,
  ratio = ratio, alpha = alpha, beta = beta,
  weight = function(x, arm0, arm1) {
    gsdmvn::wlr_weight_fh(x, arm0, arm1, rho = 0, gamma = 0.5)
  },
  upper = gs_spending_bound,
  lower = gs_spending_bound,
  upar = list(sf = gsDesign::sfLDOF, total_spend = alpha),
  lpar = list(sf = gsDesign::sfLDOF, total_spend = beta),
  analysisTimes = analysisTimes
)$bounds %>%
  mutate_if(is.numeric, round, digits = 2)
x
#> # A tibble: 6 × 11
#>   Analysis Bound  Time     N Events     Z Probability   AHR theta  info info0
#>      <dbl> <chr> <dbl> <dbl>  <dbl> <dbl>       <dbl> <dbl> <dbl> <dbl> <dbl>
#> 1        1 Upper    12  288.   61.9  6.18        0     0.78  0.63  2.08  2.09
#> 2        2 Upper    24  288.  142    2.8         0.3   0.67  0.77  8.86  9.07
#> 3        3 Upper    36  288.  191.   1.97        0.8   0.64  0.73 15.7  16.4 
#> 4        1 Lower    12  288.   61.9 -2.43        0     0.78  0.63  2.08  2.09
#> 5        2 Lower    24  288.  142    0.93        0.09  0.67  0.77  8.86  9.07
#> 6        3 Lower    36  288.  191.   1.97        0.2   0.64  0.73 15.7  16.4

10.2 Information fraction

Under the same design assumption, information fraction is different from different weight parameters in WLR.

There are two potential strategies to calculate the information fraction:

  • Option 1: Use minimal information fraction of all candidate tests (implemented in gsdmvn).
  • Option 2: Use the weighted average of information fraction.
utility <- gsdmvn:::gs_utility_combo(enrollRates, failRates, fh_test = fh_test, ratio = 1)
utility$info_all %>% mutate_if(is.numeric, round, digits = 3)
#>   test Analysis Time   N  Events   AHR  delta sigma2 theta   info  info0
#> 1    1        1   12 500 107.394 0.842 -0.009  0.054 0.172 26.841 26.899
#> 2    1        2   24 500 246.283 0.716 -0.041  0.123 0.333 61.352 62.087
#> 3    1        3   36 500 331.291 0.683 -0.062  0.164 0.381 81.918 83.944
#> 4    2        1   12 500 107.394 0.781 -0.005  0.007 0.626  3.605  3.624
#> 5    2        2   24 500 246.283 0.666 -0.024  0.031 0.765 15.371 15.737
#> 6    2        3   36 500 331.291 0.639 -0.040  0.054 0.732 27.205 28.484
#> 7    3        1   12 500 107.394 0.790 -0.004  0.006 0.676  2.898  2.910
#> 8    3        2   24 500 246.283 0.675 -0.019  0.020 0.927 10.155 10.335
#> 9    3        3   36 500 331.291 0.650 -0.029  0.030 0.973 15.070 15.526
info_frac <- tapply(utility$info_all$info0, utility$info_all$test, function(x) x / max(x))
info_frac
#> $`1`
#> [1] 0.3204395 0.7396239 1.0000000
#> 
#> $`2`
#> [1] 0.1272245 0.5525079 1.0000000
#> 
#> $`3`
#> [1] 0.1874599 0.6656420 1.0000000
min_info_frac <- apply(do.call(rbind, info_frac), 2, min)
min_info_frac
#> [1] 0.1272245 0.5525079 1.0000000

10.3 Spending function

  • gsDesign bound can not be directly used for MaxCombo test:

    • Multiple test statistics are considered in interim analysis or final analysis.
  • Design example:

    • \(\alpha = 0.025\)
    • \(\beta = 0.2\)
    • \(K=3\) total analysis.
    • test.type=4: A two-sided, asymmetric, beta-spending with non-binding lower bound.
  • Lan-DeMets spending function to approximate an O’Brien-Fleming bound Gordon Lan and DeMets (1983). (gsDesign::sfLDOF)

  • \(t\) is information fraction in the formula below.

\[f(t; \alpha)=2-2\Phi\left(\Phi^{-1}\left(\frac{1-\alpha/2}{t^{\rho/2}}\right)\right)\]

  • Other spending functions are discussed in the gsDesign technical manual and implemented in gsDesign::sf*() functions.

  • Upper bound:

    • Non-binding lower bound: lower bound are all -Inf.
alpha_spend <- gsDesign::sfLDOF(alpha = 0.025, t = min_info_frac)$spend
alpha_spend
#> [1] 3.300225e-10 2.566068e-03 2.500000e-02
  • Lower bound:
beta_spend <- gsDesign::sfLDOF(alpha = 0.2, t = min_info_frac)$spend
beta_spend
#> [1] 0.0003269604 0.0846866352 0.2000000000

10.4 Technical details

This section describes the details of calculating the sample size and events required for WLR under group sequential design. It can be skipped for the first read of this training material.

10.5 Upper bound in group sequential design

10.5.1 One test in each interim analysis

  • First interim analysis

\[\alpha_1 = \text{Pr}(Z_1 > b_1 \mid H_0)\]

qnorm(1 - alpha_spend[1])
#> [1] 6.17538
  • General formula (non-binding futility bound)

\[\alpha_k = \text{Pr}(\cap_{i=1}^{i=k-1} Z_i < b_i, Z_k > b_k \mid H_0)\]

10.5.2 MaxCombo upper bound

  • First interim analysis upper bound

\[\alpha_1 = \text{Pr}(G_1 > b_1 \mid H_0) = 1 - \text{Pr}(\cap_i Z_{i1} < b_1 \mid H_0)\]

  • General formula (non-binding futility bound)

\[\alpha_k = \text{Pr}(\cap_{i=1}^{i=k-1} G_i < b_i, G_k > b_k \mid H_0)\] \[ = \text{Pr}(\cap_{i=1}^{i=k-1} G_i < b_i \mid H_0) - \text{Pr}(\cap_{i=1}^{i=k} G_i < b_i \mid H_0)\]

  • gsdmvn implementation for MaxCombo test
gsdmvn:::gs_bound(alpha_spend, beta_spend,
  analysis = utility$info$Analysis, # Analysis indicator
  theta = rep(0, nrow(fh_test)), # Under the null hypothesis
  corr = utility$corr # Correlation
)$upper
#> [1] 6.175397 2.798651 2.097603
  • Compared with upper bound calculated from gsDesign
x <- gsDesign::gsSurv(
  k = 3, test.type = 4, alpha = 0.025,
  beta = 0.2, astar = 0, timing = c(1),
  sfu = sfLDOF, sfupar = c(0), sfl = sfLDOF,
  sflpar = c(0), lambdaC = c(0.1),
  hr = 0.6, hr0 = 1, eta = 0.01,
  gamma = c(10),
  R = c(12), S = NULL,
  T = 36, minfup = 24, ratio = 1
)
x$upper$bound
#> [1] 3.710303 2.511407 1.992970

10.6 Lower bound in group sequential design

10.6.1 One test in each interim analysis

  • First interim analysis

\[\beta_1 = \text{Pr}(Z_1 < a_1 \mid H_1)\]

n <- 400
qnorm(beta_spend[1], mean = utility$theta[1] * sqrt(n))
#> [1] -2.610665
n <- 500
qnorm(beta_spend[1], mean = utility$theta[1] * sqrt(n))
#> [1] -2.516528
  • General formula \[\beta_k = \text{Pr}(\cap_{i=1}^{i=k-1} a_i < Z_i < b_i, Z_k < a_k \mid H_1)\]

10.6.2 MaxCombo lower bound

  • First interim analysis upper bound

\[\beta_1 = \text{Pr}(G_1 < a_1 \mid H_1) = \text{Pr}(\cap_i Z_{i1} < a_1 \mid H_0)\]

  • General formula (non-binding futility bound)

\[\beta_k = \text{Pr}(\cap_{i=1}^{i=k-1} a_i <G_i < b_i, G_k < a_k \mid H_0)\]

  • gsdmvn implementation for MaxCombo test
n <- 400
utility <- gsdmvn:::gs_utility_combo(enrollRates, failRates, fh_test = fh_test, ratio = 1)
bound <- gsdmvn:::gs_bound(alpha_spend, beta_spend,
  analysis = utility$info$Analysis, # Analysis indicator
  theta = utility$theta * sqrt(n), # Under the alternative hypothesis
  corr = utility$corr # Correlation
)

bound$lower
#> [1] -2.6106678  0.9619422  2.0968098
  • Compared with lower bound calculated from gsDesign
x$lower$bound
#> [1] -0.2361874  1.1703638  1.9929702

10.7 Sample size calculation

  • Sample size and boundaries are set simultaneously using an iterative algorithm.

10.7.1 Initiate the calculation from lower bound derived at \(N = 400\)

bound
#>      upper      lower
#> 1 6.175397 -2.6106678
#> 2 2.798651  0.9619422
#> 3 2.096810  2.0968098
gs_design_combo(enrollRates,
  failRates,
  fh_test,
  alpha = 0.025,
  beta = 0.2,
  ratio = 1,
  binding = FALSE,
  upar = bound$upper,
  lpar = bound$lower,
) %>% mutate_if(is.numeric, round, digits = 2)
#>   Analysis Bound Time      N Events     Z Probability Probability_Null
#> 1        1 Upper   12 324.88  69.78  6.18        0.00             0.00
#> 3        2 Upper   24 324.88 160.03  2.80        0.24             0.00
#> 5        3 Upper   36 324.88 215.26  2.10        0.80             0.03
#> 2        1 Lower   12 324.88  69.78 -2.61        0.00               NA
#> 4        2 Lower   24 324.88 160.03  0.96        0.13               NA
#> 6        3 Lower   36 324.88 215.26  2.10        0.20               NA

10.7.2 Update bound based on newly calculated sample size

n <- 355
utility <- gsdmvn:::gs_utility_combo(enrollRates, failRates, fh_test = fh_test, ratio = 1)
bound <- gsdmvn:::gs_bound(alpha_spend, beta_spend,
  analysis = utility$info$Analysis, # Analysis indicator
  theta = utility$theta * sqrt(n), # Under the alternative hypothesis
  corr = utility$corr # Correlation
)

bound
#>      upper      lower
#> 1 6.175397 -2.6568642
#> 2 2.798651  0.8266125
#> 3 2.097617  2.0976166
  • Repeat the procedure above until the sample size and lower bound converge.