필터 지우기
필터 지우기

Issue with Cholesky decomposition and positive definiteness

조회 수: 31 (최근 30일)
James Barrett
James Barrett 2013년 9월 23일
편집: Jai 2015년 10월 28일
I'm getting the following error when using chol
Error using chol
Matrix must be positive definite.
My matrix is a squared exponential kernel matrix and by definition should always be positive definite. It's defined as
K(i,j) = exp(-(1/L)*(x_i-x_j)^2) + beta*delta(i,j)
The problem occurs during an attempt to optimise over the parameter L. When the value of L is very small (approx 1e-6) I end up with zero off diagonal elements. I have saved the matrix to a text file before calling chol and the problematic matrix looks like
10.25310436455727 0 0
0 10.25310436455727 0
0 0 10.25310436455727
Using chol on this matrix obviously doesn't lead to an error. Does anyone know why chol is giving me this error? Could there be rounding errors due to taking the exponential of very small numbers that aren't being saved to the text file?
Kind Regards, James
  댓글 수: 1
Matt J
Matt J 2013년 9월 24일
편집: Matt J 2013년 9월 24일
I have saved the matrix to a text file before calling chol and the problematic matrix looks like
Like Alan, I suspect you are not saving the matrix that is actually causing the error, even if you think you are. You should use
>> dbstop if error
to trap the error, if you aren't doing so already. You should also save the matrix K that triggers dbstop to a .mat file (not a text file) and attach it here, so that we can examine K with the full numerical precision that CHOL sees.

댓글을 달려면 로그인하십시오.

답변 (3개)

Alan Weiss
Alan Weiss 2013년 9월 24일
You say you are optimizing. Are you calling fmincon? It is possible that your optimization is stepping out of bounds, taking L to be negative or zero, during a finite differencing step.
If this is the case, perhaps you should set a lower bound on L that is a bit greater than zero, and use the sqp or interior-point algorithms because these algorithms always respect bounds.
If you are using a different optimizer, then perhaps what is happening is when L becomes small, the exponential terms become zero due to underflow, and the matrix is no longer positive definite.
Just my random thoughts.
Alan Weiss
MATLAB mathematical toolbox documentation
  댓글 수: 2
Matt J
Matt J 2013년 9월 24일
편집: Matt J 2013년 9월 24일
If this is the case, perhaps you should set a lower bound on L that is a bit greater than zero, and use the sqp or interior-point algorithms because these algorithms always respect bounds.
Or, perhaps parametrize K instead as
K(i,j) = exp( z *(x_i-x_j)^2) + beta*delta(i,j)
and optimize with respect to z. Then, later, after the optimal z is found, retransform to obtain L=-1/z. This way, K is defined for all parameter values.
James Barrett
James Barrett 2013년 9월 24일
Thanks both for your comments. I have in fact imposed a lower bound on L by writing it as
L = log(1 + LOWER_BOUND + exp(z))
and optimising over z (this is a nice way to impose constraints without using the optimiser to do so).
Alan, the matrix should be positive definite even if the exponential terms go to zero because of the diagonal beta elements (which are non negative).

댓글을 달려면 로그인하십시오.


Sean de Wolski
Sean de Wolski 2013년 9월 24일
Hi James,
When I started at MathWorks, I asked our MATLAB/Math development team a very similar question. It was the first time I was ever properly scolded by Cleve :)
Here's what's happening: from the documentation for chol
R = chol(A) produces an upper triangular matrix R from the diagonal and upper triangle of matrix A, satisfying the equation R'*R=A. The chol function assumes that A is (complex Hermitian) symmetric. If it is not, chol uses the (complex conjugate) transpose of the upper triangle as the lower triangle. Matrix A must be positive definite.
So A becomes
Ah = triu(A)+triu(A,1)'
If you look at:
eigs(Ah)
There will be a negative one.
  댓글 수: 2
James Barrett
James Barrett 2013년 9월 24일
Thanks Sean, shouldn't the matrix be symmetric (unless there are some rounding errors)? I have even done
K = (K+K')/2
before I call chol(K).
Secondly, it's not clear to me why there would any negative eigenvalues. I've constrained beta to be always greater that 1e-4 to ensure the matrix will always be positive definite.
Jai
Jai 2015년 10월 28일
편집: Jai 2015년 10월 28일
Hi Sean, It seems you are correct, indeed eigs(Ah) gives negative eigen value. Is there any workaround to make it work using chol ?

댓글을 달려면 로그인하십시오.


Leah
Leah 2013년 9월 24일
I used Spectral decomposition when this happens. Which is just a fancy name for setting the negative eigenvalues to zero before chol.

카테고리

Help CenterFile Exchange에서 Operating on Diagonal Matrices에 대해 자세히 알아보기

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by