R中OLS与固定效应模型的稳健/聚类标准误使用正确性咨询

我正在R Studio中估计无固定效应的OLS回归和含固定效应的OLS回归。了解到简单OLS模型通常使用稳健标准误，固定效应模型通常使用聚类标准误，但不确定自己的实现是否正确，恳请指导！我的数据为paneldatafinal，是包含IND（行业）和year（年份）的面板数据。

一、带稳健标准误的OLS模型 #

原代码 #

# MODEL EPS MARKET BASED:

model_mb_basic_ols <- lm(diff_log_VATFP_I ~ log_Impen_Mil + RuDintensity_vapc +
                           dummy_var*EPS_MKT_growth_MA3 + log_PRDK_growth_rate + GAP, 
                         data = paneldatafinal)


# Compute robust standard errors:
robust_se <- vcovHC(model_mb_basic, type = "HC1")

# Display the results:
summary_result <- coeftest(model_mb_basic, vcov = robust_se)

# Print summary result with robust standard errors:
print(summary_result)

代码修正与说明 #

• 变量名错误：原代码中vcovHC和coeftest调用的模型对象是model_mb_basic，但实际拟合的模型是model_mb_basic_ols，需统一变量名：

  robust_se <- vcovHC(model_mb_basic_ols, type = "HC1")
  summary_result <- coeftest(model_mb_basic_ols, vcov = robust_se)
  

• 稳健标准误选择：type = "HC1"是合理的，它经过小样本自由度调整，是应用最广泛的稳健标准误类型之一，适配OLS模型的异方差修正需求。

二、带聚类标准误的固定效应模型 #

原代码 #

## Model with clustered Standard Errors & fixed effects EPS TOTAL:

model_basicfe <- plm(diff_log_VATFP_I ~ log_Impen_Mil + RuDintensity_vapc + 
                       dummy_var*EPS_growth_MA3 + log_PRDK_growth_rate + GAP,
                    data = paneldatafinal, model = "within", index = c("IND", "year"))

# Calculate clustered standard errors:
vcov_clustered <- vcovHC(model_basicfe, type = "HC1", cluster = "group")

# Apply the clustered standard errors to the model:
coeftest(model_basicfe, vcov = vcov_clustered)

# Use stargazer to present the results:
stargazer(model_basicfe, type = "text", se = list(sqrt(diag(vcov_clustered))), header = FALSE)

代码验证与优化建议 #

• 固定效应模型拟合：仅指定model = "within"时，plm默认只吸收IND（行业）固定效应。若需同时控制年份固定效应，可二选一：

  # 方式1：公式中加入年份固定效应
  model_basicfe <- plm(diff_log_VATFP_I ~ log_Impen_Mil + RuDintensity_vapc + 
                         dummy_var*EPS_growth_MA3 + log_PRDK_growth_rate + GAP + factor(year),
                      data = paneldatafinal, model = "within", index = c("IND", "year"))
  # 方式2：使用双向固定效应参数
  model_basicfe <- plm(diff_log_VATFP_I ~ log_Impen_Mil + RuDintensity_vapc + 
                         dummy_var*EPS_growth_MA3 + log_PRDK_growth_rate + GAP,
                      data = paneldatafinal, model = "within", index = c("IND", "year"), effect = "twoways")
  

• 聚类标准误设置：cluster = "group"对应index中的第一个维度（IND），即按行业聚类，符合面板数据聚类标准误的常规逻辑（控制行业内组内相关）。也可直接写成cluster = "IND"，效果一致。
• Stargazer输出优化：原代码写法正确，也可通过vcov参数直接传入聚类方差协方差矩阵简化代码：

  stargazer(model_basicfe, type = "text", vcov = vcov_clustered, header = FALSE)
  

内容的提问来源于stack exchange，提问作者Lucas van der List