# library(smim)
library(survRM2)
#> Loading required package: survival
devtools::load_all()
#> ℹ Loading smim

Data Simulation

  • Pattern
    • 1: observed event time
    • 2: MAR imputation
    • 3: MNAR imputation
  • Delta
    • 1: no adjustment
    • larger than 1: hazard is larger than MAR
    • less than 1: hazard is smaller than MAR
n <- 1000
t.time <- rexp(n)   # event time
c.time <- runif(n, min = 0.1, max = 5)  # censoring time

time <- pmin(t.time, c.time)

status <- as.numeric(t.time < c.time)

pattern <- ifelse(status == 1, 1, ifelse(time < 0.5, 2, 3))
group <- rep(c(0,1), each = 500)
plot(survfit(Surv(time, status) ~ group), col = 1:2)

table(status, pattern, time < 0.5)
#> , ,  = FALSE
#> 
#>       pattern
#> status   1   2   3
#>      0   0   0 108
#>      1 431   0   0
#> 
#> , ,  = TRUE
#> 
#>       pattern
#> status   1   2   3
#>      0   0  61   0
#>      1 400   0   0

MAR analysis

tau <- 3

survRM2 results

rmst2(time, status, group, tau = tau)
#> 
#> The truncation time: tau = 3  was specified. 
#> 
#> Restricted Mean Survival Time (RMST) by arm 
#>               Est.   se lower .95 upper .95
#> RMST (arm=1) 0.936 0.04     0.858     1.014
#> RMST (arm=0) 0.956 0.04     0.876     1.035
#> 
#> 
#> Restricted Mean Time Lost (RMTL) by arm 
#>               Est.   se lower .95 upper .95
#> RMTL (arm=1) 2.064 0.04     1.986     2.142
#> RMTL (arm=0) 2.044 0.04     1.965     2.124
#> 
#> 
#> Between-group contrast 
#>                        Est. lower .95 upper .95     p
#> RMST (arm=1)-(arm=0) -0.019    -0.131     0.092 0.732
#> RMST (arm=1)/(arm=0)  0.980     0.871     1.102 0.732
#> RMTL (arm=1)/(arm=0)  1.010     0.956     1.066 0.732

RMST within group

delta   <- c(1,1,1)[pattern]  # the third number control delta adjustment for MNAR 
tmp <- rmst_delta(time, status, x = rep(1, length(time)), group, pattern, delta, tau, n_mi = 5, n_b = 100, seed = 123)
tmp$rmst
#>      group      rmst         sd      wb_sd
#> [1,]     0 0.9483854 0.03912381 0.04178793
#> [2,]     1 0.9377955 0.03804502 0.03616984

RMST between group

diff_rmst <- function(rmst, sd){
  diff <- diff(rmst)
  diff_sd <- sqrt( sum(sd^2) )
  p_val <- 2* (1 - pnorm( abs(diff/diff_sd) ))
  c(diff = diff, sd = diff_sd, p = p_val)
}

rbind( diff_rmst(tmp$rmst[,"rmst"], tmp$rmst[, "sd"]), 
       diff_rmst(tmp$rmst[,"rmst"], tmp$rmst[, "wb_sd"]))
#>             diff         sd         p
#> [1,] -0.01058991 0.05457194 0.8461335
#> [2,] -0.01058991 0.05526743 0.8480460