logo标签列表

Stargazer中不准确的P值

rstargazer
10

我希望在stargazer中回归输出的显著性星号与“正常输出”中的相同。

我生成数据。

library("stargazer"); library("lmtest"); library("sandwich")
set.seed(1234)
df <- data.frame(y=1001:1100)
df$x <- c(1:70,-100:-71) + rnorm(100, 0, 74.8)
model <- lm(log(y) ~ x, data=df)

并获得一些模型估计，其中x的系数具有0.1023的p值。

coeftest(model, vcov = vcovHC(model, type="HC3"))

我希望能用LaTeX展现这些结果。基于同一个函数，我计算异方差一致性标准误并让stargazer使用它们。

stderr_HC3_model <- sqrt(diag(vcovHC(model, type = "HC3")))
stargazer(model, se=list(stderr_HC3_model))

当alpha=10％时，stargazer输出的系数中有一个星号表示显着性。我希望stargazer可以给出与coeftest相同的结果。（由于与Stata的可比性，在该软件中，reg L_y x，vce（hc3）恰好给出了coeftest的结果。）

我尝试了stargazer选项p.auto、t.auto，但并没有帮助。当我执行“stargazer”时，我无法查看其他情况下可能查看的底层代码。怎么办？

Richards的答案对我有帮助。我说明了我用来发出多个回归（比如ols_a和ols_b）的步骤。

ses <- list(coeftest(ols_a, vcov = vcovHC(ols_a, type="HC3"))[,2],
        coeftest(ols_b, vcov = vcovHC(ols_b, type="HC3"))[,2])
pvals <- list(coeftest(ols_a, vcov = vcovHC(ols_a, type="HC3"))[,4],
          coeftest(ols_b, vcov = vcovHC(ols_b, type="HC3"))[,4])
stargazer(ols_a, ols_b, type="text", p=pvals, se=ses)
2
要获取源代码，请尝试使用 stargazer:::.stargazer.wrap。 - shadow
3个回答
5
您需要提供与coeftest相关联的p值。来自man page。

p是一组数字向量，将替换每个模型的默认p值。按元素名称匹配。这些将形成关于显着性星号决策的基础。

以下应该有效。
test <- coeftest(model, vcov = vcovHC(model, type="HC3"))
ses <- test[, 2]
pvals <- test[, 4]
stargazer(model, type="text", p=pvals, se=ses)

这提供了以下内容。
===============================================
                        Dependent variable:    
                    ---------------------------
                              log(y)           
-----------------------------------------------
x                            -0.00005          


Constant                     6.956***          
                              (0.003)          

-----------------------------------------------
Observations                    100            
R2                             0.026           
Adjusted R2                    0.016           
Residual Std. Error       0.027 (df = 98)      
F Statistic             2.620 (df = 1; 98)     
===============================================
Note:               *p<0.1; **p<0.05; ***p<0.01
3

这可能只是一个小问题，但Richard的回答并不完全正确 - 他的星座仪输出没有报告变量x的标准误差和潜在显着性星号。

此外，当仅报告单个模型时，在stargazer手册中必须提供系数、标准误差、p值和t值的列表。否则，stargazer将报告空列表。

稍作改正后的示例：

test <- coeftest(model, vcov = vcovHC(model, type="HC3"))
ses <- list(test[, 2])
pvals <- list(test[, 4])
stargazer(model, type="text", p=pvals, se=ses)

输出：

=======================================================================
                                         Dependent variable:           
                              -----------------------------------------
                                        Daily added investors          
                                              negative                 
                                              binomial                 
-----------------------------------------------------------------------
log(lag_raised_amount + 1)                    -0.466***                
                                               (0.124)                 

lag_target1                                   -0.661***                
                                               (0.134)                 

Constant                                      -3.480**                 
                                               (1.290)                 

-----------------------------------------------------------------------
Observations                                    6,513                  
Log Likelihood                                 -8,834                
theta                                     1.840*** (0.081)             
Akaike Inf. Crit.                              17,924                
=======================================================================
Note:                         + p<0.1; * p<0.05; ** p<0.01; *** p<0.001
2

使用se参数存在固有的风险。

在使用此方法时，用户应谨慎考虑参数t.auto和p.auto，它们都默认为TRUE。我认为将它们都设置为FALSE并手动提供t和p值是明智的选择。

如果不这样做，您可能会出现与显示的p值不同步的显著性星号。(我怀疑stargazer将简单地重用这些现在与默认值不同的se，并使用此输入重新计算显示的星号；这自然会产生意想不到的结果。)

另请参见：

网页内容由stack overflow 提供, 点击上面的
可以查看英文原文，
原文链接