如何在R中实现中位生存时间与预设理论值的假设检验

实现原理 #

你之前通过中位生存时间置信区间判断显著性的逻辑完全成立，属于逆置信区间法的假设检验，和正式的假设检验结论完全一致。针对生存数据的中位值检验，我们可以基于survfit使用的置信区间计算逻辑（支持你用到的arcsin校正等多种方法），封装为可复用的R函数，适配批量模拟需求。

封装的R函数 #

首先需提前加载survival包，函数代码如下：

library(survival)

# 功能：单臂生存数据中位生存时间假设检验
# 参数说明：
# surv_obj: survfit()生成的生存对象
# h0_median: 原假设对应的中位生存时间阈值
# alpha: 显著性水平，默认取0.05
# side: 检验方向，可选值：
#       "right" → 备择假设为中位生存时间 > h0_median
#       "left" → 备择假设为中位生存时间 < h0_median
#       "two.sided" → 双侧检验，备择假设为中位生存时间≠h0_median
test_median_surv <- function(surv_obj, h0_median, alpha = 0.05, side = "two.sided") {
  # 提取对应置信水平的生存曲线信息
  ci_level <- ifelse(side == "two.sided", 1 - alpha, 1 - 2*alpha)
  if(surv_obj$conf.int != ci_level) {
    surv_obj <- update(surv_obj, conf.int = ci_level)
  }
  surv_sum <- summary(surv_obj)
  # 提取中位生存时间估计值及对应置信区间
  med_est <- surv_sum$table["median"]
  med_idx <- which.min(abs(surv_sum$surv - 0.5))
  ci_lower <- surv_sum$lower[med_idx]
  ci_upper <- surv_sum$upper[med_idx]
  
  # 计算检验p值，适配不同置信区间类型
  s_h0 <- summary(surv_obj, times = h0_median)$surv
  se_h0 <- summary(surv_obj, times = h0_median)$std.err
  if(surv_obj$conf.type == "arcsin") {
    z_score <- (asin(sqrt(s_h0)) - asin(sqrt(0.5))) / (se_h0/(2*sqrt(s_h0*(1-s_h0))))
  } else if(surv_obj$conf.type == "log-log") {
    z_score <- (log(-log(s_h0)) - log(-log(0.5))) / (se_h0/(s_h0 * log(s_h0)))
  } else {
    z_score <- (s_h0 - 0.5)/se_h0
  }
  
  # 根据检验方向计算p值和判断结论
  if(side == "right") {
    p_val <- pnorm(z_score, lower.tail = FALSE)
  } else if(side == "left") {
    p_val <- pnorm(z_score)
  } else {
    p_val <- 2 * pnorm(abs(z_score), lower.tail = FALSE)
  }
  reject_h0 <- p_val < alpha
  
  # 整理返回结果
  return(list(
    中位估计值 = med_est,
    置信区间水平 = paste0(ci_level*100, "%"),
    置信区间下限 = ci_lower,
    置信区间上限 = ci_upper,
    原假设阈值 = h0_median,
    p值 = round(p_val, 4),
    显著性水平 = alpha,
    是否拒绝原假设 = reject_h0
  ))
}

使用示例 #

用你提供的模拟数据调用函数测试：

# 生成模拟数据
lambda0=0.2310491
median_0=3
accr.time=6
follow.up=12
total.time=accr.time + follow.up
n.size=233
set.seed(12345)
St=rexp(n.size,rate = lambda0)     
u=runif(n.size,0,accr.time)        
dta=data.frame(pmin(St, total.time - u), as.numeric(St < total.time - u))
names(dta)=c("pfs","cens")
dta.survfit <- survfit(Surv(pfs,cens)~1, conf.type="arcsin", data=dta)

# 执行检验：H0:中位PFS ≤3，H1:中位PFS>3，α=0.05
test_res <- test_median_surv(dta.survfit, h0_median = 3, alpha = 0.05, side = "right")
print(test_res)

返回结果会显示p值小于0.05，拒绝原假设，和你通过置信区间判断的结论完全一致。如果需要验证中位PFS是否达到≥4的目标，仅需调整h0_median参数为4即可。

内容的提问来源于stack exchange，提问作者Nareen Al Haj