Stargazer输出lme4包lmer/glmer聚类稳健标准误方法咨询

核心逻辑：lme4包拟合的混合效应模型默认输出无聚类调整的原生标准误，需先手动计算对应维度的聚类稳健方差，再将自定义统计量传入stargazer覆盖默认输出即可，无需修改模型拟合过程。

前置依赖包 #

需提前安装加载以下工具包：

• lme4：拟合lmer、glmer混合效应模型
• clubSandwich：计算适配混合效应模型的聚类稳健方差，小样本偏差远小于普通三明治估计
• stargazer：格式化输出回归结果表

1. lmer模型（线性混合效应）聚类标准误计算 #

拟合纳入个体随机效应的基准模型后，直接按小组维度计算聚类稳健标准误，示例代码如下：

# 加载包
library(lme4)
library(clubSandwich)
library(stargazer)

# 拟合基准lmer模型：y为因变量，treat为核心处理变量，x1/x2为控制变量，(1|individual_id)为个体层面随机截距
lmer_fit <- lmer(y ~ treat + x1 + x2 + (1 | individual_id), data = df)

# 计算小组维度的CR2类型聚类稳健方差（与Stata默认聚类SE校正逻辑一致，适配小聚类场景）
lmer_clu_vcov <- vcovCR(lmer_fit, cluster = df$group, type = "CR2")
# 提取聚类调整后标准误
lmer_clu_se <- sqrt(diag(lmer_clu_vcov))
# 提取固定效应系数
lmer_coef <- fixef(lmer_fit)
# 计算Satterthwaite自由度校正后的聚类调整p值
lmer_clu_p <- coef_test(lmer_fit, vcov = lmer_clu_vcov, test = "Satterthwaite")$p_Satt

2. glmer模型（广义线性混合效应）聚类标准误计算 #

逻辑与lmer完全一致，广义模型可直接用正态近似计算检验p值，示例代码如下：

# 拟合基准glmer模型，示例为二值因变量的logit回归，可按需替换family参数适配泊松、probit等模型
glmer_fit <- glmer(y_binary ~ treat + x1 + x2 + (1 | individual_id), 
                   family = binomial(link = "logit"), 
                   data = df)

# 计算小组维度聚类稳健方差
glmer_clu_vcov <- vcovCR(glmer_fit, cluster = df$group, type = "CR2")
# 提取聚类调整后标准误
glmer_clu_se <- sqrt(diag(glmer_clu_vcov))
# 提取固定效应系数
glmer_coef <- fixef(glmer_fit)
# 计算聚类调整z值与p值
glmer_clu_z <- glmer_coef / glmer_clu_se
glmer_clu_p <- 2 * pnorm(abs(glmer_clu_z), lower.tail = FALSE)

3. stargazer结果输出 #

将计算好的自定义系数、标准误、p值按模型顺序传入stargazer即可，示例代码如下：

stargazer(
  lmer_fit, glmer_fit,
  title = "混合效应模型回归结果",
  # 按模型顺序传入聚类调整后的统计量，覆盖默认输出
  coef = list(lmer_coef, glmer_coef),
  se = list(lmer_clu_se, glmer_clu_se),
  p = list(lmer_clu_p, glmer_clu_p),
  # 输出设置：type可选text/latex/html，out参数指定导出文件路径
  type = "text",
  out = "cluster_se_reg.html",
  star.cutoffs = c(0.1, 0.05, 0.01),
  notes = "标准误按分配小组维度聚类调整",
  notes.append = TRUE
)

注意事项 #

• 不建议用sandwich包的普通三明治方差函数计算混合效应模型聚类SE，对不平衡样本、聚类数较少的场景偏差较大，结果与主流计量软件输出一致性差。
• 若聚类数量少于20，建议在coef_test中调整小样本校正参数，避免p值估计偏误。
• 传入stargazer的系数、标准误顺序需与模型变量顺序完全一致，fixef()与vcovCR()提取的变量顺序默认匹配，无需手动调整。

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