Interpretación de la regresión logística ordinal

Question

Ejecuté esta regresión logística ordinal en R:

mtcars_ordinal <- polr(as.factor(carb) ~ mpg, mtcars)

Tengo este resumen del modelo:

summary(mtcars_ordinal)

Re-fitting to get Hessian

Call:
polr(formula = as.factor(carb) ~ mpg, data = mtcars)

Coefficients:
      Value Std. Error t value
mpg -0.2335    0.06855  -3.406

Intercepts:
    Value   Std. Error t value
1|2 -6.4706  1.6443    -3.9352
2|3 -4.4158  1.3634    -3.2388
3|4 -3.8508  1.3087    -2.9425
4|6 -1.2829  1.3254    -0.9679
6|8 -0.5544  1.5018    -0.3692

Residual Deviance: 81.36633 
AIC: 93.36633 

Puedo obtener las probabilidades logarítmicas del coeficiente para mpg así:

exp(coef(mtcars_ordinal))
 mpg 
0.7917679 

Y las probabilidades logarítmicas de los umbrales como:

exp(mtcars_ordinal$zeta)

       1|2         2|3         3|4         4|6         6|8 
0.001548286 0.012084834 0.021262900 0.277242397 0.574406353 

¿Podría alguien decirme si mi interpretación de este modelo es correcta?

En mpg aumenta en una unidad, las probabilidades de pasar de la categoría 1 de carb en cualquiera de las otras 5 categorías, disminuye en -0,23. Si las probabilidades logarítmicas cruzan el umbral de 0,0015, entonces el valor predicho para un coche será la categoría 2 de carb . Si las probabilidades logarítmicas cruzan el umbral de 0,0121, entonces el valor predicho para un coche será la categoría 3 de carb etc.

Answer 1

Has confundido perfectamente las probabilidades y las probabilidades logarítmicas. Las probabilidades logarítmicas son los coeficientes; las probabilidades son coeficientes exponenciados. Además, la interpretación de las probabilidades va al revés . (Crecí con la econometría pensando en las variables dependientes limitadas, y la interpretación de probabilidades de la regresión ordinal me resulta... uhm... divertida). Así que tu primera afirmación debería decir, "Como mpg aumenta en una unidad, las probabilidades de observar la categoría 1 de carb frente a otras 5 categorías aumentan un 21%".

En cuanto a la interpretación de los umbrales, en realidad hay que trazar todas las curvas de predicción para poder decir cuál es la predicción modal:

mpg   <- seq(from=5, to=40, by=1)
xbeta <- mpg*(-0.2335)
logistic_cdf <- function(x) {
  return( 1/(1+exp(-x) ) )
}

p1 <- logistic_cdf( -6.4706 - xbeta )
p2 <- logistic_cdf( -4.4158 - xbeta ) - logistic_cdf( -6.4706 - xbeta )
p3 <- logistic_cdf( -3.8508 - xbeta ) - logistic_cdf( -4.4158 - xbeta )
p4 <- logistic_cdf( -1.2829 - xbeta ) - logistic_cdf( -3.8508 - xbeta )
p6 <- logistic_cdf( -0.5544 - xbeta ) - logistic_cdf( -1.2829 - xbeta )
p8 <- 1 - logistic_cdf( -0.5544 - xbeta )

plot(mpg, p1, type='l', ylab='Prob')
  lines(mpg, p2, col='red')
  lines(mpg, p3, col='blue')
  lines(mpg, p4, col='green')
  lines(mpg, p6, col='purple')
  lines(mpg, p8, col='brown')
  legend("topleft", lty=1, col=c("black", "red", "blue", "green", "purple", "brown"), 
         legend=c("carb 1", "carb 2", "carb 3", "carb 4", "carb 5", "carb 6"))

La curva azul de la 3ª categoría no repuntó en ningún momento, como tampoco lo hizo la curva morada de la 6ª categoría. Así que, en todo caso, yo diría que para valores de mpg por encima de 27 tienen, la categoría más probable es 1; entre 18 y 27, categoría 2; entre 4 y 18, categoría 4; y por debajo de 4, categoría 8. (Me pregunto qué es lo que estás estudiando: ¿camiones comerciales? La mayoría de los turismos de hoy en día deberían tener mpg > 25). Quizá quieras intentar determinar los puntos de intersección con más precisión.

También me he dado cuenta de que tienes esas categorías raras que van 1, 2, 3, 4, luego 6 (saltándose 5), luego 8 (saltándose 7). Si 5 y 7 faltaban por diseño, eso está bien. Si se trata de categorías válidas que carb simplemente no cae en, esto no es bueno.

Answer 2

En el modelo logit ordenado, las probabilidades forman el cociente entre la probabilidad de estar en cualquier categoría por debajo de un umbral específico frente a la probabilidad de estar en una categoría por encima del mismo umbral (por ejemplo, con tres categorías: Probabilidad de estar en la categoría A o B frente a C, así como la probabilidad de estar en la categoría A frente a B o C).

Esto conduce al modelo logit P(Y <= k | x) = zeta_k - eta como se especifica en la descripción de polr() . Por lo tanto, las odds ratio pueden construirse para diferentes categorías o para diferentes regresores. Esta última, la más común, compara las probabilidades para las mismas categorías pero diferentes regresores e iguala a

$$\newcommand{\odds}{{\rm odds}} \frac{\odds(y_a \le k \,|\,x_a)}{\odds(y_b \le k \,|\,x_b)}~=~ \exp(-(\eta_a - \eta_b)).$$

La odds ratio para las diferentes categorías se define como

$$\frac{\odds(y_i \le k \,|\,x_i)}{\odds(y_i \le m \,|\,x_i)}~=~ \exp(\zeta_k - \zeta_m),$$

por lo que el coeficiente es independiente de los regresores. Esta propiedad da lugar a la denominación alternativa modelo de probabilidades proporcionales.

En este ejemplo sencillo, pero quizá poco intuitivo, se podría formular: Para un aumento de una unidad en el regresor mpg las probabilidades de observar la categoría 1 frente a observar cualquier categoría superior (o las probabilidades de observar cualquier categoría por debajo de un determinado umbral frente a observar cualquier categoría por encima del mismo umbral) se multiplican por 1,26 o aumentan un 26% ( exp(-(-0.233 - 0)) = 1.263 ). Si desea formular una razón de probabilidades de diferentes categorías, podría, por ejemplo, decir que las probabilidades de estar en la categoría 1 frente a cualquier categoría anterior en comparación con las probabilidades de estar en la categoría 1 o 2 frente a cualquier categoría anterior es igual a exp((-6.470) - (-4.415)) = 0.128 . Esta última interpretación no es muy útil en este caso concreto. Un ejemplo de odds ratio para diferentes categorías podría ser la probabilidad de ir a la universidad comparada con la probabilidad de ir al instituto.

Por último, podría interesarle cuánto debe cambiar una variable explicativa para alcanzar la categoría de respuesta inmediatamente superior. Para ello se compara la longitud del intervalo $(\zeta_k - \zeta_{k-1})$ con un coeficiente ajustado. Esto da una idea de lo grande que debe ser el cambio en su respectivo regresor para mover la respuesta de la categoría $k$ a la categoría superior.