Есть ли способ получить вероятность события из модели Коксме?

Question

Я хотел бы получить прогноз вероятности события из модели Кокса, для заданных значений моих объясняющих переменных и на заданном временном шаге. Мои данные структурированы так:

# Generate data
set.seed(123)
mydata <- data.frame(Site = as.factor(sample(c("SiteA", "SiteB", "SiteC"), 1000, replace = TRUE)), 
                     Block = as.factor(sample(c("A", "B", "C", "D", "E", "F"), 1000, replace = TRUE)),
                     Treatment = as.factor(sample(c("Treat.A", "Treat.B"), 1000, replace = TRUE)),
                     Origin = as.factor(sample(c("Origin.A", "Origin.B"), 1000, replace = TRUE)),
                     Time = sample(seq(3), 1000, replace = TRUE), 
                     Surv = sample(c(0, 1), 1000, replace = TRUE)) # Alive is 0, death is 1

# Coxme model
mymodel <- coxme(Surv(Time , Surv) ~ Treatment*Origin + 
                   (1|Site/Block), data = mydata)

Я хочу получить прогнозируемую вероятность выживания в момент времени 3 для каждой комбинации Лечение: Происхождение. Если бы у меня была модель Кокса (, то есть без случайных эффектов), это можно было бы легко сделать с помощью выживания из пакета выживания:

# use expand.grid to get a table with all possible combinations of Site and Treatment
newdata.surv <- with(mydata, expand.grid(Site = unique(Origin), Treatment = unique(Treatment)))

# run survfit to predict the new values
fitted <- survival::survfit(mymodel, newdata = newdata.surv)

# extract the fitted values for the time slice of interest: 3
newdata.surv$fit <- fitted$surv[3,]
newdata.surv$lower <- fitted$lower[3,] # Lower confidence interval
newdata.surv$upper <- fitted$upper[3,] # Upper confidence interval

Однако survfit не работает с coxme object. Я знаю, что predict.coxme существует в coxme package, но когда я пытаюсь его использовать, я всегда получаю сообщение об ошибке: «не могу найти функцию« Foregnet.coxme ». Я использую версию 2.2-10 coxme пакет, так что должна быть включена функцияgnast.coxme (см. https://cran.r-project.org/web/packages/coxme/news.html).

Я видел, что coxme объекты поддерживаются пакетами emmeans и lsmeans, но Я не уверен, что эти пакеты могут быть использованы для прогнозирования выживания. Заранее спасибо за вашу помощь.

Answer 1

Функция coxme.predict сама по себе не экспортируется, но вы можете вызвать ее с помощью predict(mymodel), которая затем вызовет этот метод (или вы можете напрямую вызвать coxme:::predict.coxme(mymodel) (с 3 двоеточиями)). См. ?coxme:::predict.coxme для краткого описания. Похоже, что в настоящее время он не поддерживает аргумент newdata, поэтому я не уверен, насколько он полезен для вашего случая использования.