7
$\begingroup$

I have $n$, $p$-dimensional vectors and I am constructing the $p \times p$ covariance matrix using the following formula:

$$\mathrm{Cov}(j,k) = \frac{1}{n-1} {\sum^n_{i=1}} (x_i(j) - {\mu}(j))(x_i(k) - {\mu(k)}) $$

However, the values in the $p$-dimensional vectors are between 0 and 1, as it is a normalized vector, and $n$ can be large which results in $\mathrm{Cov}(j,k)$ having rather small elements.

I am trying to calculate the determinant of the covariance matrix to determine the differential entropy of the underlying multivariate Gaussian. However, because of the small numbers, the determinant is effectively 0 (and in fact rounds to 0 in computation).

Are there any methods to avoid this? Perhaps I can avoid the division by $(n-1)$?

$\endgroup$

2 Answers 2

6
$\begingroup$

Several things.

  • The determinant is scale equivariant, so without loss of generality you can rescale your $x_i$'s (say multiply them by 10) and then factor back the scale after all the intermediate operations have been performed as shown in the following R code:

    p <- 10
    n <- 100
    a1 <- matrix(runif(n*p),n,p)
    a1 <- sweep(a1,1,sqrt(rowSums(a1*a1)),FUN="/")
    a2 <- a1*10 #you can re-scale your data
    log(det(cov(a1[1:(p+1),])))
    log(det(cov(a2[1:(p+1),])))-log(10^(2*(p))) #and compensate for the rescaling 'at the end'
    
  • You should work with the determinant rescaled to be expressed in the original units of your dataset (e.g. if $\varSigma$ is your covariance matrix and $p$ it's rank, use $|\varSigma|^{1/2p}$). This is similar to how the standard deviations are expressed in the original units but the variances in the square of the original units.

$\endgroup$
9
  • $\begingroup$ What is matrix(runif(n*p),n,p)? $\endgroup$
    – Aly
    Feb 19, 2013 at 11:28
  • $\begingroup$ @Aly: your data (a bunch of random points between 0 and 1). The next line normalizes it... $\endgroup$
    – user603
    Feb 19, 2013 at 11:40
  • $\begingroup$ so essentially $log(det(cov)) == log((1/{10^2})*det(cov*10)) == log(det(cov*10)) - 2log(10)$ ? $\endgroup$
    – Aly
    Feb 19, 2013 at 11:45
  • $\begingroup$ @Aly: you forgot the p! $\endgroup$
    – user603
    Feb 19, 2013 at 11:46
  • 1
    $\begingroup$ Ah yes, I almost forgot what I was trying to do and just focusing on the equation. I understand now $\endgroup$
    – Aly
    Feb 19, 2013 at 12:09
2
$\begingroup$

Numerical stability will be a problem for you.

One approach to make the calculation a little more stable - taking $X$ to be a mean-corrected matrix of the $x$-values, I think it makes sense to either form a QR decomposition of the $X$-matrix or failing that a Choleski decomposition of the $X'X$ matrix, and use the product of the diagonal of $R$ (or of the Choleski factor) which will be the square root of the determinant.

There are other approaches (some of which may be better than computing the QR).

You still need to scale for the $n-1$'s, though.

$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

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