CRAN Package Check Results for Package riskRegression

Last updated on 2021-12-05 07:49:21 CET.

Flavor Version Tinstall Tcheck Ttotal Status Flags
r-devel-linux-x86_64-debian-clang 2021.10.10 246.42 566.68 813.10 OK
r-devel-linux-x86_64-debian-gcc 2021.10.10 184.82 402.87 587.69 OK
r-devel-linux-x86_64-fedora-clang 2021.10.10 1065.71 NOTE
r-devel-linux-x86_64-fedora-gcc 2021.10.10 999.00 OK
r-devel-windows-x86_64-new-UL 2021.10.10 331.00 642.00 973.00 ERROR
r-devel-windows-x86_64-new-TK 2021.10.10 OK
r-devel-windows-x86_64-old 2021.10.10 337.00 503.00 840.00 OK
r-patched-linux-x86_64 2021.10.10 227.45 508.18 735.63 OK
r-release-linux-x86_64 2021.10.10 217.71 520.09 737.80 OK
r-release-macos-arm64 2021.10.10 NOTE
r-release-macos-x86_64 2021.10.10 NOTE
r-release-windows-ix86+x86_64 2021.10.10 600.00 1103.00 1703.00 OK
r-oldrel-macos-x86_64 2021.10.10 NOTE
r-oldrel-windows-ix86+x86_64 2021.10.10 621.00 1020.00 1641.00 OK

Check Details

Version: 2021.10.10
Check: installed package size
Result: NOTE
     installed size is 9.7Mb
     sub-directories of 1Mb or more:
     libs 8.3Mb
Flavors: r-devel-linux-x86_64-fedora-clang, r-release-macos-arm64, r-release-macos-x86_64, r-oldrel-macos-x86_64

Version: 2021.10.10
Check: tests
Result: ERROR
     Running 'test-AUC.R' [8s]
     Running 'test-BinomialRegression.R' [11s]
     Running 'test-BrierScore.R' [10s]
     Running 'test-CauseSpecificCoxRegression.R' [10s]
     Running 'test-FineGrayRegression.R' [8s]
     Running 'test-ate.R' [47s]
     Running 'test-boxplotScore.R' [9s]
     Running 'test-iidCox.R' [12s]
     Running 'test-phreg.R' [14s]
     Running 'test-predictCSC.R' [17s]
     Running 'test-predictCox.R' [24s]
     Running 'test-predictCoxPL.R' [10s]
     Running 'test-predictRisk-TD.R' [10s]
     Running 'test-predictRisk.R' [11s]
     Running 'test-score-bootstrap.R' [30s]
     Running 'test-score.R' [14s]
     Running 'test-summary-score.R' [16s]
    Running the tests in 'tests/test-predictCSC.R' failed.
    Complete output:
     > ## test-predictCSC.R ---
     > #----------------------------------------------------------------------
     > ## author: Brice Ozenne
     > ## created: maj 18 2017 (09:23)
     > ## Version:
     > ## last-updated: okt 7 2021 (20:02)
     > ## By: Brice Ozenne
     > ## Update #: 309
     > #----------------------------------------------------------------------
     > ##
     > ### Commentary:
     > ##
     > ### Change Log:
     > #----------------------------------------------------------------------
     > ##
     > ### Code:
     >
     > ## * Settings
     > library(riskRegression)
     riskRegression version 2021.10.10
     > library(testthat)
     > library(data.table)
     > library(mstate)
     Loading required package: survival
     > library(prodlim)
     > library(rms)
     Loading required package: Hmisc
     Loading required package: lattice
     Loading required package: Formula
     Loading required package: ggplot2
    
     Attaching package: 'Hmisc'
    
     The following object is masked from 'package:testthat':
    
     describe
    
     The following objects are masked from 'package:base':
    
     format.pval, units
    
     Loading required package: SparseM
    
     Attaching package: 'SparseM'
    
     The following object is masked from 'package:base':
    
     backsolve
    
     > tmat <- trans.comprisk(2, names = c("0", "1", "2"))
     > library(survival)
     > library(timereg); nsim.band <- 100
     >
     > ## * [predictCSC] predict.CSC vs. manual calculation and mstate on small examples
     > cat("[predictCSC] predict.CSC vs. manual calculation on small examples \n")
     [predictCSC] predict.CSC vs. manual calculation on small examples
     > ## ** Data
     > ## simulatenous events
     > df1 <- data.frame(time = rep(1:10,2),
     + event = c(rep(1,10),rep(2,10)),
     + X1 = 0:1
     + )
     > df1$event1 <- as.numeric(df1$event == 1)
     > df1$event2 <- as.numeric(df1$event == 2)
     >
     > set.seed(11)
     > dfS <- rbind(cbind(df1, grp = 1, X2 = rbinom(20,1,.4)),
     + cbind(rbind(df1,df1), grp = 2, X2 = rbinom(20,1,.4)),
     + cbind(df1, grp = 3, X2 = rbinom(20,1,.4))
     + )
     >
     > ## distinct events
     > df2 <- data.frame(time = c(1:10,-0.1+1:10),
     + event = c(rep(1,10),rep(2,10)),
     + X1 = 0:1
     + )
     > df2$event1 <- as.numeric(df2$event == 1)
     > df2$event2 <- as.numeric(df2$event == 2)
     >
     > ## ** No risk factor no strata
     > test_that("[predictCSC] (no strata, data=df1): compare to manual estimation",{
     +
     + CSC.exp <- CSC(Hist(time,event) ~ 1, data = df1)
     + pred.exp <- predict(CSC.exp, times = 1:10, newdata = data.frame(NA), cause = 1, product.limit = FALSE)
     + pred.exp2 <- predict(CSC.exp, times = 1:10, newdata = data.frame(NA), cause = 1, product.limit = TRUE)
     + CSC.prodlim <- CSC(Hist(time,event) ~ 1, data = df1, surv.type = "survival", method = "breslow")
     + pred.prodlim <- predict(CSC.prodlim, times = 1:10, newdata = data.frame(NA), cause = 1, product.limit = TRUE)
     + expect_equal(pred.prodlim,pred.exp2)
     + CSC.prodlimE <- CSC(Hist(time,event) ~ 1, data = df1, surv.type = "survival", method = "efron")
     + pred.prodlimE <- predict(CSC.prodlimE, times = 1:10, newdata = data.frame(NA), cause = 1, product.limit = TRUE)
     +
     + ## test baseline hazard
     + lambda1 <- lambda2 <- 1/seq(20,2,by = -2)
     + lambda2E <- predictCox(CSC.prodlimE$models[[2]], type = "hazard")$hazard
     + expect_equal(predictCox(CSC.exp$models[[1]], type = "hazard")$hazard,
     + lambda1)
     + expect_equal(predictCox(CSC.exp$models[[2]], type = "hazard")$hazard,
     + lambda2)
     + expect_equal(predictCox(CSC.prodlim$models[[1]], type = "hazard")$hazard,
     + lambda1)
     + expect_equal(predictCox(CSC.prodlim$models[[2]], type = "hazard")$hazard,
     + lambda1+lambda2)
     + expect_equal(basehaz(CSC.prodlimE$models[[2]])$hazard,
     + cumsum(lambda2E)
     + )
     +
     + ## test absolute risk
     + survival <- c(1,exp(cumsum(-lambda1-lambda2)))[1:10]
     + expect_equal(as.double(pred.exp$absRisk),
     + cumsum(lambda1*survival))
     + survival <- c(1,cumprod(1-(lambda1+lambda2)))[1:10]
     + expect_equal(as.double(pred.prodlim$absRisk),
     + cumsum(lambda1*survival))
     + survival <- c(1,cumprod(1-lambda2E))[1:10]
     + expect_equal(as.double(pred.prodlimE$absRisk),
     + cumsum(lambda1*survival))
     + })
     Test passed 🎉
     >
     > test_that("[predictCSC] (no strata, data=df1): compare to manual estimation",{
     + CSC.exp <- CSC(Hist(time,event) ~ 1, data = df2)
     + pred.exp <- predict(CSC.exp, times = 1:10, newdata = data.frame(NA), cause = 1, product.limit = FALSE)
     + pred.exp2 <- predict(CSC.exp, times = 1:10, newdata = data.frame(NA), cause = 1, product.limit = TRUE)
     + CSC.prodlim <- CSC(Hist(time,event) ~ 1, data = df2, surv.type = "survival", method = "breslow")
     + pred.prodlim <- predict(CSC.prodlim, times = 1:10, newdata = data.frame(NA), cause = 1, product.limit = TRUE)
     + expect_equal(pred.exp2, pred.prodlim)
     + CSC.prodlimE <- CSC(Hist(time,event) ~ 1, data = df2, surv.type = "survival", method = "efron")
     + pred.prodlimE <- predict(CSC.prodlimE, times = 1:10, newdata = data.frame(NA), cause = 1, product.limit = TRUE)
     + ## test baseline hazard
     + lambda1 <- 1/seq(19,1,by = -2)
     + lambda2 <- 1/seq(20,2,by = -2)
     + lambda2E <- as.data.table(predictCox(CSC.prodlimE$models[[2]], type = "hazard"))
     + expect_equal(predictCox(CSC.exp$models[[1]], times = 1:10, type = "hazard")$hazard,
     + lambda1)
     + expect_equal(predictCox(CSC.exp$models[[2]], times = -0.1 + 1:10, type = "hazard")$hazard,
     + lambda2)
     + expect_equal(predictCox(CSC.prodlim$models[[1]], times = 1:10, type = "hazard")$hazard,
     + lambda1)
     + expect_equal(predictCox(CSC.prodlim$models[[2]], type = "hazard")$hazard,
     + as.double(rbind(lambda2,lambda1)))
     + expect_equal(basehaz(CSC.prodlimE$models[[2]])$hazard,
     + cumsum(lambda2E$hazard)
     + )
     + ## test absolute risk
     + vec.lambda <- cumsum(as.double(rbind(lambda2,lambda1)))
     + survival <- exp(-matrix(vec.lambda, ncol = 2, byrow = TRUE)[,1])
     + expect_equal(as.double(pred.exp$absRisk),
     + cumsum(lambda1*survival))
     + vec.lambda <- as.double(rbind(lambda2,lambda1))
     + survival <- matrix(cumprod(1-vec.lambda), ncol = 2, byrow = TRUE)[,1]
     + expect_equal(as.double(pred.prodlim$absRisk),
     + cumsum(lambda1*survival))
     + vec.lambda <-lambda2E$hazard
     + survival <- matrix(cumprod(1-vec.lambda), ncol = 2, byrow = TRUE)[,1]
     + expect_equal(as.double(pred.prodlimE$absRisk),
     + cumsum(lambda1*survival))
     + })
     Test passed 😸
     >
     > test_that("[predictCSC]: compare to mstate",{
     + for(iData in 1:2){# iData <- 1
     + df <- switch(as.character(iData),
     + "1"=df1,
     + "2"=df2)
     + ## fit CSC
     + CSC.exp <- CSC(Hist(time,event) ~ 1, data = df)
     + ## fit mstate
     + dL <- msprep(time = c(NA, "time", "time"),
     + status = c(NA,"event1", "event2"),
     + data = df, keep = c("X1"),
     + trans = tmat)
     + dL.exp <- expand.covs(dL, c("X1"))
     + e.coxph <- coxph(Surv(time, status) ~ strata(trans),
     + data = dL.exp)
     + ## compare
     + newdata.L <- data.frame(trans = c(1, 2), strata = c(1, 2))
     + newdata <- data.frame(NA)
     + pred.msfit <- msfit(e.coxph, newdata = newdata.L, trans = tmat)
     + suppressWarnings(
     + pred.probtrans <- probtrans(pred.msfit,0)[[1]]
     + )
     + pred.RR <- predict(CSC.exp,
     + times = pred.probtrans[,"time"],
     + newdata = newdata,
     + cause = 1,
     + se = FALSE,
     + product.limit = TRUE)
     + expect_equal(pred.probtrans[,"pstate2"],
     + as.double(pred.RR$absRisk)
     + )
     + }
     + })
     Test passed 😀
     >
     >
     > test_that("[predict.CSC] check absolute risks add up to one",{
     + for(iData in 1:2){
     + df <- switch(as.character(iData),
     + "1"=df1,
     + "2"=df2)
     + ## fit CSC
     + CSC.exp <- CSC(Hist(time,event) ~ 1, data = df, method = "breslow")
     + ## compute all transition probabilites
     + seqTimes <- sort(unique(df$time))
     + nTimes <- length(seqTimes)
     + newdata <- data.frame(NA)
     + hazAll <- predictCox(CSC.exp$models[[1]], times = seqTimes, newdata = newdata, type = "hazard")$hazard+predictCox(CSC.exp$models[[2]], times = seqTimes, newdata = newdata, type = "hazard")$hazard
     + surv.prodlim <- cumprod(1-as.double(hazAll))
     + haz1.prodlim <- predict(CSC.exp, times = seqTimes, newdata = newdata, cause = 1, product.limit = TRUE)
     + haz2.prodlim <- predict(CSC.exp, times = seqTimes, newdata = newdata, cause = 2, product.limit = TRUE)
     + expect_equal(surv.prodlim+as.double(haz1.prodlim$absRisk)+as.double(haz2.prodlim$absRisk),rep(1,nTimes))
     + }
     + })
     Test passed 🎉
     >
     > ## ** No risk factor strata
     > test_that("[predict.CSC] (strata): compare to manual estimation",{
     + ## as.data.table(dfS)[,.N, by = "grp"]
     + CSC.exp <- CSC(Hist(time,event) ~ strata(grp), data = dfS, method = "breslow")
     + pred.exp <- predict(CSC.exp, times = 1:10, newdata = data.frame(grp = 1:3), cause = 1, product.limit = FALSE)
     + pred.exp2 <- predict(CSC.exp, times = 1:10, newdata = data.frame(grp = 1:3), cause = 1, product.limit = TRUE)
     +
     + CSC.prodlim <- CSC(Hist(time,event) ~ strata(grp), data = dfS, surv.type = "survival", method = "breslow")
     + pred.prodlim <- predict(CSC.prodlim, times = 1:10, newdata = data.frame(grp = 1:3), cause = 1, product.limit = TRUE)
     + expect_equal(pred.exp2,pred.prodlim)
     +
     + CSC.prodlimE <- CSC(Hist(time,event) ~ strata(grp), data = dfS, surv.type = "survival", method = "efron")
     + pred.prodlimE <- predict(CSC.prodlimE, times = 1:10, newdata = data.frame(grp = 1:3), cause = 1, product.limit = TRUE)
     +
     + ## baseline
     + baseline1 <- as.data.table(predictCox(CSC.exp$models[[1]], type = "hazard"))
     + lambda1 <- baseline1[strata=="grp=1",hazard]
     + expect_equal(lambda1, 1/seq(20,2,by = -2))
     +
     + baseline2 <- as.data.table(predictCox(CSC.exp$models[[1]], type = "hazard"))
     + lambda2 <- baseline2[strata=="grp=1",hazard]
     + expect_equal(lambda2, 1/seq(20,2,by = -2))
     +
     + baseline1E <- as.data.table(predictCox(CSC.prodlimE$models[[1]], type = "hazard"))
     + lambda1E.1 <- baseline1E[strata=="grp=1",hazard]
     + lambda1E.2 <- baseline1E[strata=="grp=2",hazard]
     + expect_equal(lambda1E.1, 1/seq(20,2,by = -2))
     +
     + baseline2E <- as.data.table(predictCox(CSC.prodlimE$models[[2]], type = "hazard"))
     + lambda2E.1 <- baseline2E[strata=="grp=1",hazard]
     + lambda2E.2 <- baseline2E[strata=="grp=2",hazard]
     + ##
     +
     + ## test absolute risk
     + survival <- c(1,exp(cumsum(-lambda1-lambda2)))[1:10]
     + expect_equal(as.double(pred.exp$absRisk[1,]),
     + cumsum(lambda1*survival))
     + expect_equal(as.double(pred.exp$absRisk[1,]),
     + as.double(pred.exp$absRisk[2,]))
     + expect_equal(as.double(pred.exp$absRisk[1,]),
     + as.double(pred.exp$absRisk[3,]))
     +
     + survival <- c(1,cumprod(1-(lambda1+lambda2)))[1:10]
     + expect_equal(as.double(pred.prodlim$absRisk[1,]),
     + cumsum(lambda1*survival))
     + expect_equal(as.double(pred.prodlim$absRisk[1,]),
     + as.double(pred.prodlim$absRisk[2,]))
     + expect_equal(as.double(pred.prodlim$absRisk[1,]),
     + as.double(pred.prodlim$absRisk[3,]))
     +
     + survival <- c(1,cumprod(1-lambda2E.1))[1:10]
     + expect_equal(as.double(pred.prodlimE$absRisk[pred.prodlimE$strata=="grp=1"]),
     + cumsum(lambda1E.1*survival))
     + survival <- c(1,cumprod(1-lambda2E.2))[1:10]
     + expect_equal(as.double(pred.prodlimE$absRisk[pred.prodlimE$strata=="grp=2"]),
     + cumsum(lambda1E.2*survival))
     + expect_equal(as.double(pred.prodlimE$absRisk[1,]),
     + as.double(pred.prodlimE$absRisk[3,]))
     + })
     Test passed 🌈
     >
     > ## ** risk factor no strata
     > test_that("[predict.CSC] (risk factor, strata): compare to manual estimation",{
     + CSC.exp <- CSC(Hist(time,event) ~ X1, data = df1)
     + CSC.prodlim <- CSC(Hist(time,event) ~ X1, data = df1, surv.type = "survival", method = "breslow")
     + CSC.prodlimE <- CSC(Hist(time,event) ~ X1, data = df1, surv.type = "survival", method = "efron")
     +
     + calcCIF <- function(CSC, X){
     + lambda01 <- as.data.table(predictCox(CSC$models[[1]], centered = FALSE, type = "hazard"))
     + lambda02 <- as.data.table(predictCox(CSC$models[[2]], centered = FALSE, type = "hazard"))
     + eXb.1 <- as.double(exp(X%*%coef(CSC)[[1]]))
     + eXb.2 <- as.double(exp(X%*%coef(CSC)[[2]]))
     + if(CSC$surv.type == "survival"){
     + survival <- cumprod(1-eXb.2*lambda02$hazard)
     + }else{
     + Lambda01 <- cumsum(lambda01$hazard)
     + Lambda02 <- cumsum(lambda02$hazard)
     + survival <- exp(-eXb.1*Lambda01-eXb.2*Lambda02)
     + }
     + survival <- c(1,survival[-length(survival)])
     + absRisk <- cumsum(eXb.1*lambda01$hazard*survival)
     + return(absRisk)
     + }
     + for(iX1 in 0:1){
     + newdata <- data.frame(X1=iX1)
     + pred.exp <- predict(CSC.exp, times = 1:10, newdata = newdata, cause = 1, product.limit = FALSE)
     + expect_equal(as.double(pred.exp$absRisk),
     + calcCIF(CSC.exp, X = as.matrix(newdata))
     + )
     + pred.prodlim <- predict(CSC.prodlim, times = 1:10, newdata = newdata, cause = 1, product.limit = TRUE)
     + expect_equal(as.double(pred.prodlim$absRisk),
     + calcCIF(CSC.prodlim, X = as.matrix(newdata))
     + )
     + pred.prodlimE <- predict(CSC.prodlimE, times = 1:10, newdata = newdata, cause = 1, product.limit = TRUE)
     + expect_equal(as.double(pred.prodlimE$absRisk),
     + calcCIF(CSC.prodlimE, X = as.matrix(newdata))
     + )
     + }
     + })
     Test passed 🎉
     >
     > test_that("[predict.CSC] - compare to manual estimation",{
     + CSC.exp <- CSC(Hist(time,event) ~ X1, data = df2, method = "breslow")
     +
     + CSC.prodlim <- CSC(Hist(time,event) ~ X1, data = df2, surv.type = "survival", method = "breslow")
     +
     + CSC.prodlimE <- CSC(Hist(time,event) ~ X1, data = df2, surv.type = "survival", method = "efron")
     +
     + calcCIF <- function(CSC, X){
     +
     + lambda01 <- as.data.table(predictCox(CSC$models[[1]], centered = FALSE, type = "hazard"))
     + lambda02 <- as.data.table(predictCox(CSC$models[[2]], centered = FALSE, type = "hazard"))
     + eXb.1 <- as.double(exp(X%*%coef(CSC)[[1]]))
     + eXb.2 <- as.double(exp(X%*%coef(CSC)[[2]]))
     +
     + if(CSC$surv.type == "survival"){
     + survival <- cumprod(1-eXb.2*lambda02$hazard)
     + }else{
     + Lambda01 <- cumsum(lambda01$hazard)
     + Lambda02 <- cumsum(lambda02$hazard)
     + survival <- exp(-eXb.1*Lambda01-eXb.2*Lambda02)
     + }
     +
     + survival <- c(1,survival[-length(survival)])
     + absRisk <- cumsum(eXb.1*lambda01$hazard*survival)
     + return(absRisk)
     +
     + }
     +
     + for(iX1 in 0:1){ # iX1 <- 0
     + newdata <- data.frame(X1=iX1)
     + pred.exp <- predict(CSC.exp, times = sort(unique(df2$time)), newdata = newdata, cause = 1, product.limit = FALSE)
     + expect_equal(as.double(pred.exp$absRisk),
     + calcCIF(CSC.exp, X = as.matrix(newdata))
     + )
     +
     + pred.prodlim <- predict(CSC.prodlim, times = sort(unique(df2$time)), newdata = newdata, cause = 1, product.limit = TRUE)
     + expect_equal(as.double(pred.prodlim$absRisk),
     + calcCIF(CSC.prodlim, X = as.matrix(newdata))
     + )
     +
     + pred.prodlimE <- predict(CSC.prodlimE, times = sort(unique(df2$time)), newdata = newdata, cause = 1, product.limit = TRUE)
     + expect_equal(as.double(pred.prodlimE$absRisk),
     + calcCIF(CSC.prodlimE, X = as.matrix(newdata))
     + )
     + }
     + })
     Test passed 🥇
     >
     >
     > test_that("[predict.CSC]: compare to mstate",{
     + for(iData in 1:2){
     + df <- switch(as.character(iData),
     + "1"=df1,
     + "2"=df2)
     + ## fit CSC
     + CSC.exp <- CSC(Hist(time,event) ~ X1, data = df, method = "breslow")
     + ## fit mstate
     + df$event1 <- as.numeric(df$event == 1)
     + df$event2 <- as.numeric(df$event == 2)
     + dL <- msprep(time = c(NA, "time", "time"),
     + status = c(NA,"event1", "event2"),
     + data = df, keep = c("X1"),
     + trans = tmat)
     + dL.exp <- expand.covs(dL, c("X1"))
     + e.coxph <- coxph(Surv(time, status) ~ X1.1 + X1.2 + strata(trans),
     + data = dL.exp)
     + ## compare
     + for(iX1 in 0:1){ # iX1 <- 0
     + newdata.L <- data.frame(X1.1 = c(iX1,0), X1.2 = c(0,iX1), trans = c(1, 2), strata = c(1, 2))
     + newdata <- data.frame(X1 = iX1)
     + pred.msfit <- msfit(e.coxph, newdata = newdata.L, trans = tmat)
     + suppressWarnings(
     + pred.probtrans <- probtrans(pred.msfit,0)[[1]]
     + )
     + pred.RR <- predict(CSC.exp, times = pred.probtrans[,"time"], newdata = newdata, cause = 1, product.limit = TRUE)
     + expect_equal(pred.probtrans[,"pstate2"],
     + as.double(pred.RR$absRisk)
     + )
     + }
     + }
     + })
     Test passed 🥳
     >
     > test_that("[predict.CSC]: check absolute risks add up to one",{
     + for(iData in 1:2){
     + df <- switch(as.character(iData),
     + "1"=df1,
     + "2"=df2)
     +
     + ## fit CSC
     + CSC.exp <- CSC(Hist(time,event) ~ X1, data = df, method = "breslow")
     +
     + ## compute all transition probabilites
     + seqTimes <- sort(unique(df$time))
     + nTimes <- length(seqTimes)
     +
     + for(iX1 in 0:1){
     + newdata <- data.frame(X1 = iX1)
     +
     + hazAll <- predictCox(CSC.exp$models[[1]], times = seqTimes, newdata = newdata, type = "hazard")$hazard+predictCox(CSC.exp$models[[2]], times = seqTimes, newdata = newdata, type = "hazard")$hazard
     + surv.prodlim <- cumprod(1-as.double(hazAll))
     + haz1.prodlim <- predict(CSC.exp, times = seqTimes, newdata = newdata, cause = 1, product.limit = TRUE)
     + haz2.prodlim <- predict(CSC.exp, times = seqTimes, newdata = newdata, cause = 2, product.limit = TRUE)
     +
     + expect_equal(surv.prodlim+as.double(haz1.prodlim$absRisk)+as.double(haz2.prodlim$absRisk),rep(1,nTimes))
     + }
     + }
     + })
     Test passed 🎊
     >
     > ## * [predictCSC] predict.CSC vs. mstate on simulated and real data
     > cat("[predictCSC] predict.CSC vs. mstate \n")
     [predictCSC] predict.CSC vs. mstate
     > ## ** Data
     > set.seed(10)
     > d <- sampleData(1e2, outcome = "competing.risks")[,.(time,event,X1,X2,X6)]
     > d[,X1:=as.numeric(as.character(X1))]
     > d[,X2:=as.numeric(as.character(X2))]
     > d[ , X16 := X1*X6]
     > d[ , Xcat2 := as.factor(paste0(X1,X2))]
     >
     > d <- d[event!=0]
     >
     > d[, event1 := as.numeric(event == 1)]
     > d[, event2 := as.numeric(event == 2)]
     >
     > ## ** Mstate
     > tmat <- trans.comprisk(2, names = c("0", "1", "2"))
     >
     > dL <- msprep(time = c(NA, "time", "time"),
     + status = c(NA,"event1", "event2"),
     + data = as.data.frame(d), keep = c("X1","X2","X16","Xcat2"),
     + trans = tmat)
     > dL.exp <- expand.covs(dL, c("X1","X2","X16","Xcat2"))
     >
     > ## ** No covariates
     >
     > test_that("predict.CSC (no covariates): compare to mstate",{
     + newdata <- data.frame(NA)
     + newdata.L <- data.frame(trans = c(1, 2), strata = c(1, 2))
     + # mstate
     + e.coxph <- coxph(Surv(time, status) ~1 + strata(trans),
     + data = dL.exp)
     +
     + pred.msfit <- msfit(e.coxph, newdata = newdata.L, trans = tmat)
     + pred.probtrans <- probtrans(pred.msfit,0)[[1]]
     +
     + ## riskRegression
     + CSC.RR1 <- CSC(Hist(time,event)~1, data = d, method = "breslow")
     +
     + pred.RR1a <- predict(CSC.RR1, newdata, cause = 1, time = pred.probtrans[,"time"],
     + keep.newdata = FALSE, se = TRUE, product.limit = TRUE)
     +
     + pred.RR1b <- predict(CSC.RR1, newdata, cause = 1, time = pred.probtrans[,"time"],
     + keep.newdata = FALSE, se = TRUE, product.limit = FALSE)
     +
     + pred.RR1c <- predict(CSC.RR1, newdata, cause = 2, time = pred.probtrans[,"time"],
     + keep.newdata = FALSE, se = TRUE, product.limit = TRUE)
     +
     + pred.RR1d <- predict(CSC.RR1, newdata, cause = 2, time = pred.probtrans[,"time"],
     + keep.newdata = FALSE, se = TRUE, product.limit = FALSE)
     +
     + CSC.RR2 <- CSC(Hist(time,event)~1, data = d, surv.type = "survival", method = "breslow")
     + pred.RR2a <- predict(CSC.RR2, newdata, cause = 1, time = pred.probtrans[,"time"],
     + keep.newdata = FALSE, se = TRUE, product.limit = TRUE)
     +
     + pred.RR2b <- predict(CSC.RR2, newdata, cause = 1, time = pred.probtrans[,"time"],
     + keep.newdata = FALSE, se = TRUE, product.limit = FALSE)
     +
     +
     + expect_equal(as.double(pred.RR1a$absRisk),pred.probtrans[,"pstate2"])
     + expect_equal(as.double(pred.RR1b$absRisk),pred.probtrans[,"pstate2"], tol = 5e-3)
     + expect_equal(as.double(pred.RR1c$absRisk),pred.probtrans[,"pstate3"])
     + expect_equal(as.double(pred.RR1d$absRisk),pred.probtrans[,"pstate3"], tol = 5e-2)
     + expect_equal(as.double(pred.RR2a$absRisk),pred.probtrans[,"pstate2"])
     + expect_equal(as.double(pred.RR2b$absRisk),pred.probtrans[,"pstate2"], tol = 5e-3)
     +
     + # quantile(as.double(pred.RR1a$absRisk.se) - pred.probtrans[,"se2"])
     + expect_equal(as.double(pred.RR1a$absRisk.se),pred.probtrans[,"se2"], tol = 5e-3)
     + expect_equal(as.double(pred.RR1b$absRisk.se),pred.probtrans[,"se2"], tol = 5e-3)
     + expect_equal(as.double(pred.RR1c$absRisk.se),pred.probtrans[,"se3"], tol = 5e-3)
     + expect_equal(as.double(pred.RR1d$absRisk.se),pred.probtrans[,"se3"], tol = 5e-3)
     + expect_equal(as.double(pred.RR2a$absRisk.se),pred.probtrans[,"se2"], tol = 5e-3)
     + expect_equal(as.double(pred.RR2b$absRisk.se),pred.probtrans[,"se2"], tol = 5e-3)
     + })
     Test passed 😸
     >
     > ## ** With covariates
     > test_that("predict.CSC (covariates): compare to mstate",{
     + for(iX in 0:1){
     + newdata <- data.frame(X1 = iX, X2 = 0, X16 = 0)
     + newdata.L <- data.frame(X1.1 = c(iX, 0), X1.2 = c(0, iX),
     + X2.1 = c(0, 0), X2.2 = c(0, 0),
     + X16.1 = c(0, 0), X16.2 = c(0, 0),
     + trans = c(1, 2), strata = c(1, 2))
     + # mstate
     + e.coxph <- coxph(Surv(time, status) ~ X1.1 + X1.2 + X2.1 + X2.2 + X16.1 + X16.2 + strata(trans),
     + data = dL.exp)
     + pred.msfit <- msfit(e.coxph, newdata = newdata.L, trans = tmat)
     + suppressWarnings(
     + pred.probtrans <- probtrans(pred.msfit,0)[[1]]
     + )
     + ## riskRegression
     + CSC.RR1 <- CSC(Hist(time,event)~X1+X2+X16, data = d, method = "breslow")
     + pred.RR1a <- predict(CSC.RR1, newdata, cause = 1, time = pred.probtrans[,"time"],
     + keep.newdata = FALSE, se = TRUE, product.limit = TRUE)
     + pred.RR1b <- predict(CSC.RR1, newdata, cause = 1, time = pred.probtrans[,"time"],
     + keep.newdata = FALSE, se = TRUE, product.limit = FALSE)
     + pred.RR1c <- predict(CSC.RR1, newdata, cause = 2, time = pred.probtrans[,"time"],
     + keep.newdata = FALSE, se = TRUE, product.limit = TRUE)
     + pred.RR1d <- predict(CSC.RR1, newdata, cause = 2, time = pred.probtrans[,"time"],
     + keep.newdata = FALSE, se = TRUE, product.limit = FALSE)
     + CSC.RR2 <- CSC(Hist(time,event)~X1+X2+X16, data = d, surv.type = "survival", method = "breslow")
     + pred.RR2a <- predict(CSC.RR2, newdata, cause = 1, time = pred.probtrans[,"time"],
     + keep.newdata = FALSE, se = TRUE, product.limit = TRUE)
     + pred.RR2b <- predict(CSC.RR2, newdata, cause = 1, time = pred.probtrans[,"time"],
     + keep.newdata = FALSE, se = TRUE, product.limit = FALSE)
     + expect_equal(as.double(pred.RR1a$absRisk),pred.probtrans[,"pstate2"])
     + expect_equal(as.double(pred.RR1b$absRisk),pred.probtrans[,"pstate2"], tol = 5e-3)
     + expect_equal(as.double(pred.RR1c$absRisk),pred.probtrans[,"pstate3"])
     + expect_equal(as.double(pred.RR1d$absRisk),pred.probtrans[,"pstate3"], tol = 1e-2)
     + expect_equal(as.double(pred.RR2a$absRisk),pred.probtrans[,"pstate2"], tol = 1e-2)
     + expect_equal(as.double(pred.RR2b$absRisk),pred.probtrans[,"pstate2"], tol = 1e-1)
     + #if(iX==0){
     + # quantile(as.double(pred.RR1a$absRisk.se) - pred.probtrans[,"se2"])
     + # quantile(as.double(pred.RR1c$absRisk.se) - pred.probtrans[,"se3"])
     + expect_equal(as.double(pred.RR1a$absRisk.se),pred.probtrans[,"se2"], tol = 1e-2)
     + expect_equal(as.double(pred.RR1b$absRisk.se),pred.probtrans[,"se2"], tol = 1e-2)
     + if(iX==0){
     + expect_equal(as.double(pred.RR1c$absRisk.se),pred.probtrans[,"se3"], tol = 1e-2)
     + expect_equal(as.double(pred.RR1d$absRisk.se),pred.probtrans[,"se3"], tol = 1e-2)
     + }
     + expect_equal(as.double(pred.RR2a$absRisk.se),pred.probtrans[,"se2"], tol = 1e-2)
     + expect_equal(as.double(pred.RR2b$absRisk.se),pred.probtrans[,"se2"], tol = 1e-2)
     + ## }else{
     + ## expect_equal(as.double(pred.RR1a$absRisk.se),pred.probtrans[,"se2"], tol = 1e-1)
     + ## expect_equal(as.double(pred.RR1b$absRisk.se),pred.probtrans[,"se2"], tol = 1e-1)
     + ## expect_equal(as.double(pred.RR1c$absRisk.se),pred.probtrans[,"se3"], tol = 1e-1)
     + ## expect_equal(as.double(pred.RR1d$absRisk.se),pred.probtrans[,"se3"], tol = 1e-1)
     + ## expect_equal(as.double(pred.RR2a$absRisk.se),pred.probtrans[,"se2"], tol = 1e-1)
     + ## expect_equal(as.double(pred.RR2b$absRisk.se),pred.probtrans[,"se2"], tol = 1e-1)
     + ## }
     + }
     + })
     Test passed 🥳
     >
     > ## ** Strata
     > test_that("predict.CSC (strata): compare to mstate",{
     +
     + newdata <- data.frame(X1 = 0, X2 = 0, X16 = 0)
     + newdata.L <- data.frame(X1.1 = c(0, 0), X1.2 = c(0, 0),
     + X2.1 = c(0, 0), X2.2 = c(0, 0),
     + X16.1 = c(0, 0), X16.2 = c(0, 0),
     + trans = c(1, 2), strata = c(1, 2))
     +
     + e.coxph <- coxph(Surv(time, status) ~ X1.1 + X1.2 + X2.1 + X2.2 + X16.1 + X16.2 + strata(trans),
     + data = dL.exp)
     +
     + pred.msfit <- msfit(e.coxph, newdata = newdata.L, trans = tmat)
     + suppressWarnings(
     + pred.probtrans <- probtrans(pred.msfit,0)[[1]]
     + )
     +
     + ## riskRegression no strata
     + CSC.RR1 <- CSC(Hist(time,event)~X1+X2+X16, data = d, method = "breslow")
     + pred.RR1a <- predict(CSC.RR1, newdata, cause = 1, time = pred.probtrans[,"time"],
     + keep.newdata = FALSE, se = TRUE, product.limit = TRUE)
     +
     + pred.RR1b <- predict(CSC.RR1, newdata, cause = 1, time = pred.probtrans[,"time"],
     + keep.newdata = FALSE, se = TRUE, product.limit = FALSE)
     +
     + pred.RR1c <- predict(CSC.RR1, newdata, cause = 2, time = pred.probtrans[,"time"],
     + keep.newdata = FALSE, se = TRUE, product.limit = TRUE)
     +
     + pred.RR1d <- predict(CSC.RR1, newdata, cause = 2, time = pred.probtrans[,"time"],
     + keep.newdata = FALSE, se = TRUE, product.limit = FALSE)
     +
     + CSC.RR2 <- CSC(Hist(time,event)~X1+X2+X16, data = d, surv.type = "survival", method = "breslow")
     + pred.RR2a <- predict(CSC.RR2, newdata, cause = 1, time = pred.probtrans[,"time"],
     + keep.newdata = FALSE, se = TRUE, product.limit = TRUE)
     + pred.RR2b <- predict(CSC.RR2, newdata, cause = 1, time = pred.probtrans[,"time"],
     + keep.newdata = FALSE, se = TRUE, product.limit = FALSE)
     +
     + ## riskRegression strata
     + d2 <- rbind(cbind(d,grp=1),cbind(d,grp=2))
     + CSC.RR1_strata <- CSC(Hist(time,event)~X1+X2+X16+strata(grp), data = d2, method = "breslow")
     + pred.RR1a_strata <- predict(CSC.RR1_strata, cbind(newdata,grp=1), cause = 1, time = pred.probtrans[,"time"],
     + keep.newdata = FALSE, se = TRUE, product.limit = TRUE)
     +
     + pred.RR1b_strata <- predict(CSC.RR1_strata, cbind(newdata,grp=1), cause = 1, time = pred.probtrans[,"time"],
     + keep.newdata = FALSE, se = TRUE, product.limit = FALSE)
     +
     + pred.RR1c_strata <- predict(CSC.RR1_strata, cbind(newdata,grp=1), cause = 2, time = pred.probtrans[,"time"],
     + keep.newdata = FALSE, se = TRUE, product.limit = TRUE)
     +
     + pred.RR1d_strata <- predict(CSC.RR1_strata, cbind(newdata,grp=1), cause = 2, time = pred.probtrans[,"time"],
     + keep.newdata = FALSE, se = TRUE, product.limit = FALSE)
     +
     + expect_equal(pred.RR1a_strata$absRisk,pred.RR1a$absRisk)
     + expect_equal(pred.RR1b_strata$absRisk,pred.RR1b$absRisk)
     + expect_equal(pred.RR1c_strata$absRisk,pred.RR1c$absRisk)
     + expect_equal(pred.RR1d_strata$absRisk,pred.RR1d$absRisk)
     +
     + CSC.RR2_strata<- CSC(Hist(time,event)~X1+X2+X16+strata(grp), data = d2, surv.type = "survival", method = "breslow")
     + pred.RR2a_strata <- predict(CSC.RR2_strata, cbind(newdata,grp=1), cause = 1, time = pred.probtrans[,"time"],
     + keep.newdata = FALSE, se = TRUE, product.limit = TRUE)
     + pred.RR2b_strata <- predict(CSC.RR2_strata, cbind(newdata,grp=1), cause = 1, time = pred.probtrans[,"time"],
     + keep.newdata = FALSE, se = TRUE, product.limit = FALSE)
     +
     + expect_equal(pred.RR2a_strata$absRisk,pred.RR2a$absRisk)
     + expect_equal(pred.RR2b_strata$absRisk,pred.RR2b$absRisk)
     +
     +
     + })
     Test passed 🎉
     > ## ** Vs. fixed values
     > set.seed(10)
     > n <- 300
     > df.S <- prodlim::SimCompRisk(n)
     > df.S$time <- round(df.S$time,2)
     > df.S$X3 <- rbinom(n, size = 4, prob = rep(0.25,4))
     > seqTime <- c(unique(sort(df.S$time)), max(df.S$time) + 1)[c(1,2,5,12,90,125,200,241,267,268)]
     >
     > test_that("[predictCSC] vs. fixed values",{
     + CSC.S <- CSC(Hist(time,event) ~ strata(X1) + strata(X3) + X2, data = df.S, ties = "efron", fitter = "coxph")
     +
     + ## cause 1
     + Event1.S <- predict(CSC.S, newdata = df.S[1:5,], times = seqTime, cause = 1,
     + se = TRUE, keep.newdata = FALSE)
     +
     + ## butils::object2script(as.data.table(Event1.S), digit = 6)
     + GS <- data.table("observation" = c(1, 2, 3, 4, 5, 1, 2, 3, 4, 5, 1, 2, 3, 4, 5, 1, 2, 3, 4, 5, 1, 2, 3, 4, 5, 1, 2, 3, 4, 5, 1, 2, 3, 4, 5, 1, 2, 3, 4, 5, 1, 2, 3, 4, 5, 1, 2, 3, 4, 5),
     + "times" = c( 0.14, 0.14, 0.14, 0.14, 0.14, 0.36, 0.36, 0.36, 0.36, 0.36, 0.56, 0.56, 0.56, 0.56, 0.56, 0.94, 0.94, 0.94, 0.94, 0.94, 3.23, 3.23, 3.23, 3.23, 3.23, 4.12, 4.12, 4.12, 4.12, 4.12, 6.75, 6.75, 6.75, 6.75, 6.75, 9.19, 9.19, 9.19, 9.19, 9.19, 15.90, 15.90, 15.90, 15.90, 15.90, 16.90, 16.90, 16.90, 16.90, 16.90),
     + "strata" = c("X1=1 X3=1", "X1=1 X3=1", "X1=1 X3=2", "X1=1 X3=1", "X1=0 X3=2", "X1=1 X3=1", "X1=1 X3=1", "X1=1 X3=2", "X1=1 X3=1", "X1=0 X3=2", "X1=1 X3=1", "X1=1 X3=1", "X1=1 X3=2", "X1=1 X3=1", "X1=0 X3=2", "X1=1 X3=1", "X1=1 X3=1", "X1=1 X3=2", "X1=1 X3=1", "X1=0 X3=2", "X1=1 X3=1", "X1=1 X3=1", "X1=1 X3=2", "X1=1 X3=1", "X1=0 X3=2", "X1=1 X3=1", "X1=1 X3=1", "X1=1 X3=2", "X1=1 X3=1", "X1=0 X3=2", "X1=1 X3=1", "X1=1 X3=1", "X1=1 X3=2", "X1=1 X3=1", "X1=0 X3=2", "X1=1 X3=1", "X1=1 X3=1", "X1=1 X3=2", "X1=1 X3=1", "X1=0 X3=2", "X1=1 X3=1", "X1=1 X3=1", "X1=1 X3=2", "X1=1 X3=1", "X1=0 X3=2", "X1=1 X3=1", "X1=1 X3=1", "X1=1 X3=2", "X1=1 X3=1", "X1=0 X3=2"),
     + "absRisk" = c(0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.010614, 0.008321, 0.000000, 0.005059, 0.033321, 0.010614, 0.008321, 0.000000, 0.005059, 0.033321, 0.187045, 0.149948, 0.024562, 0.094195, 0.208568, 0.285231, 0.232052, 0.047504, 0.148967, 0.320621, 0.532564, 0.455475, 0.127580, 0.315670, 0.471993, 0.592726, 0.516242, 0.341738, 0.368635, 0.694506, NA, NA, 0.804087, NA, NA, NA, NA, NA, NA, NA),
     + "absRisk.se" = c(0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.010671, 0.008392, 0.000000, 0.005145, 0.021734, 0.010671, 0.008392, 0.000000, 0.005145, 0.021734, 0.045017, 0.037950, 0.013802, 0.026526, 0.060412, 0.053118, 0.046489, 0.021607, 0.034654, 0.079323, 0.074433, 0.071882, 0.047708, 0.061923, 0.084891, 0.075101, 0.075150, 0.086033, 0.068349, 0.150281, NA, NA, 0.134410, NA, NA, NA, NA, NA, NA, NA),
     + "absRisk.lower" = c(0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000901, 0.000720, 0.000000, 0.000450, 0.007060, 0.000901, 0.000720, 0.000000, 0.000450, 0.007060, 0.108482, 0.085060, 0.006812, 0.050584, 0.105229, 0.186727, 0.147890, 0.016864, 0.088994, 0.175143, 0.377876, 0.311802, 0.052900, 0.200010, 0.300988, 0.431337, 0.361334, 0.182670, 0.237798, 0.311357, NA, NA, 0.375421, NA, NA, NA, NA, NA, NA, NA),
     + "absRisk.upper" = c(0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.052517, 0.042026, 0.000000, 0.026619, 0.096719, 0.052517, 0.042026, 0.000000, 0.026619, 0.096719, 0.282169, 0.232014, 0.063687, 0.154095, 0.335805, 0.391505, 0.327430, 0.102897, 0.223449, 0.475827, 0.665044, 0.588205, 0.236373, 0.437746, 0.625333, 0.722293, 0.650868, 0.507569, 0.499892, 0.892347, NA, NA, 0.952629, NA, NA, NA, NA, NA, NA, NA))
     +
     + ## data.table("observation" = c(1, 2, 3, 4, 5, 1, 2, 3, 4, 5, 1, 2, 3, 4, 5, 1, 2, 3, 4, 5, 1, 2, 3, 4, 5, 1, 2, 3, 4, 5, 1, 2, 3, 4, 5, 1, 2, 3, 4, 5, 1, 2, 3, 4, 5, 1, 2, 3, 4, 5),
     + ## "times" = c( 0.14, 0.14, 0.14, 0.14, 0.14, 0.36, 0.36, 0.36, 0.36, 0.36, 0.56, 0.56, 0.56, 0.56, 0.56, 0.94, 0.94, 0.94, 0.94, 0.94, 3.23, 3.23, 3.23, 3.23, 3.23, 4.12, 4.12, 4.12, 4.12, 4.12, 6.75, 6.75, 6.75, 6.75, 6.75, 9.19, 9.19, 9.19, 9.19, 9.19, 15.90, 15.90, 15.90, 15.90, 15.90, 16.90, 16.90, 16.90, 16.90, 16.90),
     + ## "strata" = c("X1=1 X3=1", "X1=1 X3=1", "X1=1 X3=2", "X1=1 X3=1", "X1=0 X3=2", "X1=1 X3=1", "X1=1 X3=1", "X1=1 X3=2", "X1=1 X3=1", "X1=0 X3=2", "X1=1 X3=1", "X1=1 X3=1", "X1=1 X3=2", "X1=1 X3=1", "X1=0 X3=2", "X1=1 X3=1", "X1=1 X3=1", "X1=1 X3=2", "X1=1 X3=1", "X1=0 X3=2", "X1=1 X3=1", "X1=1 X3=1", "X1=1 X3=2", "X1=1 X3=1", "X1=0 X3=2", "X1=1 X3=1", "X1=1 X3=1", "X1=1 X3=2", "X1=1 X3=1", "X1=0 X3=2", "X1=1 X3=1", "X1=1 X3=1", "X1=1 X3=2", "X1=1 X3=1", "X1=0 X3=2", "X1=1 X3=1", "X1=1 X3=1", "X1=1 X3=2", "X1=1 X3=1", "X1=0 X3=2", "X1=1 X3=1", "X1=1 X3=1", "X1=1 X3=2", "X1=1 X3=1", "X1=0 X3=2", "X1=1 X3=1", "X1=1 X3=1", "X1=1 X3=2", "X1=1 X3=1", "X1=0 X3=2"),
     + ## "absRisk" = c(0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.010614, 0.008321, 0.000000, 0.005059, 0.033321, 0.010614, 0.008321, 0.000000, 0.005059, 0.033321, 0.187045, 0.149948, 0.024562, 0.094195, 0.208568, 0.285231, 0.232052, 0.047504, 0.148967, 0.320621, 0.532564, 0.455475, 0.127580, 0.315670, 0.471993, 0.592726, 0.516242, 0.341738, 0.368635, 0.694506, NA, NA, 0.804087, NA, NA, NA, NA, NA, NA, NA),
     + ## "absRisk.se" = c(0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.010671, 0.008392, 0.000000, 0.005145, 0.021735, 0.010671, 0.008392, 0.000000, 0.005145, 0.021735, 0.045082, 0.037994, 0.013808, 0.026548, 0.060483, 0.053248, 0.046578, 0.021627, 0.034699, 0.079498, 0.075317, 0.072541, 0.047915, 0.062280, 0.085417, 0.076575, 0.076338, 0.087346, 0.069072, 0.159006, NA, NA, 0.145491, NA, NA, NA, NA, NA, NA, NA))
     +
     + expect_equal(as.data.table(Event1.S)$absRisk, GS$absRisk, tolerance = 1e-5)
     + expect_equal(as.data.table(Event1.S)$absRisk.se, GS$absRisk.se, tolerance = 1e-5)
     +
     + ## cause 2
     + Event2.S <- predict(CSC.S, newdata = df.S[1:5,], times = seqTime, cause = 2,
     + se = TRUE, keep.newdata = FALSE)
     + ## butils::object2script(as.data.table(Event2.S), digit = 6)
     + GS <- data.table("observation" = c(1, 2, 3, 4, 5, 1, 2, 3, 4, 5, 1, 2, 3, 4, 5, 1, 2, 3, 4, 5, 1, 2, 3, 4, 5, 1, 2, 3, 4, 5, 1, 2, 3, 4, 5, 1, 2, 3, 4, 5, 1, 2, 3, 4, 5, 1, 2, 3, 4, 5),
     + "times" = c( 0.14, 0.14, 0.14, 0.14, 0.14, 0.36, 0.36, 0.36, 0.36, 0.36, 0.56, 0.56, 0.56, 0.56, 0.56, 0.94, 0.94, 0.94, 0.94, 0.94, 3.23, 3.23, 3.23, 3.23, 3.23, 4.12, 4.12, 4.12, 4.12, 4.12, 6.75, 6.75, 6.75, 6.75, 6.75, 9.19, 9.19, 9.19, 9.19, 9.19, 15.90, 15.90, 15.90, 15.90, 15.90, 16.90, 16.90, 16.90, 16.90, 16.90),
     + "strata" = c("X1=1 X3=1", "X1=1 X3=1", "X1=1 X3=2", "X1=1 X3=1", "X1=0 X3=2", "X1=1 X3=1", "X1=1 X3=1", "X1=1 X3=2", "X1=1 X3=1", "X1=0 X3=2", "X1=1 X3=1", "X1=1 X3=1", "X1=1 X3=2", "X1=1 X3=1", "X1=0 X3=2", "X1=1 X3=1", "X1=1 X3=1", "X1=1 X3=2", "X1=1 X3=1", "X1=0 X3=2", "X1=1 X3=1", "X1=1 X3=1", "X1=1 X3=2", "X1=1 X3=1", "X1=0 X3=2", "X1=1 X3=1", "X1=1 X3=1", "X1=1 X3=2", "X1=1 X3=1", "X1=0 X3=2", "X1=1 X3=1", "X1=1 X3=1", "X1=1 X3=2", "X1=1 X3=1", "X1=0 X3=2", "X1=1 X3=1", "X1=1 X3=1", "X1=1 X3=2", "X1=1 X3=1", "X1=0 X3=2", "X1=1 X3=1", "X1=1 X3=1", "X1=1 X3=2", "X1=1 X3=1", "X1=0 X3=2", "X1=1 X3=1", "X1=1 X3=1", "X1=1 X3=2", "X1=1 X3=1", "X1=0 X3=2"),
     + "absRisk" = c(0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.015510, 0.015342, 0.028621, 0.014982, 0.000000, 0.114820, 0.115596, 0.028621, 0.115825, 0.000000, 0.133554, 0.135236, 0.065157, 0.136670, 0.000000, 0.258482, 0.280051, 0.201121, 0.314747, 0.185563, 0.301218, 0.336703, 0.368645, 0.398416, 0.185563, NA, NA, 0.368645, NA, NA, NA, NA, NA, NA, NA),
     + "absRisk.se" = c(0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.015392, 0.015249, 0.027535, 0.014981, 0.000000, 0.040606, 0.041787, 0.027535, 0.044423, 0.000000, 0.043677, 0.045256, 0.044218, 0.048715, 0.000000, 0.057897, 0.063448, 0.082190, 0.076144, 0.083832, 0.063067, 0.069936, 0.106967, 0.085867, 0.083832, NA, NA, 0.106967, NA, NA, NA, NA, NA, NA, NA),
     + "absRisk.lower" = c(0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.001301, 0.001282, 0.002379, 0.001233, 0.000000, 0.050723, 0.049970, 0.002379, 0.047117, 0.000000, 0.062787, 0.062230, 0.011742, 0.059184, 0.000000, 0.153888, 0.164616, 0.071169, 0.175141, 0.057883, 0.184670, 0.205517, 0.171289, 0.233092, 0.057883, NA, NA, 0.171289, NA, NA, NA, NA, NA, NA, NA),
     + "absRisk.upper" = c(0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.073376, 0.072786, 0.123636, 0.071741, 0.000000, 0.207778, 0.211461, 0.123636, 0.218487, 0.000000, 0.231242, 0.236571, 0.186746, 0.246346, 0.000000, 0.376051, 0.407416, 0.377807, 0.464388, 0.369464, 0.426403, 0.472880, 0.568695, 0.559050, 0.369464, NA, NA, 0.568695, NA, NA, NA, NA, NA, NA, NA))
     +
     + ## data.table("observation" = c(1, 2, 3, 4, 5, 1, 2, 3, 4, 5, 1, 2, 3, 4, 5, 1, 2, 3, 4, 5, 1, 2, 3, 4, 5, 1, 2, 3, 4, 5, 1, 2, 3, 4, 5, 1, 2, 3, 4, 5, 1, 2, 3, 4, 5, 1, 2, 3, 4, 5),
     + ## "times" = c( 0.14, 0.14, 0.14, 0.14, 0.14, 0.36, 0.36, 0.36, 0.36, 0.36, 0.56, 0.56, 0.56, 0.56, 0.56, 0.94, 0.94, 0.94, 0.94, 0.94, 3.23, 3.23, 3.23, 3.23, 3.23, 4.12, 4.12, 4.12, 4.12, 4.12, 6.75, 6.75, 6.75, 6.75, 6.75, 9.19, 9.19, 9.19, 9.19, 9.19, 15.90, 15.90, 15.90, 15.90, 15.90, 16.90, 16.90, 16.90, 16.90, 16.90),
     + ## "strata" = c("X1=1 X3=1", "X1=1 X3=1", "X1=1 X3=2", "X1=1 X3=1", "X1=0 X3=2", "X1=1 X3=1", "X1=1 X3=1", "X1=1 X3=2", "X1=1 X3=1", "X1=0 X3=2", "X1=1 X3=1", "X1=1 X3=1", "X1=1 X3=2", "X1=1 X3=1", "X1=0 X3=2", "X1=1 X3=1", "X1=1 X3=1", "X1=1 X3=2", "X1=1 X3=1", "X1=0 X3=2", "X1=1 X3=1", "X1=1 X3=1", "X1=1 X3=2", "X1=1 X3=1", "X1=0 X3=2", "X1=1 X3=1", "X1=1 X3=1", "X1=1 X3=2", "X1=1 X3=1", "X1=0 X3=2", "X1=1 X3=1", "X1=1 X3=1", "X1=1 X3=2", "X1=1 X3=1", "X1=0 X3=2", "X1=1 X3=1", "X1=1 X3=1", "X1=1 X3=2", "X1=1 X3=1", "X1=0 X3=2", "X1=1 X3=1", "X1=1 X3=1", "X1=1 X3=2", "X1=1 X3=1", "X1=0 X3=2", "X1=1 X3=1", "X1=1 X3=1", "X1=1 X3=2", "X1=1 X3=1", "X1=0 X3=2"),
     + ## "absRisk" = c(0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.015510, 0.015342, 0.028621, 0.014982, 0.000000, 0.114820, 0.115596, 0.028621, 0.115825, 0.000000, 0.133554, 0.135236, 0.065157, 0.136670, 0.000000, 0.258482, 0.280051, 0.201121, 0.314747, 0.185563, 0.301218, 0.336703, 0.368645, 0.398416, 0.185563, NA, NA, 0.368645, NA, NA, NA, NA, NA, NA, NA),
     + ## "absRisk.se" = c(0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.015393, 0.015250, 0.027535, 0.014981, 0.000000, 0.040650, 0.041823, 0.027535, 0.044450, 0.000000, 0.043741, 0.045308, 0.044233, 0.048753, 0.000000, 0.058568, 0.064056, 0.082396, 0.076642, 0.085257, 0.064560, 0.071313, 0.108398, 0.087042, 0.085257, NA, NA, 0.108398, NA, NA, NA, NA, NA, NA, NA))
     +
     + expect_equal(as.data.table(Event2.S)$absRisk, GS$absRisk, tolerance = 1e-5)
     + expect_equal(as.data.table(Event2.S)$absRisk.se, GS$absRisk.se, tolerance = 1e-5)
     +
     + })
     Test passed 🥳
     >
     > ## ** Melanoma
     > data(Melanoma, package = "riskRegression")
     > Melanoma$index <- 1:NROW(Melanoma)
     >
     > ## CSC
     > cfit1 <- CSC(formula=list(Hist(time,status)~age+logthick+epicel+sex,
     + Hist(time,status)~age+sex),
     + data=Melanoma)
     > cfit2 <- CSC(formula=list(Hist(time,status)~age+logthick+epicel+sex,
     + Hist(time,status)~age+sex),
     + data=Melanoma, surv.type = "survival")
     >
     > ## mstrate
     > Melanoma$event1 <- as.numeric(Melanoma$event == "death.malignant.melanoma")
     > Melanoma$event2 <- as.numeric(Melanoma$event == "death.other.causes")
     >
     > MelanomaL <- msprep(time = c(NA, "time", "time"),
     + status = c(NA,"event1", "event2"),
     + data = Melanoma, keep = c("age","logthick","epicel","sex","index"),
     + trans = tmat)
     > MelanomaL.exp <- expand.covs(MelanomaL, c("age","logthick","epicel","sex"))
     >
     > Melanoma.coxph <- coxph(Surv(time, status) ~ age.1 + age.2 + logthick.1 + epicelpresent.1 + sexMale.1 + sexMale.2 + strata(trans),
     + data = MelanomaL.exp)
     >
     > newdata <- Melanoma[1,,drop=FALSE]
     > newdata.exp <- MelanomaL.exp[MelanomaL.exp$index %in% newdata$index,,drop=FALSE]
     > newdata.exp$strata <- newdata.exp$trans
     >
     > pred.msfit <- msfit(Melanoma.coxph, newdata = newdata.exp, trans = tmat, variance = FALSE)
Flavor: r-devel-windows-x86_64-new-UL