Fine-Gray与病因特异性Cox回归：R中sHR和csHR的代码正确性验证

数据集编码规则 #

• status = 1：关注结局
• status = 2：竞争事件
• status = 0：截尾（未发生任一事件仍存活）

我的模型代码 #

Fine-Gray模型（亚分布风险模型） #

library(cmprsk)
library(tidycmprsk)

# 关注结局的亚分布风险比（sHR）
crr(Surv(followup, status) ~ cont1 + cont2 + cont3 + categ1 + categ2 + categ3, 
    data = d, failcode = 1)

# 竞争事件的亚分布风险比（sHR）
crr(Surv(followup, status) ~ cont1 + cont2 + cont3 + categ1 + categ2 + categ3, 
    data = d, failcode = 2)

我的理解：

• failcode = 1 输出关注结局的sHR
• failcode = 2 输出竞争事件的sHR

病因特异性Cox回归模型 #

library(survival)
# 关注结局的病因特异性风险比（csHR）
coxph(Surv(followup, status == 1) ~ cont1 + cont2 + cont3 + categ1 + categ2 + categ3, data = d)

# 竞争事件的病因特异性风险比（csHR）
coxph(Surv(followup, status == 2) ~ cont1 + cont2 + cont3 + categ1 + categ2 + categ3, data = d)

我的理解：非指定编码的事件会被视为截尾。

问题解答 #

1. tidycmprsk::crr()的failcode参数设置是否正确？ #

完全正确。tidycmprsk::crr()底层依赖cmprsk::crr()，failcode参数明确指定了要建模的目标事件：

• failcode = 1 对应关注结局的亚分布风险比（sHR），此时status=2的竞争事件和status=0的截尾数据都会被纳入风险集处理
• failcode = 2 对应竞争事件的亚分布风险比（sHR），此时status=1的关注结局会被当作风险集中的截尾情况处理

2. 病因特异性Cox回归的代码设置是否正确？ #

你的核心逻辑是对的：通过status == 1/status == 2定义目标事件，将其他事件（包括竞争事件和原始截尾）统一视为截尾，符合病因特异性模型的假设——只针对单一事件类型建模，忽略其他事件的竞争风险。

不过更规范的写法是直接用subset参数筛选目标事件和截尾数据，避免逻辑判断的隐式转换：

# 关注结局的csHR
coxph(Surv(followup, status) ~ cont1 + cont2 + cont3 + categ1 + categ2 + categ3, 
      data = d, subset = status %in% c(0,1))

# 竞争事件的csHR
coxph(Surv(followup, status) ~ cont1 + cont2 + cont3 + categ1 + categ2 + categ3, 
      data = d, subset = status %in% c(0,2))

或者用event参数明确指定：

coxph(Surv(followup, event = (status == 1)) ~ cont1 + cont2 + cont3 + categ1 + categ2 + categ3, data = d)

两种写法效果一致，但后两种可读性更强，也减少了潜在的类型转换错误。

3. sHR和csHR结果并列展示的最佳实践 #

推荐用broom包提取模型结果，再通过dplyr合并，最后用表格工具可视化，步骤如下：

步骤1：提取并整理模型结果

library(broom)
library(dplyr)

# 提取Fine-Gray模型结果
fg_interest <- tidy(crr(Surv(followup, status) ~ cont1 + cont2 + cont3 + categ1 + categ2 + categ3, 
                        data = d, failcode = 1), exponentiate = TRUE) %>%
  mutate(model_type = "Fine-Gray", target_event = "关注结局", hr_type = "sHR")

fg_competing <- tidy(crr(Surv(followup, status) ~ cont1 + cont2 + cont3 + categ1 + categ2 + categ3, 
                         data = d, failcode = 2), exponentiate = TRUE) %>%
  mutate(model_type = "Fine-Gray", target_event = "竞争事件", hr_type = "sHR")

# 提取病因特异性Cox模型结果
cs_interest <- tidy(coxph(Surv(followup, status == 1) ~ cont1 + cont2 + cont3 + categ1 + categ2 + categ3, 
                          data = d), exponentiate = TRUE) %>%
  mutate(model_type = "病因特异性Cox", target_event = "关注结局", hr_type = "csHR")

cs_competing <- tidy(coxph(Surv(followup, status == 2) ~ cont1 + cont2 + cont3 + categ1 + categ2 + categ3, 
                           data = d), exponentiate = TRUE) %>%
  mutate(model_type = "病因特异性Cox", target_event = "竞争事件", hr_type = "csHR")

# 合并所有结果
combined_results <- bind_rows(fg_interest, fg_competing, cs_interest, cs_competing) %>%
  select(term, target_event, hr_type, model_type, estimate, conf.low, conf.high, p.value)

步骤2：生成对比表格

用gt包生成美观的可交互表格：

library(gt)
combined_results %>%
  gt() %>%
  tab_header(title = "竞争风险模型结果对比") %>%
  cols_label(
    term = "协变量",
    target_event = "目标事件",
    hr_type = "风险比类型",
    model_type = "模型类型",
    estimate = "HR值",
    conf.low = "95%CI下限",
    conf.high = "95%CI上限",
    p.value = "P值"
  ) %>%
  fmt_number(columns = c(estimate, conf.low, conf.high), decimals = 2) %>%
  fmt_pvalue(columns = p.value, decimals = 3) %>%
  tab_style(style = cell_text(weight = "bold"), locations = cells_column_labels())

如果需要按变量分组对比sHR和csHR，可以用pivot_wider重塑数据：

combined_results %>%
  select(term, target_event, hr_type, estimate, conf.low, conf.high) %>%
  pivot_wider(names_from = hr_type, values_from = c(estimate, conf.low, conf.high)) %>%
  gt() %>%
  tab_header(title = "同一事件的sHR与csHR对比")

内容的提问来源于stack exchange，提问作者Konstantinos Gkirgkiris