R:观星台中强大的SE和模型诊断 [英] R: Robust SE's and model diagnostics in stargazer table
问题描述
我尝试使用stargazer程序包将通过ivreg从AER程序包中通过ivreg生成的一些2SLS回归输出放到Latex文档中.我有几个问题，但是我似乎无法解决自己.
1.我无法弄清楚ivreg摘要提供的如何插入模型诊断.即弱仪器测试，Wu-Hausmann和Sargan测试.我希望它们具有通常在表格下方报告的统计信息，例如观察数，R平方和残差. SE. stargazer函数似乎没有参数，您可以在其中提供带有其他诊断信息的列表.我没有在我的示例中提到这个问题，因为老实说我不知道从哪里开始.
2.我想将普通标准误差与健壮标准误差交换，而我发现做到这一点的唯一方法是产生具有健壮标准误差的对象，并使用se = list()将它们添加到stargazer函数中.我将其放在下面的最小工作示例中.是否可以使用一种更优雅的方式对此进行编码，或者重新估计该模型并将其保存为可靠的标准错误?帮助表示赞赏.
I try to put some 2SLS regression outputs generated via ivreg from the AER package into a Latex document using the stargazer package. I have a couple of problems however that i can't seem to solve myself.
1. I can't figur out on how to insert model diagnostics as provided by the summary of ivreg. Namely weak instruments tests, Wu-Hausmann and Sargan Test. I would like to have them with the statistics usually reported underneath the table like number of observations, R squared, and Resid. SE. The stargazer function doesn't seem to have an argument where you can provide a list with additional diagnostics. I didn't put this into my example because i honestly have no clue where to begin.
2. I want to exchange the normal standard errors with robust standard errors and the only way to do this that i found is producing objects with robust standard errors and adding them in the stargazer function with se=list(). I put this into the minimum working example below. Is there maybe a more elegant way to code this or maybe reestimate the model and save it with robust standard errors? Help appreciated.
library(AER)
library(stargazer)
y <- rnorm(100, 5, 10)
x <- rnorm(100, 3, 15)
z <- rnorm(100, 3, 7)
a <- rnorm(100, 1, 7)
b <- rnorm(100, 3, 5)
# Fitting IV models
fit1 <- ivreg(y ~ x + a |
a + z,
model = TRUE)
fit2 <- ivreg(y ~ x + a |
a + b + z,
model = TRUE)
# Here are the se's and the diagnostics i want
summary(fit1, vcov = sandwich, diagnostics=T)
summary(fit2, vcov = sandwich, diagnostics=T)
# Getting robust se's, i think HC0 is the standard
# used with "vcov=sandwich" from the above summary
cov1 <- vcovHC(fit1, type = "HC0")
robust1 <- sqrt(diag(cov1))
cov2 <- vcovHC(fit2, type = "HC0")
robust2 <- sqrt(diag(cov1))
# Create latex table
stargazer(fit1, fit2, type = "latex", se=list(robust1, robust2))
推荐答案
这是您想要做的一种方法:
Here's one way to do what you want:
require(lmtest)
rob.fit1 <- coeftest(fit1, function(x) vcovHC(x, type="HC0"))
rob.fit2 <- coeftest(fit2, function(x) vcovHC(x, type="HC0"))
summ.fit1 <- summary(fit1, vcov. = function(x) vcovHC(x, type="HC0"), diagnostics=T)
summ.fit2 <- summary(fit2, vcov. = function(x) vcovHC(x, type="HC0"), diagnostics=T)
stargazer(fit1, fit2, type = "text",
se = list(rob.fit1[,"Std. Error"], rob.fit2[,"Std. Error"]),
add.lines = list(c(rownames(summ.fit1$diagnostics)[1],
round(summ.fit1$diagnostics[1, "p-value"], 2),
round(summ.fit2$diagnostics[1, "p-value"], 2)),
c(rownames(summ.fit1$diagnostics)[2],
round(summ.fit1$diagnostics[2, "p-value"], 2),
round(summ.fit2$diagnostics[2, "p-value"], 2)) ))
哪个会产生:
==========================================================
Dependent variable:
----------------------------
y
(1) (2)
----------------------------------------------------------
x -1.222 -0.912
(1.672) (1.002)
a -0.240 -0.208
(0.301) (0.243)
Constant 9.662 8.450**
(6.912) (4.222)
----------------------------------------------------------
Weak instruments 0.45 0.56
Wu-Hausman 0.11 0.18
Observations 100 100
R2 -4.414 -2.458
Adjusted R2 -4.526 -2.529
Residual Std. Error (df = 97) 22.075 17.641
==========================================================
Note: *p<0.1; **p<0.05; ***p<0.01
如您所见，这允许在各个模型中手动包括诊断.
As you can see, this allows manually including the diagnostics in the respective models.
您可以通过创建一个函数来自动化此方法，该函数接受模型列表(例如list(summ.fit1, summ.fit2))并输出se或add.lines参数所需的对象.
You could automate this approach by creating a function that takes in a list of models (e.g. list(summ.fit1, summ.fit2)) and outputs the objects required by se or add.lines arguments.
gaze.coeft <- function(x, col="Std. Error"){
stopifnot(is.list(x))
out <- lapply(x, function(y){
y[ , col]
})
return(out)
}
gaze.coeft(list(rob.fit1, rob.fit2))
gaze.coeft(list(rob.fit1, rob.fit2), col=2)
将同时获取list个coeftest对象，并产生se期望的SE向量:
Will both take in a list of coeftest objects, and yield the SEs vector as expected by se:
[[1]]
(Intercept) x a
6.9124587 1.6716076 0.3011226
[[2]]
(Intercept) x a
4.2221491 1.0016012 0.2434801
可以对诊断进行相同的操作:
Same can be done for the diagnostics:
gaze.lines.ivreg.diagn <- function(x, col="p-value", row=1:3, digits=2){
stopifnot(is.list(x))
out <- lapply(x, function(y){
stopifnot(class(y)=="summary.ivreg")
y$diagnostics[row, col, drop=FALSE]
})
out <- as.list(data.frame(t(as.data.frame(out)), check.names = FALSE))
for(i in 1:length(out)){
out[[i]] <- c(names(out)[i], round(out[[i]], digits=digits))
}
return(out)
}
gaze.lines.ivreg.diagn(list(summ.fit1, summ.fit2), row=1:2)
gaze.lines.ivreg.diagn(list(summ.fit1, summ.fit2), col=4, row=1:2, digits=2)
两个通话都会产生:
$`Weak instruments`
[1] "Weak instruments" "0.45" "0.56"
$`Wu-Hausman`
[1] "Wu-Hausman" "0.11" "0.18"
现在stargazer()调用变得如此简单，产生与上面相同的输出:
Now the stargazer() call becomes as simple as this, yielding identical output as above:
stargazer(fit1, fit2, type = "text",
se = gaze.coeft(list(rob.fit1, rob.fit2)),
add.lines = gaze.lines.ivreg.diagn(list(summ.fit1, summ.fit2), row=1:2))
这篇关于R:观星台中强大的SE和模型诊断的文章就介绍到这了，希望我们推荐的答案对大家有所帮助，也希望大家多多支持IT屋！
