## Loading required package: Matrix
## Loaded glmnet 4.1-10
library(survival)
data(pbc)
pbc2 <- pbc[!is.na(pbc$trt), ]
pbc2$event <- as.integer(pbc$status[!is.na(pbc$trt)] == 2)
pbc2$sex_n <- as.integer(pbc2$sex == "f")
feat_cols <- c("trt","age","sex_n","ascites","hepato","spiders","edema",
"bili","chol","albumin","copper","alk.phos","ast",
"trig","platelet","protime","stage")
df <- pbc2[, c("time", "event", feat_cols)]
for (col in feat_cols)
if (any(is.na(df[[col]])))
df[[col]][is.na(df[[col]])] <- median(df[[col]], na.rm = TRUE)
set.seed(42)
idx_train <- sample(nrow(df), floor(0.75 * nrow(df)))
train <- df[idx_train, ]; test <- df[-idx_train, ]
X_tr <- as.matrix(train[, feat_cols])
X_te <- as.matrix(test[, feat_cols])
regr_lm <- function(X, y, ...) lm(y ~ ., data = data.frame(X, y = y))
fit_boost_lm <- survivalisttoo::cox_gradient_boost(X_tr, train$time, train$event,
regr_lm)
## Iter 1/100 [ ] 1% ETA 1sIter 2/100 [= ] 2% ETA 1sIter 3/100 [== ] 3% ETA 1sIter 4/100 [== ] 4% ETA 1sIter 5/100 [== ] 5% ETA 1sIter 6/100 [=== ] 6% ETA 1sIter 7/100 [==== ] 7% ETA 1sIter 8/100 [==== ] 8% ETA 1sIter 9/100 [==== ] 9% ETA 1sIter 10/100 [===== ] 10% ETA 1sIter 11/100 [====== ] 11% ETA 1sIter 12/100 [====== ] 12% ETA 1sIter 13/100 [====== ] 13% ETA 1sIter 14/100 [======= ] 14% ETA 1sIter 15/100 [======== ] 15% ETA 1sIter 16/100 [======== ] 16% ETA 1sIter 17/100 [======== ] 17% ETA 1sIter 18/100 [========= ] 18% ETA 1sIter 19/100 [========== ] 19% ETA 1sIter 20/100 [========== ] 20% ETA 1sIter 21/100 [========== ] 21% ETA 1sIter 22/100 [=========== ] 22% ETA 1sIter 23/100 [============ ] 23% ETA 1sIter 24/100 [============ ] 24% ETA 1sIter 25/100 [============ ] 25% ETA 1sIter 26/100 [============= ] 26% ETA 1sIter 27/100 [============== ] 27% ETA 1sIter 28/100 [============== ] 28% ETA 1sIter 29/100 [============== ] 29% ETA 1sIter 30/100 [=============== ] 30% ETA 1sIter 31/100 [================ ] 31% ETA 1sIter 32/100 [================ ] 32% ETA 1sIter 33/100 [================ ] 33% ETA 1sIter 34/100 [================= ] 34% ETA 1sIter 35/100 [================== ] 35% ETA 1sIter 36/100 [================== ] 36% ETA 1sIter 37/100 [================== ] 37% ETA 1sIter 38/100 [=================== ] 38% ETA 1sIter 39/100 [==================== ] 39% ETA 1sIter 40/100 [==================== ] 40% ETA 1sIter 41/100 [==================== ] 41% ETA 1sIter 42/100 [===================== ] 42% ETA 1sIter 43/100 [====================== ] 43% ETA 1sIter 44/100 [====================== ] 44% ETA 1sIter 45/100 [====================== ] 45% ETA 1sIter 46/100 [======================= ] 46% ETA 1sIter 47/100 [======================== ] 47% ETA 1sIter 48/100 [======================== ] 48% ETA 1sIter 49/100 [======================== ] 49% ETA 1sIter 50/100 [========================= ] 50% ETA 1sIter 51/100 [========================== ] 51% ETA 1sIter 52/100 [========================== ] 52% ETA 1sIter 53/100 [========================== ] 53% ETA 1sIter 54/100 [=========================== ] 54% ETA 1sIter 55/100 [============================ ] 55% ETA 1sIter 56/100 [============================ ] 56% ETA 1sIter 57/100 [============================ ] 57% ETA 1sIter 58/100 [============================= ] 58% ETA 1sIter 59/100 [============================== ] 59% ETA 0sIter 60/100 [============================== ] 60% ETA 0sIter 61/100 [============================== ] 61% ETA 0sIter 62/100 [=============================== ] 62% ETA 0sIter 63/100 [================================ ] 63% ETA 0sIter 64/100 [================================ ] 64% ETA 0sIter 65/100 [================================ ] 65% ETA 0sIter 66/100 [================================= ] 66% ETA 0sIter 67/100 [================================== ] 67% ETA 0sIter 68/100 [================================== ] 68% ETA 0sIter 69/100 [================================== ] 69% ETA 0sIter 70/100 [=================================== ] 70% ETA 0sIter 71/100 [==================================== ] 71% ETA 0sIter 72/100 [==================================== ] 72% ETA 0sIter 73/100 [==================================== ] 73% ETA 0sIter 74/100 [===================================== ] 74% ETA 0sIter 75/100 [====================================== ] 75% ETA 0sIter 76/100 [====================================== ] 76% ETA 0sIter 77/100 [====================================== ] 77% ETA 0sIter 78/100 [======================================= ] 78% ETA 0sIter 79/100 [======================================== ] 79% ETA 0sIter 80/100 [======================================== ] 80% ETA 0sIter 81/100 [======================================== ] 81% ETA 0sIter 82/100 [========================================= ] 82% ETA 0sIter 83/100 [========================================== ] 83% ETA 0sIter 84/100 [========================================== ] 84% ETA 0sIter 85/100 [========================================== ] 85% ETA 0sIter 86/100 [=========================================== ] 86% ETA 0sIter 87/100 [============================================ ] 87% ETA 0sIter 88/100 [============================================ ] 88% ETA 0sIter 89/100 [============================================ ] 89% ETA 0sIter 90/100 [============================================= ] 90% ETA 0sIter 91/100 [============================================== ] 91% ETA 0sIter 92/100 [============================================== ] 92% ETA 0sIter 93/100 [============================================== ] 93% ETA 0sIter 94/100 [=============================================== ] 94% ETA 0sIter 95/100 [================================================ ] 95% ETA 0sIter 96/100 [================================================ ] 96% ETA 0sIter 97/100 [================================================ ] 97% ETA 0sIter 98/100 [================================================= ] 98% ETA 0sIter 99/100 [==================================================] 99% ETA 0sIter 100/100 [==================================================] 100% ETA 0s
fit_cox <- coxph(Surv(time, event) ~ .,
data = train[, c("time","event",feat_cols)], x = TRUE)
y_te <- Surv(test$time, test$event)
ci_blm <- glmnet::Cindex(predict(fit_boost_lm, X_te), y_te)
ci_cox <- glmnet::Cindex(predict(fit_cox, newdata = test), y_te)
cat("\n=== Test-set C-index ===\n")
##
## === Test-set C-index ===
cat(sprintf(" CoxBoost (LM, M=100): %.4f\n", ci_blm))
## CoxBoost (LM, M=100): 0.8075
cat(sprintf(" Classical Cox PH : %.4f\n", ci_cox))
## Classical Cox PH : 0.8069