¿Cómo configurar y estimar un modelo logit multinomial en R?

20

Ejecuté un modelo logit multinomial en JMP y obtuve resultados que incluían el AIC y los valores p de chi-cuadrado para cada estimación de parámetro. El modelo tiene un resultado categórico y 7 variables explicativas categóricas.

Luego ajusté lo que pensé que construiría el mismo modelo en R, usando la multinomfunción en el paquete nnet .

El código era básicamente:

fit1 <- multinom(y ~ x1+x2+...xn,data=mydata);
summary(fit1);

Sin embargo, los dos dan resultados diferentes. Con JMP, el AIC es 2923.21, y con nnet::multinomel AIC es 3116.588.

Entonces mi primera pregunta es: ¿Está equivocado uno de los modelos?

La segunda cosa es que JMP da valores p de chi-cuadrado para cada estimación de parámetro, que necesito. Ejecutar resumen en el multinom fit1no lo hace, solo da las estimaciones, AIC y Deviance.

Mi segunda pregunta es así: ¿Hay alguna forma de obtener los valores p para el modelo y las estimaciones cuando se usa nnet::multinom?

Sé que mlogit es otro paquete R para esto y parece que su salida incluye los valores p; sin embargo, no he podido ejecutar mlogitusando mis datos. Creo que tenía los datos formateados correctamente, pero decía que tenía una fórmula no válida. Usé la misma fórmula que usé multinom, pero parece que requiere un formato diferente usando una tubería y no entiendo cómo funciona.

Gracias.

Paul
fuente
2
Puede establecer el argumento Hess = TRUE para recuperar el Hessian de multinom y luego calcular los valores p manualmente. Pero le sugiero que use la biblioteca mlogit (nnet puede tener problemas de convergencia cuando las covariables no se escalan correctamente). Las viñetas para mlogit son bastante buenas y deberían ayudarlo a configurar sus datos correctamente. Las viñetas se pueden encontrar en el lugar habitual: cran.r-project.org/web/packages/mlogit
Jason Morgan

Respuestas:

9

Estoy seguro de que ya ha encontrado sus soluciones, ya que esta publicación es muy antigua, pero para aquellos de nosotros que todavía estamos buscando soluciones, he encontrado que http://youtu.be/-Cp_KP9mq94 es una excelente fuente de instrucciones sobre cómo ejecute un modelo de regresión logística multinomial en R usando el paquete mlogit. Si va al sitio web de la academia de econonometría, ella tiene todos los scripts, datos para R y SAS y STATA, creo que o SPSS es uno de esos.

Qué tipo de explica cómo / por qué y qué hacer para transformar sus datos en el formato "largo" frente a "ancho". Lo más probable es que tenga un formato amplio, que requiere transformación.

https://sites.google.com/site/econometricsacademy/econometrics-models/multinomial-probit-and-logit-models

Kerry
fuente
3

En general, las diferencias en los valores de AIC entre dos piezas diferentes de software no son del todo sorprendentes. Calcular las probabilidades a menudo implica una constante que es la misma entre diferentes modelos de los mismos datos. Los diferentes desarrolladores pueden tomar diferentes decisiones sobre qué dejar dentro o fuera de esa constante. Cuando debe preocuparse es cuando las diferencias en los valores de AIC entre dos modelos difieren. En realidad, acabo de notar un argumento que le multinom()permite cambiar cómo se colapsan las filas con valores X idénticos, y que esto afecta la línea de base de la desviación y, por lo tanto, el AIC. Podría probar diferentes valores del argumento summ y ver si eso hace que las desviaciones estén de acuerdo. ¡No sabemos qué está haciendo JMP! :)

Si los coeficientes estimados y los errores estándar son los mismos, entonces eres bueno. Si los coeficientes no son los mismos, no olvide que JMP podría elegir un resultado de referencia diferente para calcular los coeficientes. multinom()toma diferentes decisiones mlogit(), por ejemplo.

Obtener valores p del resultado resumen () de multinom () es bastante fácil. No puedo reproducir sus modelos, así que aquí está el ejemplo de la página de ayuda en multinom ():

library("nnet")
data("Fishing", package = "mlogit")
fishing.mu <- multinom(mode ~ income, data = Fishing)
sum.fishing <- summary(fishing.mu) # gives a table of outcomes by covariates for coef and SE
str(sum.fishing)
# now get the p values by first getting the t values
pt(abs(sum.fishing$coefficients / sum.fishing$standard.errors),
  df=nrow(Fishing)-6,lower=FALSE)

¡Estoy de acuerdo en que descubrir el paquete mlogit es un desafío! Lee las viñetas con cuidado. Ellos ayudan

atiretoo
fuente
¿Cómo usaría las otras variables (genéricas) del Fishingconjunto de datos con multinom?
gregmacfarlane
@gmacfarlane Simplemente agregue las variables que desea a la fórmula en multinom (modo ~ ingreso + precio.playa, ...
atiretoo
@atiretoo Estaba buscando una manera de obtener mis pvalles, ¡así que gracias! pero solo para aclarar: ¿de dónde viene el 6 con el df? ¿Qué debo contar para obtener mi df? Tengo una variable continua y una variable categórica (4 categorías) en mi modelo. Entonces, ¿sería un df de 5? Además, la pesca es todo el conjunto de datos ¿verdad? Que es el tamaño de la muestra-grado de libertad.
Kerry
@atiretoo Pido disculpas, acabo de descubrir que si utilizo el paquete nnet para ejecutar la regresión logit, en realidad calcula el df efectivo y se almacena dentro del objeto nnet. Entonces, puedo extraer el df del modelo y usarlo. No he verificado si el objeto mlogit tendrá la misma información.
Kerry
@Kerry, me alegro de que lo hayas encontrado.
atiretoo
1

También puede intentar ejecutar un logit multinomial usando el paquete glmnet. No estoy seguro de cómo forzarlo a mantener todas las variables, pero estoy seguro de que es posible.

Zach
fuente