数据解释和可视化

Question

我正在尝试使用母亲的年龄和婴儿的性别作为协变量来拟合 Cox 比例风险模型。

我在用 R 可视化表示数据时遇到问题。

我会 post 我的代码：

这些是我的包裹：

#Getting started:

# load libraries
pkgTest <- function(pkg){
  new.pkg <- pkg[!(pkg %in% installed.packages()[,  "Package"])]
  if (length(new.pkg)) 
    install.packages(new.pkg,  dependencies = TRUE)
  sapply(pkg,  require,  character.only = TRUE)
}

lapply(c("survival", "eha", "tidyverse", "ggfortify", "stargazer"),  pkgTest)

这是我正在使用并尝试运行测试的数据：


data(infants)

imr <- with(infants, Surv(enter, exit, event))

cox <- coxph(imr ~ sex + age, data = infants)
summary(cox)
drop1(cox, test = "Chisq")
stargazer(cox, type = "text")


cox_fit <- survfit(cox)
autoplot(cox_fit)

newdat <- with(infants, 
               data.frame(
                 sex = c("male", "female"), age="Age"
               )
)

plot(survfit(cox, newdata = newdat), xscale = 12,
     conf.int = T,
     ylim = c(0.6, 1),
     col = c("red", "blue"),
     xlab = "Time",
     ylab = "Survival proportion",
     main = "")
legend("bottomleft",
       legend=c("Male", "Female"),
       lty = 1, 
       col = c("red", "blue"),
       text.col = c("red", "blue"))

# Adding an interaction
cox.int <- coxph(imr ~ sex * age, data = infants)
summary(cox.int)
drop1(cox.int, test = "Chisq")
stargazer(cox.int, type = "text")

plot(survfit(cox.int, newdata = newdat), xscale = 12,
      conf.int = T,
      ylim = c(0.6, 1),
      col = c("male", "female"),
      xlab = "Age",
      ylab = "Survival proportion",
      main = "")

对此的任何建议都将非常有帮助！我无法弄清楚我在根据最终数据生成图表时做错了什么。

Answer 1

要绘制性别对生存概率的影响，请根据 infants 中的性别级别重命名 newdat 中的性别级别。年龄固定为他们的平均值，重复不同性别水平的次数。

代码

data(infants)

imr <- with(infants, Surv(enter, exit, event))

cox <- coxph(imr ~ sex + age, data = infants)

cox_fit <- survfit(cox)

newdat <- with(infants, 
               data.frame(
                 sex = c("girl", "boy"), age=rep(mean(age, na.rm = TRUE), 2)))

plot(survfit(cox, newdata = newdat), xscale = 12,
     conf.int = T,
     ylim = c(0.6, 1),
     col = c("red", "blue"),
     xlab = "Time",
     ylab = "Survival proportion",
     main = "")
legend("bottomleft",
       legend=c("Male", "Female"),
       lty = 1, 
       col = c("red", "blue"),
       text.col = c("red", "blue"))

输出

我建议改用 ggsurvplot，这样可以进行更多自定义：

ggsurvplot(survfit(cox, newdata = newdat), data = infants,
       legend.labs = c("Male", "Female"),
       legend.title = "Sex")

数据解释和可视化

Data Interpretation and Visualisation

r

statistics

data-visualization