issues with Cholesky decomposition
19 次查看(过去 30 天)
显示 更早的评论
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 个评论
Walter Roberson
2020-2-28
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;
回答(1 个)
Christine Tobler
2020-3-2
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.
另请参阅
类别
在 Help Center 和 File Exchange 中查找有关 Matrix Decomposition 的更多信息
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!