多重插补(MI)后分类变量比例的合并标准误获取方法咨询

你当前直接合并所有插补数据集计算均值的方法，没有考虑不同插补数据集之间的变异，得到的结果是不准确的。可以用pool.scalar实现正确的合并标准误计算，因为多重插补的合并需要同时考虑插补内方差和插补间方差两部分，具体操作如下：

核心思路 #

多重插补合并估计的总方差公式为：
t = W + (1 + 1/m)*B
其中：

• W：m个插补数据集内方差的均值
• B：m个插补数据集估计值的插补间方差
• m：插补数据集数量
  合并后的标准误就是总方差t的平方根，而pool.scalar可以自动完成这些计算。

具体代码实现 #

1. 对每个插补数据集计算比例及单个数据集内的方差 #

用with()函数在mice插补对象上批量操作，计算每个插补数据集中分类变量各水平的比例，同时算出单个数据集内比例的方差（比例方差公式为p*(1-p)/n，n为单个数据集样本量）：

# 加载所需包
library(mice)
library(dplyr)

# 你的初始代码
data1 <- nhanes2
imp.data <- mice(data = data1, m = 5, maxit = 10, seed = 12345, method = "cart")

# 对每个插补数据集计算hyp各水平的比例及对应方差
imp_results <- with(imp.data, expr = {
  tab <- table(hyp)
  props <- prop.table(tab)
  n <- sum(tab)
  vars <- props * (1 - props) / n
  list(props = props, vars = vars)
})

2. 提取各水平的估计值和方差，用pool.scalar合并 #

针对每个分类水平，提取m个插补数据集的比例估计值（Q）和对应方差（U），传入pool.scalar计算合并结果：

# 处理水平0
q0 <- sapply(imp_results$analyses, function(x) x$props["0"])
u0 <- sapply(imp_results$analyses, function(x) x$vars["0"])
pool0 <- pool.scalar(Q = q0, U = u0)
cat("水平0的合并比例：", pool0$qbar, "\n")
cat("水平0的合并标准误：", sqrt(pool0$t), "\n")

# 处理水平1
q1 <- sapply(imp_results$analyses, function(x) x$props["1"])
u1 <- sapply(imp_results$analyses, function(x) x$vars["1"])
pool1 <- pool.scalar(Q = q1, U = u1)

# 处理水平2
q2 <- sapply(imp_results$analyses, function(x) x$props["2"])
u2 <- sapply(imp_results$analyses, function(x) x$vars["2"])
pool2 <- pool.scalar(Q = q2, U = u2)

结果说明 #

• pool.scalar返回的qbar是m个插补数据集比例的均值，即合并后的比例估计值
• t是合并后的总方差，对其取平方根就是你需要的合并标准误
• 这种方式同时考虑了插补内和插补间的变异，比直接合并数据集的结果更可靠

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