在R中使用boot包计算两组中位数差异的95%置信区间时返回NA结果的问题排查

你遇到的NA结果问题，核心原因是bootstrap函数里的分组取值写错了！

看一下你的数据结构：q2vsq1_baseline_cont_outcome这个因子变量的水平是"0"（对应低基线组Q1）和"1"（对应高基线组Q2），但你写的med.diff函数里却用了"2"和"1"来筛选分组——显然找不到"2"这个组，当bootstrap抽样时，自然会因为找不到对应组的观测，导致中位数计算返回NA，最终整个boot结果出现NA。

修正后的代码 #

首先，修正med.diff函数里的分组筛选条件，把"2"改成"1"，"1"改成"0"（对应Q2减Q1的中位数差异）：

library(boot)

# 修正后的中位数差异计算函数
med.diff <- function(d, i) {
  mydata_boot <- d[i,]
  # Q2（高基线组，水平"1"）的中位数减去Q1（低基线组，水平"0"）的中位数
  median(mydata_boot$change_continuous_outcome[mydata_boot$q2vsq1_baseline_cont_outcome == "1"]) - 
    median(mydata_boot$change_continuous_outcome[mydata_boot$q2vsq1_baseline_cont_outcome == "0"])
}

# 重新运行bootstrap
boot_result <- boot(data = mydata, statistic = med.diff, R = 1000)

# 查看bootstrap得到的中位数差异的中位数
median(boot_result$t)

# 计算95%百分位数置信区间
boot.ci(boot_result, type = "perc")

额外建议 #

下次写这类函数前，可以先检查因子变量的水平，避免取值错误：

levels(mydata$q2vsq1_baseline_cont_outcome)

这样能明确知道分组对应的取值是什么，减少这类低级错误~

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