user3390169
user3390169

Reputation: 1045

Cholesky decomposition failure for my correlation matrix

I am trying to use chol() to find the Cholesky decomposition of the correlation matrix below. Is there a maximum size I can use that function on? I am asking because I get the following:

d <-chol(corrMat)
Error in chol.default(corrMat) : 
  the leading minor of order 61 is not positive definite

but, I can decompose it for less than 60 elements without a problem (even when it contains the 61st element of the original):

> d  <-chol(corrMat[10:69, 10:69])
> d  <-chol(corrMat[10:70, 10:70])
Error in chol.default(corrMat[10:70, 10:70]) : 
  the leading minor of order 61 is not positive definite

Here is the matrix:

https://drive.google.com/open?id=0B0F1yWDNKi2vNkJHMDVHLWh4WjA

Upvotes: 1

Views: 2543

Answers (1)

Zheyuan Li
Zheyuan Li

Reputation: 73315

The problem is not size, but numerical rank!

d <- chol(corrMat, pivot = TRUE)

dim(corrMat)
#[1] 72 72

attr(d, "rank")
#[1] 62

corrMat is not positive-definite. Ordinary Cholesky factorization will fail, but pivoted version works.

The correct Cholesky factor here can be obtained (see Correct use of pivot in Cholesky decomposition of positive semi-definite matrix)

r <- attr(d, "rank")
reverse_piv <- order(attr(d, "pivot"))
d[-(1:r), -(1:r)] <- 0
R <- d[, reverse_piv]

Whether this is acceptable depends on your context. It might need corresponding adjustment to your other code.

Pivoted Cholesky factorization can do many things that sound impossible for a deficient, non-invertible covariance matrix, like

Upvotes: 2

Related Questions