r - 如何通过 vglm tobit 模型使用健壮的 SE 和集群 SE？

Question

我正在尝试将 tobit 模型从 Stata 迁移到 R。

稳健的 Stata 命令只需将 ,vce(robust) 添加到模型中。对于集群，它将是 ,vce(cluster idvar)。

可重现的 Stata 示例:

use http://www.ats.ucla.edu/stat/stata/dae/tobit, clear
tobit apt read math i.prog, ul(800)
tobit apt read math i.prog, ul(800) vce(cluster prog)

可重现的 R 示例:

library("VGAM")

dat <- read.csv("http://www.ats.ucla.edu/stat/data/tobit.csv")

summary(m <- vglm(apt ~ read + math + prog, tobit(Upper = 800), data = dat))

我的理解是，coeftest(m, vcov =三明治)应该给我强大的se。

但我得到以下信息:错误:未为此 S4 类定义 $ 运算符。

有人可以建议一种方法来估计 vglm 模型的鲁棒 se 并使用 vglm 聚类 se 吗？

Answer 1

我自己花了一整天的时间研究这个问题，我想我终于找到了一个合适的包:Zelig。

http://docs.zeligproject.org/en/latest/zelig-tobit.html

比较无聚类与聚类:

没有

> summary(m <- zelig(apt ~ read + math + prog,
          below=0, above=Inf, model="tobit", data = dat))


 How to cite this model in Zelig:
  Kosuke Imai, Gary King, and Olivia Lau. 2015.
  "tobit: Linear regression for Left-Censored Dependent Variable"
  in Kosuke Imai, Gary King, and Olivia Lau, "Zelig: Everyone's Statistical Software,"
  http://gking.harvard.edu/zelig


Call:
"survreg"(formula = formula, dist = "gaussian", data = data, 
    robust = robust)
                Value Std. Error     z        p
(Intercept)    242.74     29.760  8.16 3.45e-16
read             2.55      0.576  4.43 9.24e-06
math             5.38      0.651  8.27 1.31e-16
proggeneral    -13.74     11.596 -1.18 2.36e-01
progvocational -48.83     12.818 -3.81 1.39e-04
Log(scale)       4.12      0.050 82.41 0.00e+00

Scale= 61.6 

Gaussian distribution
Loglik(model)= -1107.9   Loglik(intercept only)= -1202.8
    Chisq= 189.72 on 4 degrees of freedom, p= 0 
Number of Newton-Raphson Iterations: 5 
n= 200 

与

> summary(m <- zelig(apt ~ read + math + prog, below=0,
          above=Inf, model="tobit",
          data = dat,robust=T,cluster="prog"))


 How to cite this model in Zelig:
  Kosuke Imai, Gary King, and Olivia Lau. 2015.
  "tobit: Linear regression for Left-Censored Dependent Variable"
  in Kosuke Imai, Gary King, and Olivia Lau, "Zelig: Everyone's Statistical Software,"
  http://gking.harvard.edu/zelig


Call:
"survreg"(formula = formula, dist = "gaussian", data = data, 
    robust = robust)
                Value Std. Err (Naive SE)       z        p
(Intercept)    242.74   2.8315     29.760   85.73 0.00e+00
read             2.55   0.3159      0.576    8.08 6.40e-16
math             5.38   0.2770      0.651   19.44 3.78e-84
proggeneral    -13.74   0.3252     11.596  -42.25 0.00e+00
progvocational -48.83   0.1978     12.818 -246.83 0.00e+00
Log(scale)       4.12   0.0586      0.050   70.34 0.00e+00

Scale= 61.6 

Gaussian distribution
Loglik(model)= -1107.9   Loglik(intercept only)= -1202.8
    Chisq= 189.72 on 4 degrees of freedom, p= 0 
(Loglikelihood assumes independent observations)
Number of Newton-Raphson Iterations: 5 
n= 200