1
$\begingroup$

I need to compute the log-likelihood function in a high-dimensional Gaussian time-series. I have the following model:

$\mathbf{y}_{t}\left|\mathcal{F}_{t-1}\sim\mathcal{N}\left(\mathbf{\boldsymbol{\mu}}_{t},\sigma^{2}I_{n_{t}\times n_{t}}+\boldsymbol{\Sigma}_{t}\boldsymbol{\Sigma}_{t}^{'}\right)\right.$

here $\sigma^{2}$ is a constant scaler. $\boldsymbol{\Sigma}_{t}$ is a time-varying $n_{t}\times1$ vector where $n_{t}$ is big. $I_{n_{t}\times n_{t}}$ is the identity matrix of dimension $n_{t}$.

The log-likelihood will then be given by:

$\log p\left(\mathbf{y}_{1},...,\mathbf{y}_{T}\right)=\sum_{t=1}^{T}-\frac{1}{2}\left(n_{t}\log\left(2\pi\right)+\log\left(\left|\sigma^{2}I_{n_{t}\times n_{t}}+\boldsymbol{\Sigma}_{t}\boldsymbol{\Sigma}_{t}^{'}\right|\right)+\left(\mathbf{y}_{t}-\mathbf{\boldsymbol{\mu}}_{t}\right)^{'}\left(\sigma^{2}I_{n_{t}\times n_{t}}+\boldsymbol{\Sigma}_{t}\boldsymbol{\Sigma}_{t}^{'}\right)^{-1}\left(\mathbf{y}_{t}-\mathbf{\boldsymbol{\mu}}_{t}\right)\right)$

Given the special structure of the variance matrix

$\sigma^{2}I_{n_{t}\times n_{t}}+\boldsymbol{\Sigma}_{t}\boldsymbol{\Sigma}_{t}^{'}$

is there anyway to exploit this structure to evaluate the inverse, determinant and in turn the full log-likelihood computationally fast?

$\endgroup$
1
$\begingroup$

To calculate likelihood you need to calculate the determinant $|\sigma^2 I + \Sigma \Sigma^T|$ and the quadratic function $(\mathbf{y} - \boldsymbol{\mu}_t)^{\mathrm T} (\sigma^2 I + \Sigma \Sigma^{\mathrm T})^{-1} (\mathbf{y} - \boldsymbol{\mu}_t)$.

It holds that $|A| = \prod_{i}^{n_t} a_i$, where $a_i$ are eigenvalues of $A$. Matrix $\sigma^2 I + \Sigma \Sigma^{\mathrm T}$ has $n_t - 1$ eigenvalues $\sigma^2$ and one eigenvalue $\sigma^2 + \Sigma^{\mathrm T} \Sigma$ (see for example this question). So, $$|\sigma^2 I + \Sigma \Sigma^T| = \sigma^{2(n_t - 1)} (\sigma^2 + \Sigma^{\mathrm T} \Sigma).$$ Complexity of determinant calculation is $O(n_t)$.

To get inverse of the matrix $\sigma^2 I + \Sigma \Sigma^{\mathrm T}$ you should use Sherman–Morrison formula: $$(\sigma^2 I + \Sigma \Sigma^{\mathrm T})^{-1} = \frac{1}{\sigma^2} I - \frac{1}{\sigma^4} \frac{\Sigma \Sigma^{\mathrm T}}{1 + \frac{1}{\sigma^2} \Sigma^{\mathrm T} \Sigma}.$$ To further reduce computation cost you can calculate $\xi = \Sigma^T (\mathbf{y} - \boldsymbol{\mu}_t)$, as $$(\mathbf{y} - \boldsymbol{\mu}_t)^{\mathrm T} (\sigma^2 I + \Sigma \Sigma^{\mathrm T})^{-1} (\mathbf{y} - \boldsymbol{\mu}_t) = \frac{1}{\sigma^2} (\mathbf{y} - \boldsymbol{\mu}_t)^{\mathrm T} (\mathbf{y} - \boldsymbol{\mu}_t) - \frac{1}{\sigma^4} \frac{\xi^2}{1 + \frac{1}{\sigma^2} \Sigma^{\mathrm T} \Sigma}.$$ As we need only to multiply vectors of size $n_t$ computational complexity is again $O(n_t)$.

So, we get the results for one $\boldsymbol{\mu}_t$ using $O(n_t)$ operations. For different $t$ we then get total computational cost $O(\sum_{t = 1}^T n_t)$.

$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.