Stable Log-Determinants and Gaussian Likelihood Evaluation

Why Evaluate Gaussian Likelihoods Through a Factorization?

A Gaussian likelihood combines two quantities that can be numerically delicate: a covariance determinant and a quadratic form involving the inverse covariance. Computing either by explicitly forming a determinant or an inverse is usually unnecessary and can amplify rounding error. For a symmetric positive-definite covariance matrix, a Cholesky factorization gives a direct route to both quantities using logarithms and triangular solves.

The goal here is likelihood evaluation for specified parameter values, not optimization over those values. We will use a Cholesky factor to compute the likelihood and to identify signs that a covariance estimate is close to singularity.

Likelihood Evaluation from a Cholesky Factor

Suppose \(x\in\mathbb{R}^p\) is one observation from \(N(\mu,\Sigma)\), where the specified covariance \(\Sigma\) is symmetric positive definite. Write \(r=x-\mu\). The multivariate normal log-density is

Factor \(\Sigma=LL^T\). Compute the log-determinant from the diagonal of \(L\). For the quadratic form, solve \(Ly=r\) and take the squared Euclidean length of \(y\). Neither calculation requires forming \(\Sigma^{-1}\) or \(|\Sigma|\) directly.

Proof. The determinant product rule gives \(|\Sigma|=|L||L^T|=|L|^2\). Since \(L\) is triangular, \(|L|=\prod_{j=1}^p L_{jj}\); each diagonal entry is positive, so taking logarithms yields the stated log-determinant identity. For the quadratic form, let \(u\) solve \(\Sigma u=r\). Since \(r=Ly\), the vector \(u=L^{-T}y\) satisfies \(\Sigma u=LL^TL^{-T}y=Ly=r\). Therefore \(r^T\Sigma^{-1}r=r^Tu=(Ly)^TL^{-T}y=y^Ty\). The calculation of \(y\) uses a forward triangular solve; this proof does not require computing an inverse in the numerical algorithm. \(\square\)

The word “whitened” refers here to the transformation \(y=L^{-1}r\): its squared length is the covariance-scaled squared distance of the residual. In implementation, obtain \(y\) by solving a triangular system rather than forming \(L^{-1}\).

Multiple Observations and a Fully Derived Formula

For \(n\) independent observations \(x_1,\ldots,x_n\) from the same \(N(\mu,\Sigma)\), the log-likelihood is the sum of the individual log-densities. Let \(r_i=x_i-\mu\), and let \(y_i\) solve \(Ly_i=r_i\). Summing the one-observation formula gives

Derivation. Independence means the joint density is the product of the \(n\) individual densities, so its log is their sum. Each observation contributes \(-p\log(2\pi)/2\) and \(-\log|\Sigma|/2=-\sum_j\log L_{jj}\). Its remaining contribution is \(-y_i^Ty_i/2\), by the triangular-solve identity. Summing these three terms over \(i=1,\ldots,n\) yields the displayed expression. \(\square\)

In computation, place the residuals in a matrix \(R=[r_1,\ldots,r_n]\) and solve \(LY=R\) with one triangular-system routine. The final quadratic contribution is \(-\|Y\|_F^2/2\), where \(\|Y\|_F^2\) is the sum of squares of all entries of \(Y\). This is the same calculation as summing \(y_i^Ty_i\) one observation at a time.

Near-Singularity: What the Factorization Can Tell You

A covariance matrix near singularity has one or more very small eigenvalues relative to its overall scale. In the Gaussian likelihood, this can make the log-determinant strongly negative and make quadratic forms vary sharply with the direction of the residual. Cholesky factorization is useful for evaluating a positive-definite covariance, but a successful factorization alone does not guarantee that the covariance is well-conditioned or statistically well determined.

The identity \(\log|\Sigma|=\sum_{j=1}^p\log\lambda_j\), where \(\lambda_j\) are the positive eigenvalues, shows why a very small eigenvalue matters: its logarithm can make a large negative contribution. The condition number in the Euclidean norm is \(\lambda_{\max}/\lambda_{\min}\). A large value signals that small perturbations to the matrix or arithmetic can have substantial effects on some computed quantities. Cholesky diagonal entries can help reveal a problematic factorization, but they depend on variable ordering and are not, by themselves, a condition-number estimate.

If the model uses different measurement scales, scaling can make computation and interpretation clearer. For a positive diagonal matrix \(D\), write \(\Sigma=DRD\), where \(R=D^{-1}\Sigma D^{-1}\). Then \(\log|\Sigma|=2\sum_j\log D_{jj}+\log|R|\), and the quadratic form can be evaluated by solving with the correspondingly scaled residual. This is a change of scale, not a repair for a genuinely near-singular covariance.

When a factorization fails, first check whether the input is symmetric to the accuracy expected from its construction, and inspect its eigenvalues or a condition estimate. Small asymmetries caused by rounding may sometimes be removed by replacing the matrix with \((\Sigma+\Sigma^T)/2\), but this should follow a diagnostic check rather than conceal a substantive problem. Adding a small multiple of the identity (“jitter”) changes the covariance model and therefore changes the likelihood. It should be reported and justified, not treated as an invisible numerical fix. A singular covariance requires a different, support-specific treatment; it does not have an ordinary full-dimensional Gaussian density with respect to Lebesgue measure.

Common Pitfalls and Practical Checks

• Forming an inverse to get the quadratic form. This adds needless computation and can worsen numerical error. Solve \(Ly=r\) and use \(y^Ty\).
• Taking the logarithm of a computed determinant. A determinant can overflow, underflow, or lose precision. Sum \(2\log L_{jj}\) instead.
• Ignoring scaling in the log-determinant. The formula has a factor of two: each Cholesky diagonal contributes \(2\log L_{jj}\).
• Treating Cholesky failure as harmless. Failure can indicate non-positive definiteness, material asymmetry, or a serious conditioning problem. Diagnose before modifying the matrix.
• Assuming a finite likelihood means a reliable covariance estimate. A successful factorization evaluates the supplied matrix; it does not establish that the matrix is well estimated or numerically stable.