issues with Cholesky decomposition
19 Ansichten (letzte 30 Tage)
Ältere Kommentare anzeigen
Hello guys. I have a question with Cholesky decompostion. I am running a code for MCMC sampling in high dimensional linear regression. You know that this means I have to perform Cholesky decomposition many times. In some cases I get error by non-positive definite covariance matrix. However, when I rerun the code again with same dataset, the problem can disappear. I really have no clue why this happens. Any suggestions? Thanks!
1 Kommentar
Walter Roberson
am 28 Feb. 2020
It is possible to get a small bit of round-off error in building your matrix. For any given matrix, M, the way around it is:
M = (M+M.')/2;
Antworten (1)
Christine Tobler
am 2 Mär. 2020
This can happen if your matrix is close to symmetric positive semi-definite (meaning the smallest eigenvalue is around machine epsilon compared to the largest eigenvalue). The Cholesky function does not support symmetric positive semi-definite input.
If this doesn't happen every time you run the code, it must mean that the matrices you're passing to Cholesky are not exactly equal between runs - one of the other operations used must give slightly different results between runs. The function chol guarantees that if you call it twice with the exact same input within the same MATLAB session, it will return the exact same output (unless you're changing some very specific MATLAB settings in between).
How are you generating your covariance matrix? If it's by forming C = M'*M, you could instead compute the QR decomposition of M: M = Q*R, M'*M = R'*Q'*Q*R = R'*R (using that Q'*Q is the identity matrix for the QR decomposition). This is numerically more robust; note that if C is close to positive semi-definite, the last row(s) of R will likely be close to machine epsilon relative to the other entries of R. You may need to be careful about using R depending on the next step of your algorithm.
Siehe auch
Kategorien
Mehr zu Matrix Decomposition finden Sie in Help Center und File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!