获得通过最大似然估计的系数到星形表中

Adr*_*ian 83 optimization r lm stargazer

Stargazer为lm(和其他)物体生产非常漂亮的乳胶桌.假设我按最大可能性拟合模型.我希望观星者为我的估计生成一个类似于lm的表格.我怎样才能做到这一点?

虽然它有点hacky,但有一种方法可能是创建一个包含我的估计值的"假"lm对象 - 我认为只要summary(my.fake.lm.object)有效,它就会起作用.那容易行吗?

一个例子:

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?
Run Code Online (Sandbox Code Playgroud)

更确切地说:对于lm个对象,stargazer很好地打印表格顶部的因变量,包括相应估计值下面的括号中的SE,并且在表格的底部有R ^ 2和观察数量.是否有一种(简单的)方法来获得与最大可能性估计的"自定义"模型相同的行为,如上所述?

以下是我将我的优化输出装扮为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
Run Code Online (Sandbox Code Playgroud)

小智 0

我不知道你对使用 stargazer 的决心有多大,但你可以尝试使用 broom 和 xtable 包,问题是它不会给你优化模型的标准错误

library(broom)
library(xtable)
xtable(tidy(model1))
xtable(tidy(model2))
Run Code Online (Sandbox Code Playgroud)