Reporting heteroscedasticity-robust standard errors in modelsummary for panel data (plm)

Question

I encounter a problem with the reporting of heteroscedasticity-robust standard errors for a panel data regression (plm) in the output with modelsummary. For a simple linear regression the reporting works without problems, for a linear panel data regression not.

Attached is an example syntax for the Grunfeld dataset in the plm package.

data("Grunfeld")

lm <- lm(value ~
           capital + inv,
         data = Grunfeld
           )

summary(lm)
coeftest(lm, vcov =
           vcovHC(lm, type = "HC1"))

re <- plm(value ~
            capital + inv
          , data = Grunfeld
          , model = "random"
          , effect = "twoways"
)

summary(re)
coeftest(re, vcov =
           vcovHC(re, type = "HC1"))

Both outputs with heteroscedasticity-robust standard errors should look like this:

However, when I report the model outputs in modelsummary, only for the simple linear regression the heteroscedasticity-robust standard errors are reported correctly. Most disturbingly, the output still lists the correction to the standard errors for the plm regression in the last line.

modelsummary::msummary(list(
  "lm" = lm
  , "random effects" = re
  ),
  coef_map = c(
  "(Intercept)" = "Intercept"
  , "capital" = "Capital"
  , "inv" = "Investment"
  ),
  stars = TRUE,
  statistic = c(
    "p.value"
    , "conf.int"
    , "std.error"
    ),
  vcov = "HC1"
  ) 

Can anyone please help me to fix this error? Thank you very much for your support.

Answer 1

Just replace the vcov.

re$vcov <- plm::vcovHC(re, type = "HC1")

modelsummary::msummary(list(lm = lm, `random effects` = re), 
    coef_map = c(`(Intercept)` = "Intercept", capital = "Capital", 
        inv = "Investment"), stars = TRUE, statistic = c("p.value", 
        "conf.int", "std.error"))