R equivalente a la opción de clúster cuando se usa regresión binomial negativa

10

Estoy tratando de replicar el trabajo de un colega y estoy trasladando el análisis de Stata a R. Los modelos que emplea invocan la opción "cluster" dentro de la función nbreg para agrupar los errores estándar.

Consulte http://repec.org/usug2007/crse.pdf para obtener una descripción bastante completa del qué y por qué de esta opción

Mi pregunta es cómo invocar esta misma opción para la regresión binomial negativa dentro de R?

El modelo principal en nuestro artículo se especifica en Stata de la siguiente manera

 xi: nbreg cntpd09 logpop08 pcbnkthft07 pccrunion07 urbanpop pov00 pov002 edu4yr ///
 black04 hispanic04 respop i.pdpolicy i.maxloan rollover i.region if isser4 != 1,   
 cluster(state)

y he reemplazado esto con

pday<-glm.nb(cntpd09~logpop08+pcbnkthft07+pccrunion07+urbanpop+pov00+pov002+edu4yr+
black04+hispanic04+respop+as.factor(pdpolicy)+as.factor(maxloan)+rollover+
as.factor(region),data=data[which(data$isser4 != 1),])

que obviamente carece de la pieza de errores agrupados.

¿Es posible hacer una replicación exacta? ¿Si es así, cómo? Si no, ¿cuáles son algunas alternativas razonables?

Gracias

[Editar] Como se señaló en los comentarios, esperaba una solución que no me llevara al ámbito de los modelos multinivel. Si bien mi entrenamiento me permite ver que estas cosas deberían estar relacionadas, es más un salto de lo que me siento cómodo. Como tal, seguí cavando y encontré este enlace: http://landroni.wordpress.com/2012/06/02/fama-macbeth-and-cluster-robust-by-firm-and-time-standard-errors-in- r /

eso apunta a un código bastante sencillo para hacer lo que quiero:

library(lmtest)
pday<-glm.nb(cntpd09~logpop08+pcbnkthft07+pccrunion07+urbanpop+pov00+pov002+edu4yr+
 black04+hispanic04+respop+as.factor(pdpolicy)+as.factor(maxloan)+rollover+
 as.factor(region),data=data[which(data$isser4 != 1),])
summary(pday)

coeftest(pday, vcov=function(x) vcovHC(x, cluster="state", type="HC1"))

Sin embargo, esto no replica los resultados del análisis en Stata, probablemente porque está diseñado para funcionar en OLS no binomial negativo. Entonces la búsqueda continúa. Cualquier sugerencia sobre dónde me estoy equivocando sería muy apreciada

csfowler
fuente
3
Puede encontrar útiles las notas de Ben Bolker aquí.
fmark
Y vea también esta pregunta anterior
fmark
FYI aquí es una definición de los robustos errores estándar agrupados de Stata. No parecen tan difíciles de implementar. En mi opinión, puede estar mejor con errores estándar bootstrapped o jackknifed de todos modos (consulte la ayuda en vce ). Sin embargo, no puedo sugerir ningún paquete R. ¡Buena suerte en encontrar un reemplazo!
Andy W
Gracias @fmark: comentarios muy útiles, mucho mejores que mi "respuesta" y los he actualizado en consecuencia.
Peter Ellis
Gracias a todos. Creo que la respuesta corta a mi pregunta es que no hay un reemplazo directo (por ejemplo, una función prefabricada que reemplaza exactamente la opción de clúster). Claramente, alguien con más experiencia puede ver el camino a través de las notas de Ben Bolker, pero me lleva a un nuevo territorio donde no podía estar seguro de que estaba recibiendo las declaraciones de fórmula correctas. No estoy seguro de cuál es la forma apropiada de decir "Gracias" sin aceptar una respuesta, pero sí tengo mi agradecimiento, y las deficiencias son mías.
csfowler

Respuestas:

4

Este documento muestra cómo obtener SE de clúster para una regresión glm:

http://dynaman.net/R/clrob.pdf

EddieMcGoldrick
fuente
Tendré que hacer una comparación de prueba con los resultados de stata, pero esto se ve exactamente como lo que esperaba.
csfowler
1

Esta no es una respuesta totalmente satisfactoria ...

No lo he probado, pero parece que el paquete glmmADMB podría hacer lo que quieres.

Pellizcaré descaradamente el comentario de @fmark sobre la pregunta y estaré de acuerdo con él en que las notas de Ben Bolker son útiles, como es esta pregunta anterior , que no es un duplicado exacto pero cubre temas muy similares.

Peter Ellis
fuente