Running 'test.kmi.R' [167s]
Running the tests in 'tests/test.kmi.R' failed.
Complete output:
  > require(survival)
  Loading required package: survival
  > require(kmi)
  Loading required package: kmi
  > 
  > ### test 1
  > 
  > id <- 1:10
  > time <- 1:10
  > ev <- c(0, 0, 0, rep(1, 3), rep(2, 4))
  > data <- data.frame(id, time, ev)
  > 
  > aa <- kmi(Surv(time, ev != 0) ~ 1, data, etype = ev)
  > 
  > bb <- aa$imputed.data
  > test1 <- lapply(bb, function(x) {
  +     x[aa$original.data$id %in% 7:10, ]
  + })
  > test1 <- do.call(rbind, test1)
  > 
  > all(test1[, 1] == 11)
  [1] TRUE
  > all(test1[, 2] == 0)
  [1] TRUE
  > 
  > aa <- kmi(Surv(time, ev != 0) ~ 1, data, etype = ev, failcode = 2)
  > 
  > bb <- aa$imputed.data
  > test1 <- lapply(bb, function(x) {
  +     x[aa$original.data$id %in% 4:6, ]
  + })
  > test1 <- do.call(rbind, test1)
  > 
  > all(test1[, 1] == 11)
  [1] TRUE
  > all(test1[, 2] == 0)
  [1] TRUE
  > 
  > aa <- kmi(Surv(time, ev != 0) ~ 1, data, etype = ev, epsilon = 2)
  > 
  > bb <- aa$imputed.data
  > test1 <- lapply(bb, function(x) {
  +     x[aa$original.data$id %in% 7:10, ]
  + })
  > test1 <- do.call(rbind, test1)
  > 
  > all(test1[, 1] == 12)
  [1] TRUE
  > all(test1[, 2] == 0)
  [1] TRUE
  > 
  > aa <- kmi(Surv(time, ev != 0) ~ 1, data, etype = ev, nimp = 13)
  > length(aa$imputed.data) == 13
  [1] TRUE
  > 
  > 
  > ### test 2
  > 
  > set.seed(198740)
  > time <- rexp(100)
  > ev <- sample(c(0, 1, 2), 100, replace = TRUE)
  > cov <- rbinom(100, 1, 0.5)
  > dd <- data.frame(time, ev, cov)
  > 
  > ## add a test when etype is a factor
  > dd$status <- factor(ifelse(dd$ev == 0, "cens", ifelse(dd$ev == 1, "rel", "dc")))
  > 
  > set.seed(1440293)
  > dat.kmi <- kmi(Surv(time, ev != 0) ~ 1, dd, etype = ev, nimp = 5)
  > set.seed(1440293)
  > dat.kmi.fact <- kmi(Surv(time, status != "cens") ~ 1, dd, etype = status,
  +                     nimp = 5, failcode = "rel")
  > set.seed(1440293)
  > dat.kmi.mixed1 <- kmi(Surv(time, ev != 0) ~ 1, dd, etype = status,
  +                       nimp = 5, failcode = "rel")
  > set.seed(1440293)
  > dat.kmi.mixed2 <- kmi(Surv(time, status != "cens") ~ 1, dd, etype = ev,
  +                       nimp = 5, failcode = 1)
  > 
  > 
  > fit.kmi <- cox.kmi(Surv(time, ev == 1) ~ cov, dat.kmi)
  > fit.kmi.fact <- cox.kmi(Surv(time, status == "rel") ~ cov, dat.kmi)
  > fit.kmi.mixed1 <- cox.kmi(Surv(time, status == "rel") ~ cov, dat.kmi.mixed1)
  > fit.kmi.mixed2 <- cox.kmi(Surv(time, ev == 1) ~ cov, dat.kmi.mixed2)
  > 
  > all.equal(coef(fit.kmi), coef(fit.kmi.fact))
  [1] TRUE
  > all.equal(coef(fit.kmi), coef(fit.kmi.mixed1))
  [1] TRUE
  > all.equal(coef(fit.kmi), coef(fit.kmi.mixed2))
  [1] TRUE
  > all.equal(fit.kmi$variance, fit.kmi.fact$variance)
  [1] TRUE
  > all.equal(fit.kmi$variance, fit.kmi.mixed1$variance)
  [1] TRUE
  > all.equal(fit.kmi$variance, fit.kmi.mixed2$variance)
  [1] TRUE
  > 
  > ## fit.kmi
  > 
  > ## summary(fit.kmi)
  > 
  > ## fit.kmi.fact
  > 
  > ## summary(fit.kmi.fact)
  > 
  > ## avec bootstrap
  > 
  > set.seed(867988)
  > dat.kmib <- kmi(Surv(time, ev != 0) ~ 1, dd, etype = ev, nimp = 5,
  +                 boot = TRUE, nboot = 5)
  > set.seed(867988)
  > dat.kmib.fact <- kmi(Surv(time, ev != 0) ~ 1, dd, etype = status, nimp = 5,
  +                      boot = TRUE, nboot = 5, failcode = "rel")
  > 
  > fit.kmib <- cox.kmi(Surv(time, ev == 1) ~ cov, dat.kmib)
  > fit.kmib.fact <- cox.kmi(Surv(time, status == "rel") ~ cov, dat.kmib.fact)
  > 
  > all.equal(coef(fit.kmib), coef(fit.kmib.fact))
  [1] TRUE
  > 
  > fit.kmib
  Call:
  cox.kmi(formula = Surv(time, ev == 1) ~ cov, imp.data = dat.kmib)
  
  *****************
  Pooled estimates:
  *****************
             coef exp(coef)  se(coef)          t         p
  cov -0.07328941 0.9293318 0.3126316 -0.2344274 0.8146571
  
  > fit.kmib.fact
  Call:
  cox.kmi(formula = Surv(time, status == "rel") ~ cov, imp.data = dat.kmib.fact)
  
  *****************
  Pooled estimates:
  *****************
             coef exp(coef)  se(coef)          t         p
  cov -0.07328941 0.9293318 0.3126316 -0.2344274 0.8146571
  
  > 
  > summary(fit.kmib)
  Call:
  cox.kmi(formula = Surv(time, ev == 1) ~ cov, imp.data = dat.kmib)
  
  
  *****************
  Pooled estimates:
  *****************
          coef exp(coef) se(coef)      t Pr(>|t|)
  cov -0.07329   0.92933  0.31263 -0.234    0.815
  
      exp(coef) exp(-coef) lower .95 upper .95
  cov    0.9293      1.076    0.5035     1.715
  
  > 
  > 
  > ### test 4
  > 
  > data(icu.pneu)
  > 
  > set.seed(1313)
  > dat <- kmi(Surv(start, stop, status) ~ 1, data = icu.pneu,
  +            etype = event, id= id, failcode = 2, nimp = 5)
  > 
  > ### add a factor for testing purposes
  > icu.pneu$ev.fact <- factor(ifelse(icu.pneu$event == 3, "disch", "death"))
  > 
  > icu.pneu$ev <- icu.pneu$event
  > icu.pneu$ev[icu.pneu$status == 0] <- 0
  > 
  > set.seed(1313)
  > dat2 <- kmi(Surv(start, stop, ev != 0) ~ 1, data = icu.pneu,
  +            etype = ev, id= id, failcode = 2, nimp = 5)
  > 
  > set.seed(1313)
  > dat3 <- kmi(Surv(start, stop, status) ~ 1, data = icu.pneu,
  +             etype = ev.fact, id = id, failcode = "death", nimp = 5)
  > 
  > a <- logical(5)
  > for (i in 1:5) a[i] <- all.equal(dat$imputed.data[[i]][, 1], dat2$imputed.data[[i]][, 1])
  > a
  [1] TRUE TRUE TRUE TRUE TRUE
  > 
  > fit.kmi <- cox.kmi(Surv(start, stop, event == 2) ~ pneu, dat)
  > 
  > fit.kmi2 <- cox.kmi(Surv(start, stop, ev == 2) ~ pneu, dat2)
  > 
  > fit.kmi3 <- cox.kmi(Surv(start, stop, ev.fact == "death") ~ pneu, dat3)
  > 
  > all.equal(fit.kmi$coefficients, fit.kmi2$coefficients)
  [1] TRUE
  > all.equal(coef(fit.kmi), coef(fit.kmi3))
  [1] TRUE
  > all.equal(fit.kmi$variance, fit.kmi2$variance)
  [1] TRUE
  > 
  > fit.kmi
  Call:
  cox.kmi(formula = Surv(start, stop, event == 2) ~ pneu, imp.data = dat)
  
  *****************
  Pooled estimates:
  *****************
           coef exp(coef)  se(coef)       t            p
  pneu1 1.10382  3.015664 0.2398608 4.60192 4.186147e-06
  
  > 
  > fit.kmi2
  Call:
  cox.kmi(formula = Surv(start, stop, ev == 2) ~ pneu, imp.data = dat2)
  
  *****************
  Pooled estimates:
  *****************
           coef exp(coef)  se(coef)       t            p
  pneu1 1.10382  3.015664 0.2398608 4.60192 4.186147e-06
  
  > 
  > fit.kmi3
  Call:
  cox.kmi(formula = Surv(start, stop, ev.fact == "death") ~ pneu, 
      imp.data = dat3)
  
  *****************
  Pooled estimates:
  *****************
           coef exp(coef)  se(coef)       t            p
  pneu1 1.10382  3.015664 0.2398608 4.60192 4.186147e-06
  
  > 
  > ## avec bootstrap
  > 
  > set.seed(598085)
  > dat <- kmi(Surv(start, stop, status) ~ 1, data = icu.pneu,
  +            etype = event, id= id, failcode = 2, nimp = 5,
  +            boot = TRUE, nboot = 5)
  > 
  > set.seed(598085)
  > dat2 <- kmi(Surv(start, stop, ev != 0) ~ 1, data = icu.pneu,
  +             etype = ev, id= id, failcode = 2, nimp = 5,
  +             boot = TRUE, nboot = 5)
  > 
  > set.seed(598085)
  > dat3 <- kmi(Surv(start, stop, status) ~ 1, data = icu.pneu,
  +             etype = ev.fact, id = id, failcode = "death", nimp = 5,
  +             boot = TRUE, nboot = 5)
  > 
  > 
  > a <- logical(5)
  > for (i in 1:5) a[i] <- all.equal(dat$imputed.data[[i]][, 1], dat2$imputed.data[[i]][, 1])
  > a
  [1] TRUE TRUE TRUE TRUE TRUE
  > 
  > fit.kmi <- cox.kmi(Surv(start, stop, event == 2) ~ pneu, dat)
  > 
  > fit.kmi2 <- cox.kmi(Surv(start, stop, ev == 2) ~ pneu, dat2)
  > 
  > fit.kmi3 <- cox.kmi(Surv(start, stop, ev.fact == "death") ~ pneu, dat3)
  > 
  > all.equal(fit.kmi$coefficients, fit.kmi2$coefficients)
  [1] TRUE
  > all.equal(coef(fit.kmi), coef(fit.kmi3))
  [1] TRUE
  > all.equal(fit.kmi$variance, fit.kmi2$variance)
  [1] TRUE
  > 
  > fit.kmi
  Call:
  cox.kmi(formula = Surv(start, stop, event == 2) ~ pneu, imp.data = dat)
  
  *****************
  Pooled estimates:
  *****************
            coef exp(coef)  se(coef)        t            p
  pneu1 1.104148  3.016652 0.2399367 4.601829 4.188064e-06
  
  > 
  > fit.kmi2
  Call:
  cox.kmi(formula = Surv(start, stop, ev == 2) ~ pneu, imp.data = dat2)
  
  *****************
  Pooled estimates:
  *****************
            coef exp(coef)  se(coef)        t            p
  pneu1 1.104148  3.016652 0.2399367 4.601829 4.188064e-06
  
  > 
  > fit.kmi3
  Call:
  cox.kmi(formula = Surv(start, stop, ev.fact == "death") ~ pneu, 
      imp.data = dat3)
  
  *****************
  Pooled estimates:
  *****************
            coef exp(coef)  se(coef)        t            p
  pneu1 1.104148  3.016652 0.2399367 4.601829 4.188064e-06
  
  > 
  > ### with covariates
  > ## classic
  > set.seed(1)
  > dd$juhu <- rnorm(nrow(dd))
  > 
  > set.seed(78223)
  > imp.dd <- kmi(Surv(time, ev != 0) ~ juhu, dd,
  +               etype = ev, nimp = 5)
  > set.seed(44889)
  > imp.ddb <- kmi(Surv(time, ev != 0) ~ juhu, dd, nimp = 5,
  +               etype = ev, boot = TRUE, nboot = 5)
  > 
  > summary(cox.kmi(Surv(time, ev == 1) ~ cov, imp.dd))
  Call:
  cox.kmi(formula = Surv(time, ev == 1) ~ cov, imp.data = imp.dd)
  
  
  *****************
  Pooled estimates:
  *****************
         coef exp(coef) se(coef)      t Pr(>|t|)
  cov -0.0485    0.9527   0.3155 -0.154    0.878
  
      exp(coef) exp(-coef) lower .95 upper .95
  cov    0.9527       1.05    0.5132     1.768
  
  > summary(cox.kmi(Surv(time, ev == 1) ~ cov, imp.ddb))
  Call:
  cox.kmi(formula = Surv(time, ev == 1) ~ cov, imp.data = imp.ddb)
  
  
  *****************
  Pooled estimates:
  *****************
          coef exp(coef) se(coef)      t Pr(>|t|)
  cov -0.08706   0.91662  0.31606 -0.275    0.783
  
      exp(coef) exp(-coef) lower .95 upper .95
  cov    0.9166      1.091    0.4932     1.704
  
  > 
  > ## time-dependent covariates
  > set.seed(9763)
  > imp.dat.c <- kmi(Surv(start, stop, status) ~ age + sex,
  +                  data = icu.pneu, etype = event, id = id,
  +                  failcode = 2, nimp = 5)
  > 
  > set.seed(19832)
  > imp.dat.cb <- kmi(Surv(start, stop, status) ~ age + sex,
  +                   data = icu.pneu, etype = event, id = id,
  +                   failcode = 2, nimp = 5, boot = TRUE,
  +                   nboot = 5)
  > 
  > summary(cox.kmi(Surv(start, stop, event == 2) ~ pneu,
  +                 imp.dat.c))
  Call:
  cox.kmi(formula = Surv(start, stop, event == 2) ~ pneu, imp.data = imp.dat.c)
  
  
  *****************
  Pooled estimates:
  *****************
         coef exp(coef) se(coef)     t Pr(>|t|)    
  pneu1 1.117     3.055    0.240 4.654 3.25e-06 ***
  ---
  Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
  
        exp(coef) exp(-coef) lower .95 upper .95
  pneu1     3.055     0.3273     1.909      4.89
  
  > summary(cox.kmi(Surv(start, stop, event == 2) ~ pneu,
  +                 imp.dat.cb))
  Call:
  cox.kmi(formula = Surv(start, stop, event == 2) ~ pneu, imp.data = imp.dat.cb)
  
  
  *****************
  Pooled estimates:
  *****************
          coef exp(coef) se(coef)     t Pr(>|t|)    
  pneu1 1.1077    3.0274   0.2403 4.609 4.04e-06 ***
  ---
  Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
  
        exp(coef) exp(-coef) lower .95 upper .95
  pneu1     3.027     0.3303      1.89     4.849
  
  > 
  > 
  > ## Test ordering of the data set
  > 
  > data(icu.pneu)
  > icu.pneu$outcome <- with(icu.pneu, status * event)
  > 
  > mydata1 <- icu.pneu[1:2, ]
  > mydata2 <- icu.pneu[3:1421, ]
  > 
  > ## page 221
  > set.seed(453)
  > imp.dat <- kmi(Surv(start, stop, outcome != 0) ~ 1,
  +     data = icu.pneu, etype = outcome,
  +     id = id, failcode = 2, nimp = 10)
  > 
  > ## page 222
  > kmi.sh.hap <- cox.kmi(Surv(start, stop, outcome == 2) ~ pneu, imp.dat)
  > summary(kmi.sh.hap)
  Call:
  cox.kmi(formula = Surv(start, stop, outcome == 2) ~ pneu, imp.data = imp.dat)
  
  
  *****************
  Pooled estimates:
  *****************
          coef exp(coef) se(coef)     t Pr(>|t|)    
  pneu1 1.1057    3.0213   0.2399 4.608 4.06e-06 ***
  ---
  Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
  
        exp(coef) exp(-coef) lower .95 upper .95
  pneu1     3.021      0.331     1.888     4.835
  
  > 
  > newicu.pneu <- as.data.frame(rbind( mydata2, mydata1))
  > tail(newicu.pneu)
            id start stop status event pneu adm.cens.exit      age sex outcome
  1418 2020485     0   21      1     3    0           445 69.85106   M       3
  1419 2020508     0   49      1     3    0           537 60.23257   M       3
  1420 2020634     0   35      1     3    0           555 59.52186   M       3
  1421 2020811     0    3      1     3    0           465 73.95806   F       3
  1         12     0    3      1     3    0           421 62.53307   F       3
  2         41     0    5      1     3    0           545 75.34153   F       3
  > 
  > set.seed(453)
  > imp.datnew <- kmi(Surv(start, stop, outcome != 0) ~ 1,
  +     data = newicu.pneu, etype = outcome,
  +     id = id, failcode = 2, nimp = 10)
  > 
  > # page 222
  > kmi.sh.hap <- cox.kmi(Surv(start, stop, outcome == 2) ~ pneu, imp.datnew)
  > summary(kmi.sh.hap)
  Call:
  cox.kmi(formula = Surv(start, stop, outcome == 2) ~ pneu, imp.data = imp.datnew)
  
  
  *****************
  Pooled estimates:
  *****************
          coef exp(coef) se(coef)     t Pr(>|t|)    
  pneu1 1.1057    3.0213   0.2399 4.608 4.06e-06 ***
  ---
  Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
  
        exp(coef) exp(-coef) lower .95 upper .95
  pneu1     3.021      0.331     1.888     4.835
  
  > 
  > proc.time()
     user  system elapsed 
     6.03    0.53    6.54