Derivación de la estimación de máxima verosimilitud (MLE) de un parámetro para una distribución gaussiana inversa

Question

Dada la siguiente función de probabilidad

$$f(y|x,\tau) = \prod_{i=0}^Nf_T(u_i-x_i-\tau) \tag{1}$$

donde, $f_T(t)$ es la función de densidad de probabilidad de un Distribución gaussiana inversa dado por

$$f_T(t) = \sqrt\frac{\lambda}{2\pi t^3} \exp\Bigl(- \frac{\lambda (t-\mu)^2}{2\mu^2t}\Bigr)\tag{2}$$

El objetivo aquí es determinar la MLE del parámetro $\tau$

$$ \hat{\tau}_{MLE} := \mathop{argmax}\limits_\tau f(y|x,\tau) \tag{3}$$

Según el principio de MLE y sustituyendo $(2)$ en $(1)$ obtendremos lo siguiente

\begin {align}L( \tau ) & = \prod_ {i=1}^N \sqrt\frac { \lambda }{2 \pi (u_i-x_i- \tau )^3} \exp\Bigl (- \frac { \lambda (u_i-x_i- \tau - \mu )^2}{2 \mu ^2(u_i-x_i- \tau )} \Bigr ) \\\\ & = \Bigl ( \frac { \lambda }{2 \pi } \Bigr )^{N/2} \prod_ {i=1}^N(u_i-x_i- \tau )^{-3/2} \exp\Bigl (- \frac { \lambda }{2 \mu ^2} \sum_ {i=1}^N \frac {(u_i-x_i- \tau - \mu )^2}{u_i-x_i- \tau } \Bigr ) \tag {4} \end {align}

Tomando el logaritmo, obtenemos

\begin {alinear} logL( \tau ) & = \frac {N}{2} log \Bigl ( \frac { \lambda }{2 \pi } \Bigr ) - \frac {3}{2} \sum_ {i=1}^N \log (u_i-x_i- \tau ) - \frac { \lambda }{2 \mu ^2} \sum_ {i=1}^N \frac {(u_i-x_i- \tau - \mu )^2}{u_i-x_i- \tau } \tag {5} \end {align}

Ahora tomando la deriativa con respecto a. $\tau$

\begin {align} \frac {d(logL( \tau ))}{d \tau }& = 0 - \frac {3}{2} \sum_ {i=1}^N \frac {1} {(u_i-x_i- \tau )}(-1) - \frac { \lambda }{2 \mu ^2} \sum_ {i=1}^N \left ( \frac {2(u_i-x_i- \tau - \mu )}{u_i-x_i- \tau }(-1) - \frac {(u_i-x_i- \tau - \mu )^2}{(u_i-x_i- \tau )^2}(-1) \right ) \\\\ & = \frac {3}{2} \sum_ {i=1}^N \frac {1} {(u_i-x_i- \tau )}- \frac { \lambda }{2 \mu ^2} \sum_ {i=1}^N \left ( \frac {-2(u_i-x_i- \tau - \mu )}{u_i-x_i- \tau } + \frac {(u_i-x_i- \tau - \mu )^2}{(u_i-x_i- \tau )^2} \right ) \tag {6} \end {align}

Ecuación de ajuste $6$ a $0$

\begin {align} \frac {3}{2} \sum_ {i=1}^N \frac {1} {(u_i-x_i- \tau )}- \frac { \lambda }{2 \mu ^2} \sum_ {i=1}^N \left ( \frac {-2(u_i-x_i- \tau - \mu )}{u_i-x_i- \tau } + \frac {(u_i-x_i- \tau - \mu )^2}{(u_i-x_i- \tau )^2} \right )= 0 \tag {7} \end {align}

Antes de llegar al problema en cuestión, ¿son correctas las derivaciones realizadas hasta ahora?

Aquí está el cuello de botella:

¿Cómo debo proceder a partir de aquí? El segundo término de la suma se ha vuelto muy complicado y no puedo averiguar cómo derivar $\tau$ .

[ACTUALIZACIÓN 2] según las aportaciones de @gunes

$=>\frac{3}{2}\sum_{i=1}^N \frac{1} {(u_i-x_i-\tau)}- \frac{\lambda }{2\mu^2} \sum_{i=1}^N \left(-1+1^2 - 2*1*\left(\frac{u_i-x_i-\tau-\mu}{u_i-x_i-\tau}\right) + \left(\frac{u_i-x_i-\tau-\mu}{u_i-x_i-\tau} \right)^2\right)= 0 $

$=>\frac{3}{2}\sum_{i=1}^N \frac{1} {(u_i-x_i-\tau)}- \frac{\lambda }{2\mu^2} \sum_{i=1}^N -1+\left(\frac{\require{cancel} \cancel{u_i}-\require{cancel} \cancel{x_i}-\require{cancel} \cancel{\tau} -\require{cancel} \cancel{u_i}+\require{cancel} \cancel{x_i}+\require{cancel} \cancel{\tau}+\mu}{u_i-x_i-\tau} \right)^2= 0 $

$=>\frac{3}{2}\sum_{i=1}^N \frac{1} {(u_i-x_i-\tau)}- \frac{\lambda }{2\mu^2} \sum_{i=1}^N -1+\left(\frac{\mu}{u_i-x_i-\tau} \right)^2= 0 $

$=>\frac{3}{2}\sum_{i=1}^N \frac{1} {(u_i-x_i-\tau)}+\frac{N\lambda}{2\mu^2 }- \frac{\lambda N\require{cancel} \cancel{\mu^2} }{2\require{cancel} \cancel{\mu^2}} \sum_{i=1}^N \frac{1}{(u_i-x_i-\tau)^2} = 0 $

$=>\frac{3}{2}\sum_{i=1}^N \frac{1} {(u_i-x_i-\tau)}- \frac{\lambda N }{2} \sum_{i=1}^N \frac{1}{(u_i-x_i-\tau)^2}= -\frac{N\lambda}{2\mu^2 } $

$=>\sum_{i=1}^N \frac{3(u_i-x_i-\tau) - \lambda}{(u_i-x_i-\tau)^2} = -\frac{N\lambda}{\mu^2 } $

$=>\sum_{i=1}^N \frac{ \lambda-3(u_i-x_i-\tau)}{(u_i-x_i-\tau)^2} = \frac{N\lambda}{\mu^2 } $

[ACTUALIZACIÓN 3] Según la derivación proporcionada por @Ben

$$1 + 3 H_{-1}(\tau)^2 H_1(\tau)^2 - 5 H_{-1}(\tau) H_1(\tau) + H_1(\tau)^2 H_{-2}(\tau) = 0.$$

Como dice @Ben más abajo, nos quedamos con la ecuación mencionada, que obviamente no es fácil de estimar $\tau$ .

Ahora nos quedan las siguientes preguntas: ¿Cómo podemos resolver esto numéricamente? ¿Existen paquetes de software que puedan realizar este tipo de soluciones numéricas? ¿O es mejor escribir una nosotros mismos?

Answer 1

La derivación completa de los MLEs para los datos de IID de un distribución gaussiana inversa se puede encontrar en la respuesta a esta pregunta relacionada . En tu caso has añadido una capa adicional de complicación al tener valores de datos observables $t_i = u_i - x_i - \tau$ que dependen de algunas covariables condicionantes y de un parámetro adicional. A partir de esta formulación, su densidad de muestreo es:

$$f(\mathbf{u} | \mathbf{x}, \tau, \mu, \lambda) = \prod_{i=1}^n \Big( \frac{\lambda}{2 \pi (u_i-x_i-\tau)^3} \Big)^{1/2} \exp \Big( - \sum_{i=1}^n \frac{\lambda (u_i-x_i-\tau - \mu)^2}{2 \mu^2 (u_i-x_i-\tau)} \Big)$$

sobre el soporte $\mathbf{u} \geqslant \mathbf{x} + \tau \mathbf{1}$ . La función de probabilidad logarítmica se define sobre $\tau \leqslant \min (u_i-x_i)$ y se da en este rango por:

$$\ell_{\mathbf{u},\mathbf{x}}(\tau, \mu, \lambda) = \text{const} + \frac{n}{2} \ln (\lambda) - \frac{3}{2} \sum_{i=1}^n \ln (u_i-x_i-\tau) - \frac{\lambda}{2 \mu^2 } \sum_{i=1}^n \frac{(u_i-x_i-\tau - \mu)^2}{(u_i-x_i-\tau)}.$$

---

Encontrar la MLE: Para facilitar nuestro análisis, definimos las funciones

$$H_k(\tau) \equiv \frac{1}{n} \sum_{i=1}^n (u_i-x_i-\tau)^k.$$

Entonces tenemos:

$$\begin{equation} \begin{aligned} \frac{\partial \ell_{\mathbf{u},\mathbf{x}}}{\partial \tau}(\tau, \mu, \lambda) &= \frac{3}{2} \sum_{i=1}^n \frac{1}{u_i-x_i-\tau} + \frac{\lambda}{2 \mu^2 } \sum_{i=1}^n \frac{(u_i - x_i - \tau + \mu)(u_i-x_i-\tau - \mu)}{(u_i-x_i-\tau)^2} \\[10pt] &= \frac{3}{2} \sum_{i=1}^n \frac{1}{u_i-x_i-\tau} + \frac{\lambda}{2 \mu^2 } \sum_{i=1}^n \frac{(u_i - x_i - \tau)^2 -2 \mu (u_i-x_i-\tau) + \mu^2}{(u_i-x_i-\tau)^2} \\[10pt] &= \frac{3}{2} \sum_{i=1}^n \frac{1}{u_i-x_i-\tau} + \frac{\lambda}{2 \mu^2 } \Big[ n - 2\mu \sum_{i=1}^n \frac{1}{u_i-x_i-\tau} + \mu^2 \sum_{i=1}^n \frac{1}{(u_i-x_i-\tau)^2} \Big] \\[10pt] &= \frac{3n}{2} H_{-1}(\tau) + \frac{n \lambda}{2 \mu^2 } \Big[ 1 - 2 \mu H_{-1}(\tau) + \mu^2 H_{-2}(\tau) \Big]. \\[10pt] \end{aligned} \end{equation}$$

Tomando $\tau$ para fijarse por el momento, el MLEs de la distribución gaussiana inversa son:

$$\hat{\mu}(\tau) = H_1(\tau) \quad \quad \quad \frac{1}{\hat{\lambda}(\tau)} = H_{-1}(\tau) - \frac{1}{H_1(\tau)}.$$

Sustituyendo estas funciones se obtiene:

$$\begin{equation} \begin{aligned} \frac{\partial \ell_{\mathbf{u},\mathbf{x}}}{\partial \tau}(\tau, \hat{\mu}(\tau), \hat{\lambda}(\tau)) &= \frac{3n}{2} H_{-1}(\tau) + \frac{n}{2 H_1(\tau)^2 } \frac{1 - 2 H_1(\tau) H_{-1}(\tau) + H_1(\tau)^2 H_{-2}(\tau)}{H_{-1}(\tau) - H_1(\tau)^{-1}} \\[10pt] &= \frac{n}{2} \cdot \frac{1}{H_1(\tau)^2} \Big[ 3 H_{-1}(\tau) H_1(\tau)^2 - \frac{2 H_1(\tau) H_{-1}(\tau) - H_1(\tau)^2 H_{-2}(\tau) - 1}{H_{-1}(\tau) - H_1(\tau)^{-1}} \Big]. \\[10pt] \end{aligned} \end{equation}$$

Al fijar esta derivada parcial en cero se obtiene la ecuación del punto crítico:

$$1 + 3 H_{-1}(\tau)^2 H_1(\tau)^2 - 5 H_{-1}(\tau) H_1(\tau) + H_1(\tau)^2 H_{-2}(\tau) = 0.$$

Esta ecuación del punto crítico deberá resolverse numéricamente, ya que no existe una expresión sencilla para la solución.