如何在R中为plm对象获取非面板时间/组ID层级的多维度聚类标准误？

没错，plm包里的vcovHC确实只能按面板自带的个体（组）或时间ID聚类，要是咱们需要在其他层级——比如行业、地区这种跨面板单位的维度，或者同时在多个维度聚类，其实可以通过sandwich和lmtest这两个包来搞定，甚至能自定义聚类变量。下面给你一步步拆解操作方法：

1. 单维度非面板层级聚类 #

假设你的数据里有一个额外的聚类变量（比如industry行业），不是面板的个体或时间ID，操作如下：

• 先拟合你的面板模型：

library(plm)
library(sandwich)
library(lmtest)

# 用plm自带的Grunfeld数据做示例
data("Grunfeld")
# 模拟一个非面板的行业聚类变量
Grunfeld$industry <- rep(c(1,2), length.out = nrow(Grunfeld))

# 拟合随机效应面板模型（你也可以换成固定效应等其他模型）
model <- plm(inv ~ value + capital, data = Grunfeld, model = "random", index = c("firm", "year"))

• 用sandwich的vcovCL函数计算自定义聚类的稳健标准误：

# 指定按industry变量聚类，标准误类型选HC1（小样本表现更稳定）
vcov_custom <- vcovCL(model, cluster = ~ industry, type = "HC1")

• 最后用lmtest的coeftest输出带聚类标准误的回归结果：

coeftest(model, vcov = vcov_custom)

2. 多维度聚类（含非面板维度） #

如果需要同时在多个维度聚类——比如同时按行业和时间，或者个体+行业这种混合维度，vcovCL同样支持多变量聚类：

# 同时按industry和year两个维度聚类
vcov_multi <- vcovCL(model, cluster = ~ industry + year, type = "HC1")
# 输出结果
coeftest(model, vcov = vcov_multi)

3. 几个关键注意点 #

• 聚类变量可以是因子、数值型的分组标识，vcovCL会自动识别分组，不需要额外转换（但如果是字符型，建议先转成因子）。
• 如果拟合的是固定效应模型，plm已经吸收了个体/时间固定效应，聚类时不用重复控制，但多维度聚类依然适用。
• 标准误类型可以调整：sandwich支持HC1、HC2、HC3等多种类型，默认是HC3，HC1在小样本中更常用，你可以根据需求替换type参数。

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