chol() returns error even the matrix is positive definite.

7 次查看（过去 30 天）
Hung Dao
Hung Dao 2022-12-6
编辑： Bruno Luong 2022-12-9
I have a matrix , where I is a identity matrix and B is a matrix. It is clear that M is a positive definite matrix. However, I am facing a situation in which the Matlab Cholesky decomposition function
chol(M);
returns an error stating that M must be positive definite. I examine matrix B and realise that it's elements are pretty small, resulting in the elements in matrix are extreme large in magnitude (the magnitudes are about e+15 to e+17).
Any advice/suggestion how to handle this issue is much appreciated !
Update on 9 Dec 22:
To help explaning the issue more cleary, I attached the data containing the matrix. Please load the data 'example.mat' and run the code below
[p, r] = size(B);
B1 = d.\B;
z1 = z./d;
z2 = B1'*z1;
R = chol(eye(r) + B1'*B1);
  2 个评论
Walter Roberson
Walter Roberson 2022-12-6
Taking B'*B cannot result in large values if the entries are all small... not unless there are a huge number of entries. Each element of B'*B is the dot product between a row and another row of B. The maximum possible value would be max(abs(B))^2 * number of columns
Perhaps you are looking at the pseudoinverse of the matrix instead of the transpose
Hung Dao
Hung Dao 2022-12-9
Thank you for your comment. My apologies for my late reply as I am travelling. You are right. Matrix B is quite a big matrix (it has the size of ). Many of the elements are small (23,036 entries are smaller than 1e-1) but I do realise that B doeve has some very large elements (26 elements are greater than 1e+6 with the maximum is about 1e+9).
I attached the data and the code. Hopefully it helps explain my issue better.
Thank you for answering my question. I appreciate it !

请先登录，再进行评论。

请先登录，再回答此问题。

回答（2 个）

Bruno Luong
Bruno Luong 2022-12-9
编辑：Bruno Luong 2022-12-9
This statement
B1 = d.\B;
looks very suspecious to me (the .\ operator), if you remove the dot it will work.
Here you do element-wise of d that oscillate to 0, create a huge number in B1 (it can goes to infinity indeed s such noisy data)
[p, r] = size(B);
B1 = d\B;
R = chol(eye(r) + B1'*B1);
Now wht you trying to do and implement only you know, but I suspect your code is incorrect for whatever the goal you want to achieve.
  2 个评论
Hung Dao
Hung Dao 2022-12-9
Thank you for your answer. I appreciate your support. However, I respectfully disagree with you.
I am using
B1 = d.\B
as I want to calculate matrix , where , a diagonal matrix with the main diagonal is vector d. Hence, replacing my formula using the .\ operator by
B1 = d\B
is incorrect.
Bruno Luong
Bruno Luong 2022-12-9
编辑：Bruno Luong 2022-12-9
Only you know what is correct, but your d looks like this
It looks like it contains noise close to 0 and it cross randomly 0 mny time. 1./ d.\1 is then randomy big.
It that make sense for you then OK, but don't expect getting anything meaningful from that.
B1'*B1 have condition number and large element of the order of probably 1e20. Add eye(r) to it is meaning less, doing cholesky doesn't make any sense numerically. Your code won't produce anything usuful but gardebage with such data.

请先登录，再进行评论。

Askic V
Askic V 2022-12-6
编辑：Walter Roberson 2022-12-7
Ran in:
Without actual code it is very hard to know for sure, but using chol() function in Matlab is actually a preferred way to test if the matrix is positive definite.
For me, it is very hard to understand how values of B'*B can be extremly large if values of B are small. Something is wrong. Please have a look at the example:
B = randn(3,3)*1e-6;
C = B'*B
C = 3×3
1.0e-11 * 0.3916 0.0891 -0.0314 0.0891 0.0303 -0.0417 -0.0314 -0.0417 0.1399
max(C(:))
ans = 3.9156e-12
min(C(:))
ans = -4.1702e-13
Regarding your problem, I suggest you to read these advices:
https://www.mathworks.com/matlabcentral/answers/225177-chol-gives-error-for-a-barely-positive-definite-matrix?s_tid=answers_rc1-1_p1_MLT
  4 个评论
显示 2更早的评论
Bruno Luong
Bruno Luong 2022-12-9
编辑：Bruno Luong 2022-12-9
@Askic V "In this particular example, you do have symmetric but not positive definite matrix"
I claim is wrong
C = eye(r) + B1'*B1
is always spd theorically since
x'*C*x is equal to |x|^2+|B1*x|^ 2 >= |x|^2 > 0 for any vector x ~= 0. Which is the requirement of spd (symmetric is easy to see).
Even B1'*B1 is non negative, eig fails to estimate the smallest eigen value.
OP simply has numerical error due to bad scaling issue.

请先登录，再进行评论。

类别

在 Help Center 和 File Exchange 中查找有关 Linear Algebra 的更多信息

标签

Community Treasure Hunt

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

Start Hunting!

Translated by