Betareg in Stargazer (Error to allocate vector of size __)

Question

betareg default residuals are heavy, which may cause an error to allocate vector due to its high size. This can be solved by changing the type of residuals in the summary call as explained here.

However, when presenting regression tables with stargazer, the type of residuals can't be explicitly set.

Is there any way to make (large) betareg objects work in stargazer?

The potential solutions that I can think of, but don't know how to implement, are:

• Indicating the residual type in the original betareg call (type = "pearson" (or any other type) doesn't work).
• Explicitly indicating the argument that stargazer should include when calling summary on the betareg object.
• Any other?

Example:

set.seed(1)
df <- data.frame(y=runif(100000), x=runif(100000))

library(betareg)
beta <- betareg(y ~ x, data=df)

library(stargazer)
stargazer(beta)
# Error: cannot allocate vector of size 74.5 Gb

Answer 1

The way I finally solved this issue was by modifying the stargazer function. Particularly, in .stargazer.wrap, which code you can find here, I added the following code:

else if (class(object.name)[1] == "betareg") { #betareg
          .summary.object <<- summary(object.name, type = "person")

inside .add.model, i.e. just after these lines of code:

.add.model <-
      function(object.name, user.coef=NULL, user.se=NULL, user.t=NULL, user.p=NULL, auto.t=TRUE, auto.p=TRUE, user.ci.lb=NULL, user.ci.rb=NULL) {
        
        if (class(object.name)[1] == "Glm") {
          .summary.object <<- summary.glm(object.name)
        }

Note that you need to create a new stargazer function with the new .stargazer.wrap. I just called it stargazer2:

stargazer2 <-
  function(..., type = "latex", title="", style="default", summary=NULL, out=NULL, out.header=FALSE,
           column.labels=NULL, column.separate = NULL, covariate.labels=NULL, dep.var.caption=NULL, 
           dep.var.labels=NULL, dep.var.labels.include=TRUE, align=FALSE, coef=NULL, se=NULL, t=NULL, 
           p=NULL,  t.auto=TRUE, p.auto=TRUE, ci=FALSE, ci.custom=NULL, ci.level=0.95, ci.separator=NULL, 
           add.lines=NULL, apply.coef=NULL, apply.se=NULL, apply.t=NULL, apply.p=NULL, apply.ci=NULL, 
           colnames = NULL, column.sep.width = "5pt", 
           decimal.mark=NULL, df=TRUE, digit.separate=NULL, digit.separator=NULL, digits=NULL, digits.extra=NULL, 
           flip=FALSE,
           float=TRUE, float.env="table", font.size=NULL, header=TRUE, initial.zero=NULL, intercept.bottom=TRUE, 
           intercept.top=FALSE, keep=NULL, keep.stat=NULL, label="", model.names=NULL, model.numbers=NULL, 
           multicolumn=TRUE, no.space=NULL, notes=NULL, notes.align=NULL, notes.append=TRUE, notes.label=NULL, 
           object.names=FALSE,
           omit=NULL, omit.labels=NULL, omit.stat=NULL, omit.summary.stat=NULL, omit.table.layout=NULL,
           omit.yes.no=c("Yes","No"), order=NULL, ord.intercepts=FALSE, 
           perl=FALSE, report=NULL, rownames = NULL,
           rq.se = "nid", selection.equation=FALSE, single.row=FALSE, star.char=NULL, 
           star.cutoffs=NULL, suppress.errors=FALSE, table.layout=NULL, table.placement = "!htbp", zero.component=FALSE, 
           summary.logical=TRUE, summary.stat=NULL, nobs=TRUE, mean.sd=TRUE, min.max=TRUE, median=FALSE, 
           iqr=FALSE) {
    
    save.warn.option <- getOption("warn") 
    options(warn=-1)
    return(.stargazer.wrap(..., type=type, title=title, style=style, summary=summary, out=out, out.header=out.header,
                           column.labels=column.labels, column.separate = column.separate,
                           covariate.labels=covariate.labels, dep.var.caption = dep.var.caption,
                           dep.var.labels=dep.var.labels, dep.var.labels.include=dep.var.labels.include,
                           align=align, coef=coef, se=se, t=t, p=p, t.auto=t.auto, p.auto=p.auto, 
                           ci=ci, ci.custom=ci.custom, ci.level=ci.level, ci.separator = ci.separator,
                           add.lines=add.lines, apply.coef=apply.coef, apply.se=apply.se, apply.t=apply.t, 
                           apply.p=apply.p, apply.ci=apply.ci, colnames=colnames,
                           column.sep.width=column.sep.width, decimal.mark=decimal.mark, df=df,
                           digit.separate=digit.separate, digit.separator=digit.separator, 
                           digits=digits, digits.extra=digits.extra, 
                           flip=flip,
                           float=float, float.env=float.env,
                           font.size=font.size, header=header, 
                           initial.zero=initial.zero, 
                           intercept.bottom=intercept.bottom, intercept.top=intercept.top, 
                           keep = keep, keep.stat = keep.stat,
                           label = label,
                           model.names=model.names, model.numbers=model.numbers,
                           multicolumn = multicolumn,
                           no.space=no.space, notes=notes, notes.align=notes.align, 
                           notes.append=notes.append, notes.label=notes.label, object.names=object.names,
                           omit=omit, omit.labels=omit.labels, omit.stat=omit.stat, omit.summary.stat=omit.summary.stat,
                           omit.table.layout=omit.table.layout,omit.yes.no=omit.yes.no, 
                           order=order, ord.intercepts=ord.intercepts, perl=perl, 
                           report=report, rownames=rownames, rq.se=rq.se, selection.equation=selection.equation,
                           single.row=single.row, star.char=star.char, 
                           star.cutoffs=star.cutoffs, suppress.errors=suppress.errors,
                           table.layout=table.layout,
                           table.placement = table.placement, zero.component=zero.component,
                           summary.logical = summary.logical,
                           summary.stat = summary.stat,
                           nobs=nobs, mean.sd=mean.sd, 
                           min.max=min.max, median=median, iqr=iqr, warn=save.warn.option))
  }

Answer 2

You could take stargazer::.stargazer.wrap to task and change how the model summary is calculated, but it's about 8k lines long.

As a workaround, it's probably much easier to first compute a smaller regression that can be digested by the table-generating function, then replace the coefficients and statistics before compiling.

I show you the method using the texreg package, it may also be possible with stargazer.

foo <- y ~ x
beta0 <- betareg::betareg(foo, data=df[sample(nrow(df), 1e3), ])  ## small version
beta <- betareg::betareg(foo, data=df)  ## actual version

(beta.sum <- betareg:::summary.betareg(beta, type='pearson')$coe)  ## actual summary
# $mean
#                 Estimate  Std. Error    z value  Pr(>|z|)
# (Intercept)  0.007933124 0.006967576  1.1385773 0.2548795
# x           -0.010845820 0.012080681 -0.8977822 0.3693017
# 
# $precision
#       Estimate  Std. Error z value Pr(>|z|)
# (phi) 1.999218 0.007501926 266.494        0

Table

Screen:

texreg::screenreg(beta0,
                  digits=3,
                  override.coef=unlist(lapply(beta.sum, \(x) x[, 1])),
                  override.se=unlist(lapply(beta.sum, \(x) x[, 2])),
                  override.pvalues=unlist(lapply(beta.sum, \(x) x[, 4]))
)
# ==============================
#                   Model 1     
# ------------------------------
# (Intercept)          0.008    
#                     (0.007)   
# x                   -0.011    
#                     (0.012)   
# Precision: (phi)     1.999 ***
#                     (0.008)   
# ------------------------------
# Pseudo R^2           0.000    
# Log Likelihood       2.137    
# Num. obs.         1000        
# ==============================
# *** p < 0.001; ** p < 0.01; * p < 0.05

LaTeX:

texreg::texreg(beta0,
                  digits=3,
                  override.coef=unlist(lapply(beta.sum, \(x) x[, 1])),
                  override.se=unlist(lapply(beta.sum, \(x) x[, 2])),
                  override.pvalues=unlist(lapply(beta.sum, \(x) x[, 4]))
)
# \begin{table}
# \begin{center}
# \begin{tabular}{l c}
# \hline
# & Model 1 \\
# \hline
# (Intercept)      & $0.008$       \\
# & $(0.007)$     \\
# x                & $-0.011$      \\
# & $(0.012)$     \\
# Precision: (phi) & $1.999^{***}$ \\
# & $(0.008)$     \\
# \hline
# Pseudo R$^2$     & $0.000$       \\
# Log Likelihood   & $2.137$       \\
# Num. obs.        & $1000$        \\
# \hline
# \multicolumn{2}{l}{\scriptsize{$^{***}p<0.001$; $^{**}p<0.01$; $^{*}p<0.05$}}
# \end{tabular}
# \caption{Statistical models}
# \label{table:coefficients}
# \end{center}
# \end{table}

For getting HTML code, use:

texreg::htmlreg

---

Data:

set.seed(1)
n <- 1e5
df <- data.frame(y=runif(n), x=runif(n))