Извлечь результаты теста ранг (балл) с p-значением для модели Кокша

У меня есть 100 копий модели Coxph, встроенной в петлю. Я пытаюсь извлечь результаты теста оценки лог-ранга с p-значениями для каждой реплики во фрейме данных или списке. Я использую следующее. Но это дает мне только оценку рейтинга журнала, а не р-значение. Любая помощь будет очень ценится.

Я могу поделиться набором данных, но не уверен, как прикрепить здесь.

спасибо крина

Repl_List <- unique(dat3$Repl)
doLogRank = function(sel_name) {
dum <- dat3[dat3$Repl == sel_name,]
reg <- with(dum, coxph(Surv(TIME_day, STATUS) ~ Treatment, ties = "breslow"))
LogRank <- with(reg, reg$score) 
}
LogRank <- t(as.data.frame(lapply(Repl_List, doLogRank)))

1 ответ

Решение

Вот пример, который я взял со страницы помощи функции coxph. Я просто повторил набор данных 100 раз, чтобы создать сценарий. Я настоятельно рекомендую начать использовать tidyverse пакеты для такой работы. broom является отличным дополнением вместе с dplyr а также tidyr,

library(survival)
library(tidyverse)
library(broom)
  test <- data.frame(time=c(4,3,1,1,2,2,3), 
              status=c(1,1,1,0,1,1,0), 
              x=c(0,2,1,1,1,0,0), 
              sex=c(0,0,0,0,1,1,1))

Ниже я повторяю набор данных 100 раз, используя replicate функция.

r <- replicate(test,n = 100,simplify = FALSE) %>% bind_rows %>% 
  mutate(rep = rep(seq(1,100,1),each=7))

Я настроил модель Кокса как небольшую функцию, которую я могу передать каждой копии кадра данных.

cxph_mod <- function(df) {
  coxph(Surv(time, status) ~ x + strata(sex), df)
}

Ниже приведен пошаговый процесс подгонки модели и извлечения значений.

tidyr::nest кадр данных purrr::map модель в каждое гнездоnest функция в library(tidyr)map функция похожа на lapply в library(purrr)

nested <- r %>% 
  group_by(rep) %>% 
  nest %>% 
  mutate(model = data %>% map(cxph_mod))

загляните в первый повтор, чтобы увидеть вывод coxph. Вы увидите объект модели, хранящийся в ячейках кадра данных, что облегчает доступ.

nested %>% filter(rep==1)

Для каждого объекта модели теперь используйте метлу для получения оценок параметров и прогноза из модели во вложенный набор данных.

nested <- nested %>% 
  mutate(
    ests = model %>% map(broom::tidy)
  )

tidyr::unnest чтобы просмотреть ваши прогнозы по подгонке каждого набора данных

ests <- unnest(nested,ests,.drop=TRUE) %>% dplyr::select(rep,estimate:conf.high)

В этом случае, поскольку я повторяю один и тот же набор данных 100 раз, значение будет одинаковым, но в вашем случае у вас будет 100 различных наборов данных и, следовательно, 100 различных значений p.value.

ggplot(data=ests,aes(y=p.value,x=rep))+geom_point()

Виджай

Другие вопросы по тегам