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

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

# 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 для каждой комбинации лечение: происхождение. Если бы у меня была модель coxph (т.е. без случайных эффектов), это можно было бы легко сделать с помощью выживания из пакета Survival:

# 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, но когда я пытаюсь его использовать, я всегда получаю сообщение об ошибке: «не удалось найти функцию» pred.coxme. Я использую версию 2.2-10 пакета coxme, поэтому прогнозирую Должна быть включена функция .coxme (см. https://cran.r-project.org/web/packages/coxme/news.html).

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


r survival lsmeans
person Julien Barrère    schedule 29.02.2020    source источник

Ответы (1)


arrow_upward
2
arrow_downward

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

person user12728748    schedule 01.03.2020
comment
Действительно, я использовал только два двоеточия, поэтому я не смог заставить функцию работать (coxme::predict.coxme), поэтому теперь она работает. Тем не менее, единственный тип прогноза, при котором функция, возвращаемая функцией, является либо линейным предиктором, либо абсолютным риском, и ни одна из этих переменных не может использоваться для вычисления вероятности выживания (насколько мне известно). Есть другие идеи? - person Julien Barrère; 01.03.2020
comment
См. stats.stackexchange.com/questions/45538/ для обсуждения и потенциального обходного пути. Кажется, с тех пор ничего не изменилось. - person user12728748; 02.03.2020