通过最大可能性将系数估计到一个观星表中
Stargazer为lm（和其他）物品生产非常好的乳胶桌子。 假设我已经最大可能地拟合了一个模型。 我想要观星者为我的估计制作一个类似于lm的桌子。 我怎样才能做到这一点？
虽然这有点不好意思，但是有一种方法可能是创build一个包含我的估计的“虚假”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? 更精确地说，对于lm对象，stargazer很好地在表格的顶部打印因variables，在相应的估计值下面包括括号内的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
我不知道你是如何承诺使用stargazer，但你可以尝试使用扫帚和xtable包，问题是它不会给你的优化模型的标准错误
library(broom) library(xtable) xtable(tidy(model1)) xtable(tidy(model2))
我只是有这个问题，并通过使用coef克服了这一点，并omit星空中的function…例如
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="")
你需要首先实例化一个虚拟lm对象，然后打扮：
#... 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
（然后当然要确保其余的汇总统计是正确的）
