12.3 Estimación e inferencia en los modelos logit y probit

Hasta ahora no se ha dicho nada sobre cómo se estiman los modelos Logit y Probit mediante software estadístico. La razón por la que esto es interesante es que ambos modelos son no lineales en los parámetros y, por lo tanto, no pueden estimarse usando MCO. En su lugar, uno se basa en estimación de máxima verosimilitud (EMV). Otro enfoque es la estimación por mínimos cuadrados no lineales (MCNL).

Mínimos cuadrados no lineales

Considere el modelo Probit de regresión múltiple:

\[\begin{align} E(Y_i\vert X_{1i}, \dots, X_{ki}) = P(Y_i=1\vert X_{1i}, \dots, X_{ki}) = \Phi(\beta_0 + \beta_1 X_{1i} + \dots + \beta_k X_{ki}). \tag{12.8} \end{align}\]

De manera similar a MCO, MCNL estima los parámetros \(\beta_0,\beta_1,\dots,\beta_k\) minimizando la suma de errores al cuadrado

\[\sum_{i=1}^n\left[ Y_i - \Phi(b_0 + b_1 X_{1i} + \dots + b_k X_{ki}) \right]^2.\]

La estimación de MCNL es un enfoque coherente que produce estimaciones que normalmente se distribuyen en muestras grandes. En R existen funciones como NLM() del paquete stats que proporcionan algoritmos para resolver problemas de mínimos cuadrados no lineales.

Sin embargo, MCNL es ineficiente, lo que significa que existen técnicas de estimación que tienen una varianza más pequeña.

Estimación de máxima verosimilitud

En la EMV se busca estimar los parámetros desconocidos eligiéndolos de manera que se maximice la probabilidad de extraer la muestra observada. Esta probabilidad se mide mediante la función de verosimilitud, la distribución de probabilidad conjunta de los datos tratados en función de los parámetros desconocidos. Dicho de otra manera, las estimaciones de máxima verosimilitud de los parámetros desconocidos son los valores que dan como resultado un modelo que es más probable que produzca los datos observados. Resulta que EMV es más eficiente que MCNL.

Como las estimaciones de máxima verosimilitud se distribuyen normalmente en muestras grandes, la inferencia estadística de los coeficientes en modelos no lineales como la regresión Logit y Probit se puede realizar utilizando las mismas herramientas que se utilizan para los modelos de regresión lineal: Se puede calcular el estadístico \(t\) e intervalos de confianza.

Muchos paquetes de software utilizan un algoritmo EMV para la estimación de modelos no lineales. La función glm() utiliza un algoritmo llamado mínimos cuadrados iterativamente reponderados.

Medidas de ajuste

Es importante tener en cuenta que los valores habituales \(R^2\) y \(\bar{R}^2\) son inválidos para los modelos de regresión no lineal. La razón de esto es simple: Ambas medidas asumen que la relación entre la variable dependiente y explicativa(s) es lineal. Obviamente, esto no es válido para los modelos Probit y Logit. Por tanto, \(R^2\) no necesita estar entre \(0\) y \(1\) y no existe una interpretación significativa. Sin embargo, el software estadístico a veces informa estas medidas de todos modos.

Existen muchas medidas de ajuste para los modelos de regresión no lineal y no existe consenso sobre cuál debe informarse. La situación es aún más complicada porque no existe una medida de ajuste que sea significativa en general. Para modelos con una variable de respuesta binaria como \(deny\), uno podría usar la siguiente regla: Si \(Y_i = 1\) y \(\widehat{P(Y_i|X_{i1}, \dots, X_{ik})} > 0.5\) o si \(Y_i = 0\) y \(\widehat{P(Y_i|X_{i1}, \dots, X_{ik})} < 0.5\), considere el \(Y_i\) como se predijo correctamente. De lo contrario, se dice que \(Y_i\) se predice incorrectamente. La medida de ajuste es la proporción de observaciones predichas correctamente. La desventaja de este enfoque es que no refleja la calidad de la predicción: Si \(\widehat{P(Y_i = 1|X_{i1}, \dots, X_{ik}) = 0.51}\) o \(\widehat{P(Y_i =1|X_{i1}, \dots, X_{ik}) = 0.99}\) no se refleja, solo se predice \(Y_i = 1\).¹

Una alternativa a este último son las llamadas medidas pseudo-\(R^2\). Para medir la calidad del ajuste, estas medidas comparan el valor de la probabilidad maximizada (log-) del modelo con todos los regresores (el modelo completo) con la probabilidad de un modelo sin regresores (modelo nulo, regresión sobre una constante).

Por ejemplo, considere una regresión Probit. El \(\text{pseudo-}R^2\) viene dado por \[\text{pseudo-}R^2 = 1 - \frac{\ln(f^{max}_{full})}{\ln(f^{max}_{null})}\] donde \(f^{max}_j \in [0,1]\) denota la probabilidad maximizada para el modelo \(j\).

El razonamiento detrás de esto es que la probabilidad maximizada aumenta a medida que se agregan regresores adicionales al modelo, de manera similar a la disminución en \(SSR\) cuando se agregan regresores en un modelo de regresión lineal. Si el modelo completo tiene una probabilidad maximizada similar a la del modelo nulo, el modelo completo no mejora realmente sobre un modelo que usa solo la información en la variable dependiente, entonces \(\text{pseudo-}R^2 \approx 0\). Si el modelo completo se ajusta muy bien a los datos, la probabilidad maximizada debe estar cerca de \(1\), de modo que \(\ln(f^{max}_{full}) \approx 0\) y \(\text{pseudo-}R^2 \approx 1\). Consultar en internet para obtener más información sobre las medidas EMV y pseudo-\(R^2\).

summary() no reporta \(\text{pseudo-}R^2\) para modelos estimados por glm(), pero se pueden usar las entradas desviación residual (desviación) y desviación nula (desviación nula) en su lugar. Estos se calculan como

\[\text{deviance} = -2 \times \left[\ln(f^{max}_{saturated}) - \ln(f^{max}_{full}) \right]\]

y

\[\text{null deviance} = -2 \times \left[\ln(f^{max}_{saturated}) - \ln(f^{max}_{null}) \right]\]

donde \(f^{max}_{saturated}\) es la probabilidad maximizada para un modelo que asume que cada observación tiene su propio parámetro (existen \(n + 1\) parámetros para estimar que conducen a un ajuste perfecto). Para modelos con una variable dependiente binaria, se sostiene que \[\text{pseudo-}R^2 = 1 - \frac{\text{deviance}}{\text{null deviance}} = 1- \frac{\ln(f^{max}_{full})}{\ln(f^{max}_{null})}.\]

Ahora se calcula \(\text{pseudo-}R^2\) para el modelo Probit aumentado de denegación de hipoteca.

# calcular pseudo-R2 para el modelo probit de denegación de hipoteca
pseudoR2 <- 1 - (denyprobit2$deviance) / (denyprobit2$null.deviance)
pseudoR2
#> [1] 0.08594259

Otra forma de obtener \(\text{pseudo-}R^2\) es estimar el modelo nulo usando glm() y extraer las probabilidades logarítmicas maximizadas tanto para el modelo nulo como para el modelo completo usando la función logLik().

# calcular el modelo nulo
denyprobit_null <- glm(formula = deny ~ 1, 
                       family = binomial(link = "probit"), 
                       data = HMDA)

# calcular la pseudo-R2 usando 'logLik'
1 - logLik(denyprobit2)[1]/logLik(denyprobit_null)[1]
#> [1] 0.08594259

---

1. Esto contrasta con el caso de una variable dependiente numérica en la que se utilizan los errores al cuadrado para evaluar la calidad de la predicción.↩︎