Различные надежные стандартные ошибки регрессии логита в Stata и R

Я пытаюсь реплицировать регрессию logit из Stata в R. В Stata я использую опцию "robust", чтобы иметь надежную стандартную ошибку (стандартную ошибку, совместимую с гетероседастичностью). Я могу копировать точно такие же коэффициенты из Stata, но я не могу иметь такую ​​же надежную стандартную ошибку с пакетом "сэндвич".

Я пробовал некоторые примеры линейной регрессии OLS; похоже, что сэндвич-оценки R и Stata дают мне такую ​​же надежную стандартную ошибку для OLS. Кто-нибудь знает, как Stata вычисляет сэндвич-оценку для нелинейной регрессии, в моем случае логит-регрессия?

Спасибо!

Прикрепленные коды: в R:

library(sandwich)
library(lmtest) 
mydata <- read.csv("http://www.ats.ucla.edu/stat/data/binary.csv") 
mydata$rank<-factor(mydata$rank) 
myfit<-glm(admit~gre+gpa+rank,data=mydata,family=binomial(link="logit")) 
summary(myfit) 
coeftest(myfit, vcov = sandwich) 
coeftest(myfit, vcov = vcovHC(myfit, "HC0")) 
coeftest(myfit, vcov = vcovHC(myfit)) 
coeftest(myfit, vcov = vcovHC(myfit, "HC3")) 
coeftest(myfit, vcov = vcovHC(myfit, "HC1")) 
coeftest(myfit, vcov = vcovHC(myfit, "HC2")) 
coeftest(myfit, vcov = vcovHC(myfit, "HC")) 
coeftest(myfit, vcov = vcovHC(myfit, "const")) 
coeftest(myfit, vcov = vcovHC(myfit, "HC4")) 
coeftest(myfit, vcov = vcovHC(myfit, "HC4m")) 
coeftest(myfit, vcov = vcovHC(myfit, "HC5"))

Stata:

use http://www.ats.ucla.edu/stat/stata/dae/binary.dta, clear 
logit admit gre gpa i.rank, robust
1 ответ

По умолчанию так называемые "надежные" стандартные ошибки в Stata соответствуют тому, что вычисляет sandwich() из пакета с одним и тем же именем. Единственное различие заключается в том, как выполняется корректировка с конечной выборкой. В функции sandwich(...) по умолчанию не выполняется настройка конечных выборок, т.е. Сэндвич делится на 1/n, где n - количество наблюдений. В качестве альтернативы можно использовать sandwich(..., adjust = TRUE), которая делит на 1/(n - k), где k - число регрессоров. И Stata делит на 1/(n - 1).

Конечно, асимптотически они вообще не различаются. И за исключением нескольких особых случаев (например, линейной регрессии OLS) нет аргументов для 1/(n - k) или 1/(n - 1), чтобы работать "правильно" в конечных выборках (например, несмещенность). По крайней мере, насколько мне известно.

Итак, чтобы получить те же результаты, что и в Stata, вы можете сделать:

sandwich1 <- function(object, ...) sandwich(object) * nobs(object) / (nobs(object) - 1)
coeftest(myfit, vcov = sandwich1)

Это дает

z test of coefficients:
 Estimate Std. Error z value Pr(>|z|) 
(Intercept) -3.9899791 1.1380890 -3.5059 0.0004551 ***
gre 0.0022644 0.0011027 2.0536 0.0400192 * 
gpa 0.8040375 0.3451359 2.3296 0.0198259 * 
rank2 -0.6754429 0.3144686 -2.1479 0.0317228 * 
rank3 -1.3402039 0.3445257 -3.8900 0.0001002 ***
rank4 -1.5514637 0.4160544 -3.7290 0.0001922 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

И только для записи: в случае двоичного ответа эти "надежные" стандартные ошибки не устойчивы ни к чему. При условии, что модель правильно указана, они согласуются, и их можно использовать, но они не защищают от какой-либо ошибки в модели. Поскольку основное допущение для стандартных ошибок сэндвича заключается в том, что модельное уравнение (или, точнее, соответствующая функция оценки) правильно задано, в то время как остальная часть модели может быть неопределенной. Однако в двоичной регрессии нет места для спецификации неопределенности, потому что модельное уравнение состоит только из среднего (= вероятность), а вероятность - среднее и среднее значение соответственно. Это контрастирует с линейной или подсчетной регрессией данных, где могут быть гетероскедастичность, сверхдисперсия и т.д.

licensed under cc by-sa 3.0 with attribution.