How to have output from lm() include std. error and others without using summary() for stargazer

Question

I'm fitting several linear models in r in the following way:

set.seed(12345)
n = 100
x1 = rnorm(n)
x2 = rnorm(n)+0.1
y = x + rnorm(n)

df <- data.frame(x1, x2, y)
x_str <- c("x1", "x1+x2")

regf_lm <- function(df,y_var, x_str ) {
  frmla <- formula(paste0(y_var," ~ ", x_str ))
  fit <- lm(frmla, data = df ) 
  summary(fit) #fit
}

gbind_lm <- function(vv) {
  n <-  vv %>% length()
  fits <- list()
  coefs <- list()
  ses <- list()
  for (i in 1:n ) {
    coefs[[i]] <- vv[[i]]$coefficients[,1]
    ses[[i]] <- vv[[i]]$coefficients[,2]
    fits[[i]] <- vv[[i]]
  }
  
  list("fits" = fits, "coefs" = coefs, "ses" = ses)
}
  
stargazer_lm <- function(mylist, fname, title_str,m_type = "html",...) {
    stargazer(mylist$fits, coef =  mylist$coefs,
              se = mylist$ses,
              type = m_type, title = title_str, 
              out = paste0("~/projects/outputs",fname),  single.row = T ,...)
}

p_2 <- map(x_str, 
           ~ regf_lm (df = df ,
                      y_var = "y", x_str = .))
m_all <- do.call(c, list(p_2)) %>% gbind_lm()
stargazer_lm(m_all,"name.html","My model", m_type = "html")

In regf_lm, if I use summary(fit) on the last line, I'm able to generate reg output with columns for estimated coefficients, std. error, etc. But Stargazer() does not work with summary(lm()) (returns error $ operator is invalid for atomic vectors). However, if I just use "fit" on the last line in regf_lm, the output shows only the estimated coefficients and not std error, R sq...and gbind_lm() won't work because I cannot extract ses or fit.

Any advice is greatly appreciated.

Answer 1

You can directly export model statistics in tidy format with the package broom

library(broom)

set.seed(12345)
n = 100
x1 = rnorm(n)
x2 = rnorm(n)+0.1
y = x1 + rnorm(n)

df <- data.frame(x1, x2, y)
x_str <- c("x1", "x1+x2")

regf_lm <- function(df,y_var, x_str ) {
  frmla <- formula(paste0(y_var," ~ ", x_str ))
  fit <- lm(frmla, data = df ) 
  return(list(fit,select(broom::tidy(fit),std.error))) #fit
}

exm_model <- regf_lm(iris,'Sepal.Width','Sepal.Length')

stargazer(exm_model[[1]], coef =  exm_model[[2]], title = 'x_model', 
              out ='abc',  single.row = T)

This piece of code worked on my local with no problem, I think you can apply this in your workflow.