varias - Obtener coeficientes estimados por máxima verosimilitud en una tabla de estrellas

superponer graficas en r

Question

Stargazer produce muy buenas tablas de látex para objetos de lm (y otros). Supongamos que he ajustado un modelo por máxima verosimilitud. Me gustaría que Stargazer produzca una tabla tipo lm para mis estimaciones. ¿Cómo puedo hacer esto?

Aunque es un poco hacky, una forma podría ser crear un objeto lm "falso" que contenga mis estimaciones: creo que funcionaría siempre que funcione el resumen (my.fake.lm.object). ¿Es eso factible?

Un ejemplo:

library(stargazer)

N <- 200
df <- data.frame(x=runif(N, 0, 50))
df$y <- 10 + 2 * df$x + 4 * rt(N, 4)  # True params
plot(df$x, df$y)

model1 <- lm(y ~ x, data=df)
stargazer(model1, title="A Model")  # I''d like to produce a similar table for the model below

ll <- function(params) {
    ## Log likelihood for y ~ x + student''s t errors
    params <- as.list(params)
    return(sum(dt((df$y - params$const - params$beta*df$x) / params$scale, df=params$degrees.freedom, log=TRUE) -
               log(params$scale)))
}

model2 <- optim(par=c(const=5, beta=1, scale=3, degrees.freedom=5), lower=c(-Inf, -Inf, 0.1, 0.1),
                fn=ll, method="L-BFGS-B", control=list(fnscale=-1), hessian=TRUE)
model2.coefs <- data.frame(coefficient=names(model2$par), value=as.numeric(model2$par),
                           se=as.numeric(sqrt(diag(solve(-model2$hessian)))))

stargazer(model2.coefs, title="Another Model", summary=FALSE)  # Works, but how can I mimic what stargazer does with lm objects?

Para ser más precisos: con objetos lm, Stargazer imprime muy bien la variable dependiente en la parte superior de la tabla, incluye SE entre paréntesis debajo de las estimaciones correspondientes, y tiene el R ^ 2 y el número de observaciones en la parte inferior de la tabla. ¿Hay una (n fácil) forma de obtener el mismo comportamiento con un modelo "personalizado" estimado por máxima verosimilitud, como se indicó anteriormente?

Estos son mis débiles intentos de disfrazar mi salida optimizada como un objeto lm:

model2.lm <- list()  # Mimic an lm object
class(model2.lm) <- c(class(model2.lm), "lm")
model2.lm$rank <- model1$rank  # Problematic?
model2.lm$coefficients <- model2$par
names(model2.lm$coefficients)[1:2] <- names(model1$coefficients)
model2.lm$fitted.values <- model2$par["const"] + model2$par["beta"]*df$x
model2.lm$residuals <- df$y - model2.lm$fitted.values
model2.lm$model <- df
model2.lm$terms <- model1$terms  # Problematic?
summary(model2.lm)  # Not working

Answer 1

No sé qué tan comprometido está usted con el uso de Stargazer, pero puede intentar usar los paquetes de escoba y xtable, el problema es que no le dará los errores estándar para el modelo optim.

library(broom)
library(xtable)
xtable(tidy(model1))
xtable(tidy(model2))

Answer 2

Primero necesita crear una instancia de un objeto ficticio lm y luego vestirlo:

#...
model2.lm = lm(y ~ ., data.frame(y=runif(5), beta=runif(5), scale=runif(5), degrees.freedom=runif(5)))
model2.lm$coefficients <- model2$par
model2.lm$fitted.values <- model2$par["const"] + model2$par["beta"]*df$x
model2.lm$residuals <- df$y - model2.lm$fitted.values
stargazer(model2.lm, se = list(model2.coefs$se), summary=FALSE, type=''text'')

# ===============================================
#                         Dependent variable:    
#                     ---------------------------
#                                  y             
# -----------------------------------------------
# const                        10.127***         
#                               (0.680)          
#                                                
# beta                         1.995***          
#                               (0.024)          
#                                                
# scale                        3.836***          
#                               (0.393)          
#                                                
# degrees.freedom              3.682***          
#                               (1.187)          
#                                                
# -----------------------------------------------
# Observations                    200            
# R2                             0.965           
# Adjusted R2                    0.858           
# Residual Std. Error       75.581 (df = 1)      
# F Statistic              9.076 (df = 3; 1)     
# ===============================================
# Note:               *p<0.1; **p<0.05; ***p<0.01

(y luego, por supuesto, asegúrese de que las estadísticas resumidas restantes sean correctas)

Answer 3

Solo estaba teniendo este problema y lo supero mediante el uso del coef , y omit funciones dentro de stargazer ... por ejemplo

stargazer(regressions, ...
                     coef = list(... list of coefs...),
                     se = list(... list of standard errors...),
                     omit = c(sequence),
                     covariate.labels = c("new names"),
                     dep.var.labels.include = FALSE,
                     notes.append=FALSE), file="")