¿Cómo garantizar propiedades de la matriz de covarianzas al modelo normal multivariante usando máxima verosimilitud?

Question

Supongamos que tengo el siguiente modelo

$$y_i=f(x_i,\theta)+\varepsilon_i$$

where $y_i\in \mathbb{R}^K$ , $x_i$ is a vector of explanatory variables, $\theta$ is the parameters of non-linear function $f$ and $\varepsilon_i\sim N(0,\Sigma)$, where $\Sigma$ naturally is $K\times K$ matrix.

The goal is the usual to estimate $\theta$ and $\Sigma$. The obvious choice is maximum likelihood method. Log-likelihood for this model (assuming we have a sample $(y_i,x_i),i=1,...,n$) looks like

$$l(\theta,\Sigma)=-\frac{n}{2}\log(2\pi)-\frac{n}{2} \log\det\Sigma-\sum_{i=1}^n(y_i-f(x_i,\theta))'\Sigma^{-1}(y-f(x_i,\theta)))$$

Now this seems simple, the log-likelihood is specified, put in data, and use some algorithm for non-linear optimisation. The problem is how to ensure that $\Sigma$ is positive definite. Using for example optim in R (or any other non-linear optimisation algorithm) will not guarantee me that $\Sigma$ is positive definite.

So the question is how to ensure that $\Sigma$ stays positive definite? I see two possible solutions:

1. Reparametrise $\Sigma$ as $RR'$ where $R$ is upper-triangular or symmetric matrix. Then $\Sigma$ will always be positive-definite and $R$ puede ser sin restricciones.

2. Perfil de uso de la probabilidad. Derivar las fórmulas para$\hat\theta(\Sigma)$$\hat{\Sigma}(\theta)$. Empezar con algo de $\theta_0$ e iterar $\hat{\Sigma}_j=\hat\Sigma(\hat\theta_{j-1})$, $\hat{\theta}_j=\hat\theta(\hat\Sigma_{j-1})$ hasta la convergencia.

¿Hay alguna otra manera y lo que acerca de estos 2 métodos, van a trabajar, son estándar? Esto parece bastante estándar problema, pero una búsqueda rápida no me dan los punteros. Sé que la estimación Bayesiana también sería posible, pero por el momento no quisiera participar en ella.

Answer 1

¿Suponiendo que en la construcción de la matriz de covarianza, automáticamente está cuidando de la cuestión de la simetría, la log-verosimilitud será $-\infty$ cuando $\Sigma$ no es positiva definida por el término de $\log {\rm det} \ \Sigma$ en el modelo adecuado? Para evitar un error numérico si ${\rm det} \ \Sigma < 0$ quisiera precalculate ${\rm det} \ \Sigma$ y, si no es positivo, entonces la probabilidad de registro igual -Inf, lo contrario. Tienes que calcular el determinante de todos modos, así que esto no le cuesta ningún cálculo extra.

Answer 2

Resulta que usted puede utilizar el perfil de máxima verosimilitud para garantizar las propiedades necesarias. Se puede demostrar que para un determinado $\hat\theta$, $l(\hat\theta,\Sigma)$ es maximizada por

$$\hat\Sigma=\frac{1}{n}\sum_{i=1}^n\hat{\varepsilon}_i\hat{\varepsilon}_i',$$

donde

$$\hat{\varepsilon}_i=y_i-f(x_i,\hat\theta)$$

Entonces es posible demostrar que

$$\sum_{i=1}^n(y_i-f(x_i,\hat\theta))'\hat\Sigma^{-1}(y-f(x_i,\hat\theta)))=const,$$

por lo tanto sólo tenemos que maximizar

$$l_R(\theta,\Sigma)=-\frac{n}{2} \log\det\hat\Sigma.$$

Naturalmente, en este caso $\Sigma$ va a satisfacer todas las propiedades necesarias. Las pruebas son las mismas para el caso de al $f$ es lineal, el cual puede ser encontrado en el Análisis de Series de Tiempo por J. D. Hamilton página 295, de ahí que se omite.

Answer 3

Una alternativa parametrización de la matriz de covarianza es en términos de autovalores $\lambda_1,...,\lambda_p$ $p(p-1)/2$ "Dados" ángulos $\theta_ij$.

Es decir, podemos escribir

$$\Sigma = G^T \Lambda G$$

donde $G$ es ortonormales, y

$$\Lambda = diag(\lambda_1, ..., \lambda_p)$$

con $\lambda_1 \geq ... \geq \lambda_p \geq 0$.

Mientras tanto, $G$ puede ser parametrizado de forma única en términos de $p(p-1)/2$ ángulos, $\theta_{ij}$ donde$i = 1,2,...,p-1$$j = i, ..., p-1$.[1]

(detalles)

[1]: Hoffman, Raffenetti, Ruedenberg. "La generalización de los Ángulos de Euler con N‐Dimensional Ortogonal de Matrices". J. Math. Phys. 13, 528 (1972)

Answer 4

A lo largo de las líneas de solución de charles.y.zheng, usted tal vez desee modelo $\Sigma = \Lambda + C C^{\top}$, donde $\Lambda$ es una matriz diagonal, y $C$ es una factorización de Cholesky de una actualización de fila a $\Lambda$. Sólo deberá mantener la diagonal de $\Lambda$ positiva para mantener $\Sigma$ positiva definida. Es decir, debe calcular la diagonal de $\Lambda$ y los elementos de $C$ en vez de estimar $\Sigma$.