Skip to contents

Two-arm two-stage: example

two_arm_two_stage_example.Rmd
library(pwb)
nsim = 10^4
N = 2428
RD_interim = RD_final = SE_interim = SE_final = pval = rep(NA, N)
set.seed(7)

Simulate the adaptive design trial only: 


for(i in 1:nsim){
  # Stage 1:
  X0 = rbinom(N/4, 1, 0.2)
  X1 = rbinom(N/4, 1, 0.15)
  
  EE = sum(X1) # Events on Experimental arm at interim
  EN = N/4 - EE  # Non-events on Experimental arm at interim
  
  CE = sum(X0) # Events on Control arm at interim
  CN = N/4 - CE # Non-events on Control arm at interim
  
  RD_interim[i] = (EE/(EE+EN)) - (CE/(CE+CN))
  
  SE_interim[i] = sqrt(EE*EN/(EE+EN)^3 + CE*CN/(CE+CN)^3)
  
  pval[i] = pnorm(RD_interim[i]/SE_interim[i])
  
  # Stage 2:
  Y0 = rbinom(N/4, 1, 0.2)
  Y1 = rbinom(N/4, 1, 0.15)
  
  EE = sum(X1+Y1) # Events on Experimental arm at end
  EN = N/2 - EE # Non-events on Experimental arm at end
  
  CE = sum(X0+Y0) # Events on Control arm at end
  CN = N/2 - CE # Non-events on Control arm at end
  
  RD_final[i] = (EE/(EE+EN)) - (CE/(CE+CN))
  
  SE_final[i] = sqrt(EE*EN/(EE+EN)^3 + CE*CN/(CE+CN)^3)
  
}

stop_early = (pval <= 0.001)
continue = (pval > 0.001)

Plot the results:

Table of results: 

#>                         nsims Bias in RD Mean SE of RD
#> Stopped early            2090  -0.029800        0.0217
#> Not stopped early        7910   0.004240        0.0154
#> All                     10000  -0.002870        0.0167
#> All, precision weighted 10000   0.000289        0.0161

Meta-analysis

Simulate data: 

########### 4 more studies with fixed sample size
RD_fixed = SE_fixed = matrix(nrow = nsim, ncol = 4)

for(j in 1:4){
  
  for(i in 1:nsim){
    
    X0 = rbinom(N/2, 1, 0.2)
    X1 = rbinom(N/2, 1, 0.15)
    
    EE = sum(X1)
    EN = N/2 - EE
    
    CE = sum(X0)
    CN = N/2 - CE
    
    RD_fixed[i,j] = (EE/(EE+EN)) - (CE/(CE+CN))
    
    SE_fixed[i,j] = sqrt(EE*EN/(EE+EN)^3 + CE*CN/(CE+CN)^3)
  }
}

Do the meta-analysis: 

# Meta-analysis

MA_all = MA_exclude = rep(NA, nsim)


for(i in 1:nsim){
  
  if(pval[i] <= 0.001){
    
    MA_all[i] = weighted.mean(c(RD_interim[i], RD_fixed[i,]),
                           c(1/SE_interim[i]^2, 1/SE_fixed[i,]^2))
    
    MA_exclude[i] = weighted.mean(c(RD_fixed[i,]), c(1/SE_fixed[i,]^2))
    
  } else {
    
    MA_all[i] = MA_exclude[i] = weighted.mean(c(RD_final[i], RD_fixed[i,]),
                           c(1/SE_final[i]^2, 1/SE_fixed[i,]^2))
    
  }
  
}

Meta-analysis results: 

#>                           nsims Mean no. studies Bias in RD Mean SE of RD
#> All trials                10000             5.00  -0.000127       0.00697
#> Exc. stopped early trials 10000             4.79   0.000561       0.00694