Modelo probit con ruido gaussiano

Question

Supongamos que tenemos el siguiente modelo $$\Phi^{-1}(D)=\alpha+\beta X+\epsilon$$ donde $\epsilon\sim N(0,\sigma^{2})$ y $D_{i}=\{0,1\}$ . Esto implica que $$\text{Pr}(D_{i}=1\,|\,X,\epsilon)=\Phi(\alpha+\beta X+\epsilon)$$ pero tenemos que integrar el $\epsilon$ $$\begin{align} \text{Pr}(D_{i}=1\,|\,X)&=\int_{\epsilon}\text{Pr}(D_{i}=1\,|\,X,\epsilon)\,f_{\epsilon}(\epsilon)\,d\epsilon\\ &=\Phi\bigg(\frac{\alpha+\beta x}{\sqrt{1+\sigma^{2}}}\bigg) \end{align}$$

Así que si tuviéramos que estimar los parámetros de este modelo $(\hat{\alpha},\hat{\beta},\hat{\sigma})$ podríamos hacerlo cambiando la función de verosimilitud del modelo probit simple de $$\begin{align} L=\sum_{i=1}^{n}D_{i}\log(\Phi(\alpha+\beta X+\epsilon))+(1-D_{i})\log(\Phi(\alpha+\beta X+\epsilon)) \end{align}$$ a $$\begin{align} L^{*}=\sum_{i=1}^{n}D_{i}\log\bigg(\Phi\bigg(\frac{\alpha+\beta X}{\sqrt{1+\sigma^{2}}}\bigg)\bigg)+(1-D_{i})\log\bigg(1-\Phi\bigg(\frac{\alpha+\beta X}{\sqrt{1+\sigma^{2}}}\bigg)\bigg) \end{align}$$ Sin embargo, he observado que realizar una optimización fiable de esta probabilidad es difícil. Dado el ejemplo de juguete

n = 10000
a = -2
b = 0.01
x = runif(n, min = 1, max = 5) + rnorm(n, 0, 0.15)
p = pnorm(a + b*x)
d = rbinom(n, size = 1, prob = p)
y = tibble::as_tibble(data.frame(x, p, d))

y la probabilidad definida como

fn = function(par, x, d) {
  return(-sum(
    d*log(pmax(10^-23, pnorm((par[1] + par[2]*x)/sqrt(1 + par[3]^2)))) +
      (1-d)*log(pmax(10^-23, 1 - pnorm((par[1] + par[2]*x)/sqrt(1 + par[3]^2))))
  ))
}

y utilizando métodos cuasi-Newton

optim(par = c(0, 0, 0.5),
      fn = fn,
      x = y$x,
      d = y$d,
      method = "L-BFGS-B",
      lower = c(-Inf, -Inf, 0),
      upper = c(Inf, Inf, Inf),
      hessian = TRUE)

normalmente no se comporta muy bien. De hecho, el $\hat{\sigma}$ normalmente converge a un punto cercano al valor inicial.

¿Hay algún cambio obvio (elección del algoritmo, aproximaciones a la función de verosimilitud, mejor elección de los valores de partida) que pueda hacerse para que la estimación de $(\hat{\alpha},\hat{\beta},\hat{\sigma})$ ¿más fiable?

Answer 1

El modelo no está identificado, lo que significa que no existe una solución única para el problema de optimización. Hay infinitos valores de los parámetros que darán la misma probabilidad. Por ejemplo, $\alpha = .5$ , $\beta = 1$ y $\sigma = 2$ dará exactamente la misma probabilidad que $\alpha = 1$ , $\beta = 2$ y $\sigma = \sqrt{19}$ .

De forma más general, considere el máximo de la probabilidad, $L^*$ que se encuentra cuando $(\alpha, \beta, \sigma) = (\alpha^*, \beta^*, \sigma^*)$ . Para cualquier $k$ , $$\left(k\alpha^*, k\beta^*, \sqrt{k^2+1+(k\sigma^*)^2}\right)$$ dará exactamente la misma probabilidad. Por lo tanto, no existe un valor único de los parámetros que maximice la probabilidad. Por eso la optimización es inestable; cualquier solución específica a la que llegue se deberá puramente a la inestabilidad numérica.

Obsérvese que en una regresión probit típica, suponemos que $\sigma=0$ es decir, que no existe una variable latente $\epsilon$ que está en paradero desconocido. Esto difiere de la formulación de variable latente de la regresión probit, en la que suponemos que $$Y^*=X\beta+\epsilon$$ donde $\epsilon \sim N(0, 1)$ y $$P(D=1|X) = P(Y^*>0|X)=P(X\beta + \epsilon > 0)$$ lo que implica $P(D=1|X) = \Phi(X\beta)$ .