product of nonsingular matrices is singular?
1 回表示 (過去 30 日間)
古いコメントを表示
I am generating this 1000x1000 real matrix A that is nonsingular. Linear algebra claims that det(A) = det(A') and that det(AA') = det(A'A) = det(A)det(A'). So if A is nonsingular, so are A' and A'A. Therefore if rank(A) = 1000 (full rank), then rank(A') = 1000 and rank(A'A) = 1000. Yet:
rank(A)
ans =
1000
>> rank(A')
ans =
1000
>> rank(A' * A)
ans =
927
Of course I want it to be rank(A'*A) = 1000. What can I do?
0 件のコメント
回答 (2 件)
Matt J
2017 年 12 月 2 日
編集済み: Matt J
2017 年 12 月 2 日
Therefore if rank(A) = 1000 (full rank), then rank(A') = 1000 and rank(A'A) = 1000
But just how strongly non-singular is A? Have you checked with cond(A)?
The linear algebra theory that you are relying on does not take into account the floating point computer errors involved in calculating A'*A nor the floating point errors required to evaluate the rank, nor in fact the inexact computations that produced A itself, if there were any. This is why the rank command takes a tolerance argument rank(A,tol). The default tolerance, tol, that you are using may not be appropriate for your matrix.
0 件のコメント
John D'Errico
2017 年 12 月 2 日
編集済み: John D'Errico
2017 年 12 月 2 日
Trivial.
However, NEVER. NEVER. NEVER use det to check if a matrix is singular. There, I said it three times so it must be true.
So what can you do? DON'T USE DET!!!!!!!!!!!! The simple example of why not is...
Is the matrix eye(1000) singular? DUH. Of course not. It is as non-singular as a matrix can be. What does det tell us?
A = eye(1000);
det(A)
ans =
1
WOW. DET even got that one right.
But how about A/3? Surely if A is as non-singular as possible, multiplying it by a scalar non-zero constant will not change the singularity.
det(A/3)
ans =
0
det(A*3)
ans =
Inf
Hmm. Just be scaling A by some small constant in either direction, I can get ANY number I want out.
The point is, DET is the worst thing you can use here. Yes, I know, your teacher told you to use the determinant to check for a singular matrix. But guess what, they were wrong!
The problem is, people learn from their peers, their teachers, their mentors. But if your various mentors never learned the truth, they will propagate falsehoods. ANd then so will you. So learn this now:
DO NOT USE DET. At least don't use it unless you know enough about the question to know why and when you should not use DET! (Yes, there are some circumstances where use of det can be reasonable. Again, unless you understand DET well enough to know what they are, then don't use DET.)
So now to answer your real question. Can the product of two square matrices have lesser rank? Consider this simple example:
A = vander(1:10);
rank(A)
ans =
10
cond(A)
ans =
2.1063e+12
A is not singular. But it is not that far away from singularity.
What does det say?
det(A)
ans =
-1.8349e+21
Well, det is a bit confused here, as we showed above. But, what happens if we just form AA= A*A?
AA = A*A;
rank(AA)
ans =
9
cond(AA)
ans =
1.1464e+18
Rank thinks it is singular. So does cond.
det(AA)
ans =
3.3373e+42
DET is still confused.
So AA is NUMERICALLY singular, despite the fact that AA is not truly singular. We know that to be true.
But remember that all work here was done using double precision arithmetic.
In the end, you asked what you can do. Well, don't use det. Use rank, or cond, or svd. Learn enough about double precision arithmetic to know why those tools are reasonable.
0 件のコメント
参考
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!