как создать один результат с несколькими моделями линейной регрессии, включая кластерные стандартные ошибки в Stargazer

Я знаю, как создать модель линейной регрессии lm и как использовать функцию summary для получения кластеризованных стандартных ошибок и добавления их к выводу stargazer:

# estimate models
ols1 <- lm(y ~ x) 

# summary with cluster-robust SEs
summary(ols1, cluster="cluster_id") 

# create table in stargazer
stargazer(ols1, se=list(coef(summary(ols1,cluster = c("cluster_id")))[, 2]), type = "text")

Кто-нибудь знает, как должен выглядеть код, если я хочу создать один вывод звездочета с несколькими моделями регрессии и кластеризованными стандартными ошибками?

Логика кода следующая:

Шаг 1: создание моделей lm

ols1 <- lm(y ~ x) 
ols2 <- lm(y ~ x + z) 
ols3 <- lm(y ~ x + z + a) 
ols2 <- lm(y ~ x + z + a + b) 

Шаг 2: включите стандартные ошибки

summary(ols1, cluster="cluster_id")
summary(ols2, cluster="cluster_id")
summary(ols3, cluster="cluster_id")
summary(ols4, cluster="cluster_id")

Шаг 3: создайте один результат с 4 разными моделями

stargazer(ols1,ols2,ols3,ols4, type="html", dep.var.labels=c("ROA"), intercept.bottom = FALSE,
          out="OLS1")

Я думаю, что шаг 1 и шаг 2 не критичны, но я не знаю, как настроить код для шага 3.

Я не знаю, как реализовать следующий код на шаге 3:

# create table in stargazer
stargazer(ols1, se=list(coef(summary(ols1,cluster = c("cluster_id")))[, 2]), type = "text")

Кто-нибудь может помочь?

Большое спасибо!!!


person Yaron Nolan    schedule 28.12.2019    source источник


Ответы (1)


Вы уже видели, что опция se может обрабатывать list. Так что просто сохраните надежный SE в другом месте и легко доступном, например:

library("sandwich")
library("plm")
library("stargazer")

data("Produc", package = "plm")

# Regression    
model1 <- plm(log(gsp) ~ log(pcap) + log(pc) + log(emp) + unemp,
             data = Produc, 
             index = c("state","year"),
             method="pooling")

model2 <- plm(log(gsp) ~ log(pcap) + log(pc) + log(emp),
              data = Produc, 
              index = c("state","year"),
              method="pooling")

# Adjust standard errors
cov1         <- vcovHC(model1, type = "HC1")
cov2         <- vcovHC(model2, type = "HC1")
robust_se1    <- sqrt(diag(cov1))
robust_se2    <- sqrt(diag(cov2))

# Stargazer output (with and without RSE)
stargazer(model1, model2, type = "text",
          se = list(robust_se1, robust_se2))

Теперь у вас может быть как любые модели рядом друг с другом с надежным SE в stargazer.

person Marco    schedule 14.01.2020