2

对于我的stargazer表，我想包括弱仪器的诊断和 Wu Hausman 测试。我阅读了以下答案：R: Robust SE's and model diagnostics in stargazer table。

虽然，我真的不明白这一点。如何将其应用于具有 9 或 10 列的更大模型？如何插入您在编码时获得的诊断信息summary(ivreg1, diagnostics=TRUE)？我对观星者表有以下代码：

ivreg1 <- ivreg(budget ~ policyfactor1 + control1 + control2 + control3 | iv + control1 + control2 + control3, data=df)
ivreg2 <- ivreg(budget ~ policyfactor2 + control1 + control2 + control3 | iv + control1 + control2 + control3, data=df)
ivreg3 <- ivreg(budget ~ policyfactor3 + control1 + control2 + control3 | iv + control1 + control2 + control3, data=df)
ivreg4 <- ivreg(budget ~ policyfactor4 + control1 + control2 + control3 | iv + control1 + control2 + control3, data=df)
ivreg5 <- ivreg(budget ~ policyfactor5 + control1 + control2 + control3 | iv + control1 + control2 + control3, data=df)
ivreg6 <- ivreg(budget ~ policyfactor6 + control1 + control2 + control3 | iv + control1 + control2 + control3, data=df)
ivreg7 <- ivreg(budget ~ policyfactor7 + control1 + control2 + control3 | iv + control1 + control2 + control3, data=df)
ivreg8 <- ivreg(budget ~ policyfactor8 + control1 + control2 + control3 | iv + control1 + control2 + control3, data=df)
ivreg9 <- ivreg(budget ~ policyfactor9 + control1 + control2 + control3 | iv + control1 + control2 + control3, data=df)
stargazer(ivreg1, ivreg2, ivreg3, ivreg4, ivreg5, ivreg6, ivreg7, ivreg8, ivreg9, title="Regression results", align=TRUE, column.sep.width = "-15pt", font.size = "small", type="latex")

我有以下数据框示例：

data.frame(
         budget = c(10, 8, -7, 8, 3, 2, 0.5, 1.5, -2, -15, 30, -0.5, 12, 18),
  policyfactor1 = c(1L, 1L, 0L, 0L, 0L, 1L, 1L, 0L, 0L, 0L, 1L, 0L, 1L, 1L),
  policyfactor2 = c(0L, 1L, 0L, 1L, 0L, 0L, 1L, 1L, 0L, 0L, 1L, 0L, 1L, 0L),
  policyfactor3 = c(0L, 1L, 1L, 0L, 1L, 1L, 1L, 0L, 0L, 1L, 1L, 1L, 0L, 0L),
  policyfactor4 = c(1L, 0L, 1L, 0L, 0L, 1L, NA, NA, 1L, 1L, 0L, 1L, 1L, 0L),
  policyfactor5 = c(1L, 0L, 1L, 0L, 0L, 1L, 1L, 0L, 0L, 1L, 1L, 0L, 1L, 1L),
  policyfactor6 = c(0L, 1L, 1L, 0L, 1L, 1L, 0L, 1L, 1L, 1L, 0L, 1L, NA, 0L),
  policyfactor7 = c(0L, 1L, 1L, 1L, 1L, 0L, 1L, 0L, 1L, 0L, 0L, 1L, 0L, 0L),
  policyfactor8 = c(1L, 0L, 1L, 1L, NA, 1L, 1L, 0L, 0L, 1L, 0L, 0L, 1L, 0L),
  policyfactor9 = c(0L, 0L, 1L, 0L, 1L, 0L, 1L, 1L, 0L, 1L, NA, 1L, 1L, 1L),
             IV = c(17L,46L,18L,23L,35L,10L,0L,19L,
                    15L,12L,21L,6L,27L,10L),
       control1 = c(29.4,27.4,33.2,27.1,27.4,34.2,
                    26.8,32.9,26,26.3,27.3,33.4,23.5,31.3),
       control2 = c(1.3,20,-0.2,3,3.4,0.3,-1.1,1.9,
                    -2,-1.6,2.8,1.9,2,1.8)
)
#>    budget policyfactor1 policyfactor2 policyfactor3 policyfactor4 policyfactor5
#> 1    10.0             1             0             0             1             1
#> 2     8.0             1             1             1             0             0
#> 3    -7.0             0             0             1             1             1
#> 4     8.0             0             1             0             0             0
#> 5     3.0             0             0             1             0             0
#> 6     2.0             1             0             1             1             1
#> 7     0.5             1             1             1            NA             1
#> 8     1.5             0             1             0            NA             0
#> 9    -2.0             0             0             0             1             0
#> 10  -15.0             0             0             1             1             1
#> 11   30.0             1             1             1             0             1
#> 12   -0.5             0             0             1             1             0
#> 13   12.0             1             1             0             1             1
#> 14   18.0             1             0             0             0             1
#>    policyfactor6 policyfactor7 policyfactor8 policyfactor9 IV control1 control2
#> 1              0             0             1             0 17     29.4      1.3
#> 2              1             1             0             0 46     27.4     20.0
#> 3              1             1             1             1 18     33.2     -0.2
#> 4              0             1             1             0 23     27.1      3.0
#> 5              1             1            NA             1 35     27.4      3.4
#> 6              1             0             1             0 10     34.2      0.3
#> 7              0             1             1             1  0     26.8     -1.1
#> 8              1             0             0             1 19     32.9      1.9
#> 9              1             1             0             0 15     26.0     -2.0
#> 10             1             0             1             1 12     26.3     -1.6
#> 11             0             0             0            NA 21     27.3      2.8
#> 12             1             1             0             1  6     33.4      1.9
#> 13            NA             0             1             1 27     23.5      2.0
#> 14             0             0             0             1 10     31.3      1.8

由reprex 包（v0.3.0）于 2020-06-20 创建

4

1 回答 1

3

如果您的全局环境中已经有对象，您只需要list模型适合您可能使用的位置mget（不过，请参阅答案底部的替代方案）。

iv.fit <- mget(paste0("ivreg", 1:9))

然后在stargazer::stargazer函数中gaze.lines.ivreg.diagn添加包含诊断的摘要，使用lapply.

library(stargazer)
stargazer(iv.fit[1:3], title="Regression results", align=TRUE, 
          type="text",  ## to get LaTeX output use `type="latex"`
          column.sep.width="-15pt", font.size="small",
          add.lines=
            gaze.lines.ivreg.diagn(lapply(iv.fit[1:3], summary, diagnostics=TRUE), row=1:2))

结果stargazer

# Regression results
# ===========================================================
#                                    Dependent variable:     
#                               -----------------------------
#                                          budget            
#                                  (1)       (2)       (3)   
# -----------------------------------------------------------
# policyfactor1                   1.217                      
#                               (15.134)                     
# 
# policyfactor2                             3.209            
#                                         (39.816)           
# 
# policyfactor3                                       3.302  
#                                                   (44.912) 
# 
# control1                       -0.394    -0.237    -0.513  
#                                (0.984)   (2.407)   (1.687) 
# 
# control2                        0.517     0.439     0.486  
#                                (0.721)   (1.496)   (1.079) 
# 
# Constant                       14.493     9.342    16.727  
#                               (31.255)  (82.795)  (33.982) 
# 
# -----------------------------------------------------------
# Weak instruments                0.17      0.57      0.63   
# Wu-Hausman                      0.36      0.91      0.83   
# Observations                     14        14        14    
# R2                              0.154     0.159    -0.011  
# Adjusted R2                    -0.099    -0.094    -0.315  
# Residual Std. Error (df = 10)  11.465    11.437    12.539  
# ===========================================================
# Note:                           *p<0.1; **p<0.05; ***p<0.01

 

这也是可能的texreg::texreg。同样，我们将摘要放入列表中。

iv.gof2 <- lapply(iv.fit, function(x) 
  summary(x, diagnostics=T)$diagnostics[c("Weak instruments", "Wu-Hausman"), "p-value"])

然后我们使用texreg::extract提取相关模型细节来编辑GOF信息。这给出了iv.fit.ex包含"texreg"对象的列表。

library(texreg)
iv.fit.ex <- lapply(iv.fit, extract)

通过一个小函数，我们可以使用Map.

addGof<- function(x, y) {
  x@gof.names <- c(x@gof.names, names(y))
  x@gof <- c(x@gof, y)
  x@gof.decimal <- c(x@gof.decimal, rep(TRUE, length(y)))
  x
}

iv.fit.ex <- Map(addGof, iv.fit.ex, iv.gof2)

然后我们继续texreg::texreg运行iv.fit.ex。（我在这里显示控制台版本texreg::screenreg）。

screenreg(iv.fit.ex[1:3], digits=3, column.spacing=1,
          custom.model.names=sprintf("   (%s)", 1:length(iv.fit.ex[1:3])))

结果texreg

# ===========================================
#                     (1)      (2)      (3)  
# -------------------------------------------
# (Intercept)       14.493    9.342   16.727 
#                  (31.255) (82.795) (33.982)
# policyfactor1      1.217                   
#                  (15.134)                  
# control1          -0.394   -0.237   -0.513 
#                   (0.984)  (2.407)  (1.687)
# control2           0.517    0.439    0.486 
#                   (0.721)  (1.496)  (1.079)
# policyfactor2               3.209          
#                           (39.816)         
# policyfactor3                        3.302 
#                                    (44.912)
# -------------------------------------------
# R^2                0.154    0.159   -0.011 
# Adj. R^2          -0.099   -0.094   -0.315 
# Num. obs.         14       14       14     
# Weak instruments   0.168    0.567    0.628 
# Wu-Hausman         0.357    0.911    0.825 
# ===========================================
# *** p < 0.001; ** p < 0.01; * p < 0.05

 

顺便说一句，您还可以以一种不那么乏味的方式列出您的回归模型，而不必单独编码每个回归：

fo <- sapply(paste0("policyfactor", 1:9), function(i) 
  as.formula(sprintf("budget ~ %s + control1 + control2 | iv + control1 + control2", i)))
iv.fit <- lapply(fo, function(fo) do.call(ivreg, list(formula=fo, data=quote(dat))))

数据：

dat <- structure(list(budget = c(10, 8, -7, 8, 3, 2, 0.5, 1.5, -2, -15, 
30, -0.5, 12, 18), policyfactor1 = c(1L, 1L, 0L, 0L, 0L, 1L, 
1L, 0L, 0L, 0L, 1L, 0L, 1L, 1L), policyfactor2 = c(0L, 1L, 0L, 
1L, 0L, 0L, 1L, 1L, 0L, 0L, 1L, 0L, 1L, 0L), policyfactor3 = c(0L, 
1L, 1L, 0L, 1L, 1L, 1L, 0L, 0L, 1L, 1L, 1L, 0L, 0L), policyfactor4 = c(1L, 
0L, 1L, 0L, 0L, 1L, NA, NA, 1L, 1L, 0L, 1L, 1L, 0L), policyfactor5 = c(1L, 
0L, 1L, 0L, 0L, 1L, 1L, 0L, 0L, 1L, 1L, 0L, 1L, 1L), policyfactor6 = c(0L, 
1L, 1L, 0L, 1L, 1L, 0L, 1L, 1L, 1L, 0L, 1L, NA, 0L), policyfactor7 = c(0L, 
1L, 1L, 1L, 1L, 0L, 1L, 0L, 1L, 0L, 0L, 1L, 0L, 0L), policyfactor8 = c(1L, 
0L, 1L, 1L, NA, 1L, 1L, 0L, 0L, 1L, 0L, 0L, 1L, 0L), policyfactor9 = c(0L, 
0L, 1L, 0L, 1L, 0L, 1L, 1L, 0L, 1L, NA, 1L, 1L, 1L), iv = c(17L, 
46L, 18L, 23L, 35L, 10L, 0L, 19L, 15L, 12L, 21L, 6L, 27L, 10L
), control1 = c(29.4, 27.4, 33.2, 27.1, 27.4, 34.2, 26.8, 32.9, 
26, 26.3, 27.3, 33.4, 23.5, 31.3), control2 = c(1.3, 20, -0.2, 
3, 3.4, 0.3, -1.1, 1.9, -2, -1.6, 2.8, 1.9, 2, 1.8)), class = "data.frame", row.names = c(NA, 
-14L))
于 2020-06-21T08:24:22.130 回答