如何为变量的每个四分位数拟合生存模型？

Question

我尝试为数据集中一个变量的每个四分位数拟合一个生存模型。以survival包中的肺癌数据集为例

library(survival)
datalung <- lung
attach(datalung)
fit<- survfit(Surv(time,status) ~ quantile(age)[2:5],type = "kaplan-meier")

但是我得到一个关于变量长度不同的错误。我想拟合四个模型，每个四分位数一个。

小于或等于 25%

大于 25 且小于等于 50%

大于 50% 小于等于 75%

大于 75%

我该怎么做？

Answer 1

默认情况下，quantile 返回 5 个值，位于 prob = seq(0, 1, 0.25)。我认为您想使用cut 来获取因子变量：

library(survival)
datalung <- lung
datalung$fage <- with(datalung, cut(age, quantile(age), include = TRUE))

## don't use `attach()`; use the `data` argument of model fitting routine
fit <- survfit(Surv(time,status) ~ fage, data = datalung, type="kaplan-meier")

#Call: survfit(formula = Surv(time, status) ~ fage, data = datalung, 
#    type = "kaplan-meier")
#
#              n events median 0.95LCL 0.95UCL
#fage=[39,56] 58     39    337     239     457
#fage=(56,63] 59     41    348     245     574
#fage=(63,69] 55     39    329     285     477
#fage=(69,82] 56     46    283     222     361

---

跟进

@42-也用过quantile，得到的是“左闭右开”区间。

你的问题是：

• 小于或等于 25%
• 大于 25 且小于等于 50%
• 大于 50% 小于等于 75%
• 大于 75%

很明显你想要“左开右闭”的间隔。因此，我的代码正是你想要的。

What is the meaning of include.lowest in reclassify raster package 详细解释了 include.lowest 和 right 内部的 cut 和 raster::reclassify 参数。现在让我们比较一下：

## my factor
table(with(datalung, cut(age, quantile(age), include.lowest = TRUE)))
#[39,56] (56,63] (63,69] (69,82] 
#     58      59      55      56 

## 42-'s factor
table(with(datalung, cut(age, quantile(age), include.lowest = TRUE, right = FALSE)))
#[39,56) [56,63) [63,69) [69,82] 
#     49      57      55      67 

Answer 2

我尝试使用我喜欢的方法来创建四分位指标：

library(survival)
datalung <- lung
detach(datalung)  # Agree with Zheyuan Li that attach()-ing is dangerous practice.
fit3<- survfit(Surv(time,status) ~ findInterval(age, quantile(age)[-5]), 
                     data=datalung, type = "kaplan-meier")

需要删除向量中的第五个项目是拆分值，因为 findInterval 的拆分在左侧关闭，并且会得到第五个只有最大年龄的组。请注意，我们的四分位数结果是不同的。他的方法丢失了案例，而不仅仅是最小或最大组。他们去了哪里……我还不确定：

> fit3
Call: survfit(formula = Surv(time, status) ~ findInterval(age, quantile(age)[-5]), 
    data = datalung, type = "kaplan-meier")

                                        n events median 0.95LCL 0.95UCL
findInterval(age, quantile(age)[-5])=1 49     32    320     226     533
findInterval(age, quantile(age)[-5])=2 57     41    340     245     433
findInterval(age, quantile(age)[-5])=3 55     39    310     267     524
findInterval(age, quantile(age)[-5])=4 67     53    285     229     363

您向李哲元提出的关于 ggplot 中级别顺序的问题暴露了使用 cut 的另一个陷阱，至少如果不提供带有“标签”参数的名称的话。级别是按词法排序的，“[”是>而不是“（”：

> levels(datalung$fage)
[1] "[39,56]" "(56,63]" "(63,69]" "(69,82]"
> "[" < "("
[1] FALSE

要解决我使用分位数与@ZheyuanLi 使用以及他对我的方法的错误描述的问题，只需检查一下：

> quantile(datalung$age)
  0%  25%  50%  75% 100% 
  39   56   63   69   82 

> with( datalung, table( findInterval(age, quantile(datalung$age)[-5] )))

 1  2  3  4 
49 57 55 67 

所以大部分区别在于 56 岁的处理方式：

>  sum(lung$age==56)
[1] 9

在使用 cut() 时试图解决标签问题（这不是我的责任，对吗？）：

> library(ggplot2)  # checked to make sure I have the most recent version per CRAN
> autoplot(fit2)
Error: Objects of type survfit not supported by autoplot.