0

我在循环中安装了 100 个 coxph 模型的副本。我正在尝试为数据框或列表中的每个复制提取具有 p 值的对数秩分数测试结果。我正在使用以下内容。但是，它只给我对数排名分数，而不是 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)))
4

1 回答 1

2

这是我从 coxph 函数的帮助页面中获取的一个模拟示例。我刚刚复制了 100 次数据集来创建您的场景。我强烈建议开始使用这些tidyverse包来做这样的工作。与和broom一起是一个很好的补充。dplyrtidyr

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))

下面我使用该replicate函数复制数据集 100 次。

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

我将 cox 模型设置为一个小函数，我可以将它们传递给数据帧的每个副本。

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

下面是拟合模型和提取值的分步过程。

tidyr::nestpurrr::map模型进入每个嵌套 的数据框 nest是函数 inlibrary(tidyr) map是类似于lapplyin的函数library(purrr)

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

查看第一个代表以查看 coxph 输出。您将看到模型对象存储在数据框的单元格中，从而更容易访问。

nested %>% filter(rep==1)

对于每个模型对象，现在使用 broom 将模型中的参数估计和预测放入嵌套数据集中

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

tidyr::unnest查看您对拟合每个重采样数据集的预测

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

在这种情况下，由于我重复相同的数据集 100 次，pvalue 将是相同的，但在您的情况下，您将有 100 个不同的数据集，因此有 100 个不同的 p.values。

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

维杰

于 2016-12-04T18:40:35.923 回答