Dett

Computes the determinant of non-square matrices.

Description

Computes the determinant of a matrix of any size using the QR decomposition. Version 2, Nov 2006.
See my updated Adj also.

• MATLAB release: MATLAB 6.5 (R13)

[Jos x@y.z]

This is mathematical nonsense.
For non-square matrices, there are always rows/columns that are not indepenent of other rows/columns. So, by definition, the determinant of non-square matrices should be zero!
That is why textbooks always say that the determinant is only defined for square matrices.

[Paul Godfrey]

Not nonsense at all.

Consider that for non-square matrices, a left or right inverse may exist.
The inverse may be defined as Adjoint(A)/Det(A)
This program computes a meaningful value for Det(A) consistant with that definition.
Taker a look at the Adj.m program to put all the pieces together.

[Paul Godfrey]

Regarding the comment of Jos x@y.z
"For non-square matrices, there are always rows/columns that are not indepenent of other rows/columns."
Not true at all. Consider orth(rand(3,2))
You have 2 orthogonal vectors in a non-square matrix.

[John D'Errico]

A right or left pseudo-inverse MAY exist. But this "determinant" seems to tell you nothing that rank or cond will not say better. (NEVER use any variety of determinant to identify the singularity status of a matrix.)

If your goal is merely to compute a pseudo-inverse, then use pinv - a far faster solution than one based on the adjoint.

Looking more closely at the author's claim, I'll note that adj fails when applied to non-square matrices. I just download the latest version of adj. Try this:

adj(rand(5,3))
??? Error using ==> eig
Matrix must be square.

Error in ==> adj at 16
ce = poly(eig(A));

The other way fails too:

adj(rand(3,5))
??? Error using ==> eig
Matrix must be square.

Error in ==> adj at 16
ce = poly(eig(A));

This does not even touch the question of complex (particularly self adjoint) matrices as brought up by a review of adj. I've not checked that issue.

For square matrices, this code is also far slower than Matlab's existing det function anyway. So at the very least users should be cautioned away from using this code on square systems. Since apparently non-square systems are not handled as the author claims, I see no use at all, unless you wish to compute the pseudo-inverse of a scalar.

Just use pinv for your pseudo-inverses. Pinv will work for rank deficient problems too, whereas this code will fail due to the singularity.

[Jos x@y.z]

1) The adjugate is also only defined for square matrices.
2) orth(rand(3,2)) produces three rows that are NOT linearly independent, e.g., its rank is 2.

[Paul Godfrey]

Regarding some comments by Jos x@y.z:
1) The adjugate is also only defined for square matrices.
The real question to ask is "Can an Adjugate matrix be defined in a meaningful manner for non-square matrices?"
For many non-square matrices the Matlab Pinv command can generate a left, right, or true inverse for the matrix.
The matrix inverse is typically defined as Adj(A)/det(A).
So the problem then becomes one of finding a useful value for Adj(A) and det(A) when A is non-square that "matches" in Pinv(A) in some sense.
My Adj(A) and dett(A) programs do exactly that.

2) orth(rand(3,2)) produces three rows that are NOT linearly independent, e.g., its rank is 2.
And that tells me that a left inverse exists but not a right inverse.
The left inverse should generare Eye(2)
and the right inverse would (if it existed) generate Eye(3).
Therefore, that inverse should be Adj(A)/dett(A).
Or, in terms of dett(A), Adj(A)*A = I*dett(A)
To see more clearly what the value of dett(A) should be try this:
A=rand(3,2)
B=Adj(A)
B*A
eye(2)*dett(A)

They are equal, as expected.

I think that the real problem here is that most Linear ALgebra books don't completely cover the Adj() function.
The reason for that is, probably, that the Adj() was computed using co-factors and such which is horribly inefficient.
I don't think many people looked at computing the Adj() using the SVD.
Try pinv(A), Adj(A)/dett(A)

[Paul Godfrey]

Regarding the comments of John D'Errico:
A right or left pseudo-inverse MAY exist. But this "determinant" seems to tell you nothing that rank or cond will not say better.
<I agree. But the purpose of my Adj() and Dett() programs is to provide meaningful values for non-square matrices.>

If your goal is merely to compute a pseudo-inverse, then use pinv - a far faster solution than one based on the adjoint.
<No, my goal is to extend the definition in a meaningfull manner>

Looking more closely at the author's claim, I'll note that adj fails when applied to non-square matrices. I just download the latest version of adj.
<Sigh, you don't have the latest version.
Is your version dated April 1998 or October 2006?>

Try this:
adj(rand(5,3))
<Works fine with the latest version.
Give the Mathworks people a day or to to make the latest version available.>

This does not even touch the question of complex (particularly self adjoint) matrices as brought up by a review of adj.
<It works fine for complex matrices as well.
There are multiple definitions of Adj().
My definition is based on Inv(A) = Adj(A)/Det(A)
as given on page 169 of "Applied Linear Algebra, 3rd edition" by B. Noble and J. Daniel

For square matrices, this code is also far slower than Matlab's existing det function anyway.
<True. And is also true the Matlab's own Det() function doesn't handle non-square matrices like mine does>

So at the very least users should be cautioned away from using this code on square systems.
<True>

Since apparently non-square systems are not handled as the author claims,
<Try again in a day or two.>

I see no use at all,
<Other than providing a meaningful extension of the definition of Adj() and Det()>

[John D'Errico]

Cogent arguments by Paul in favor of dett. If one has the goal of an extension of the determinant in a way that is consistent with the pseudoinverse, then adj(A)/dett(A) succeeds, WHEN A is nonsingular. Obviously this fails for singular arguments.

Were there other uses proposed for this extension of dett, I might up my rating further. For example, det is often used to compute the volume of an n-d simplex. Is there a similar extension for dett? (Consider a simplex that lives in a planar subset of a higher dimensional space. Dett should indeed apply here.) Another use for the determinant is as applied to the Jacobian in variable transformations. Is dett consistent as an extension there?

The new version of adj is up there now, so I could test them both further. I was wondering if LU would offer a more efficient approach to its computation. This would work if a standard property of the determinant also applied:

A = triu(rand(3,4));
dett(A) ~= prod(diag(A))

An LU based dett would allow it to apply to sparse matrices - which might offer an improvement over pinv.

So, along those lines, why not consider a QR version of dett? QR applies to sparse systems. It is numerically well behaved, especially if the column pivoted form is used. It will be faster than the svd solution. The only flaw is a QR based dett will work only for non-square systems with more rows than columns. So when there are more columns than rows, apply the QR to the transpose.

[n,m] = size(A);
if n<m
  A = A';
end
[Q,R] = qr(A);
D = det(Q)*prod(diag(R));

This can be hacked easily enough to use a column pivoted QR, though the column pivoted version does not work for sparse problems. You would need to revert to the simpler code above when sparse.

[Paul Godfrey]

Regarding the comments of John D'Errico:
"...WHEN A is nonsingular. Obviously this fails for singular arguments."
And that "failure" is the correct answer for a singular matrix.
For any matrix,(square or non-square), with one or more singular values equal to 0, no inverse (left, unique, or right) is going to exist.
In that case we resort to the SVD to provide the "best" answer under those circumstances.

"Were there other uses proposed for this extension of dett, I might up my rating further..."
Well, some of us consider a consistant extended definition of Det() to be "use" enough... ;)

[John D'Errico]

So why should anyone bother to use your code, when they can copy a better version out of this review? Can you show ANY reason why the svd version is clearly the best answer as you claim?

function [d] = dett2(A)
%Dett2 Determinant of any sized matrix
%
%usage: d = dett2(A)
%
%see also: Det, Cond, Eig, SVD

[r,c]=size(A);
if r>=c
  [Q,R] = qr(A);
else
  [Q,R] = qr(A');
end
if r==1 | c==1
   R=R(1);
else
   R=prod(diag(R));
end
d=det(Q)*R;

Comparison:

>> A = rand(500,200);
>> tic,p=dett(A);toc
Elapsed time is 1.351050 seconds.
>> tic,p=dett2(A);toc
Elapsed time is 0.688662 seconds.

So, dett, based on the svd is slower than it needs to be.
Is it more capable? Dett2, as I've just supplied, works
nicely on sparse matrices. Can that be said for dett?

>> A = sprand(500,200,.1);
>> tic,p=dett2(A);toc
Elapsed time is 3.930455 seconds.

>> tic,p=dett(A);toc
??? Error using ==> svd
Use svds for sparse singular values and vectors.

Error in ==> dett at 20
[u,s,v]=svd(A);

So, no, I don't consider a poor implementation of any tool to be worth an excellent rating. Not when a superior one is so trivial to write.

[Paul Godfrey]

Regarding the comments of John D'Errico:
"So why should anyone bother to use your code, when they can copy a better version out of this review?"
Because it's publicly available as a submitted .m file.
Why don't you go ahead and submit your own version?
Won't that make it easy for all Matlab users to find and use?

"Can you show ANY reason why the svd version is clearly the best answer as you claim?"
I never claimed that the SVD algorithm was the best one to use for computing DETT.
I claimed that the SVD algorithm is the best one to use when "solving" a singular set of equations.
It is common knowledge that the SVD gives you the Least squares solution for over determined systems, the exact answer for full rank systems, and the minimum norm answer for under determined systems.

Why did you take a discussion about solutions of linear systems and twist it into a claim about the best algorithm to use?

IMHO, you appear to want to try and win some kind of argument or contest of some sort.
Your postings have been in a rather belligerent and "in your face" style.
I prefer a more thought based rational discussion instead.

"dett, based on the svd is slower than it needs to be"
Ok, good for you.

"nicely on sparse matrices. Can that be said for dett?"
Well, I never claimed that it was designed for sparse matrices.

"Error in ==> dett at 20"
So, go ahead and submit it.

"So, no, I don't consider a poor implementation of any tool to be worth an excellent rating"
Ok, but I never asked you to rate it.
Why do you think that I care?

"Not when a superior one is so trivial to write."
Sigh.
You really want to win something here don't you?
I can't help but ask:
Before you saw this submission, how did you calculate the Determinant and Adjoint of non-square matrices?

[Comment Responsibly]

Mr. Godfrey I thank you. There are a few bullies running around on the Matlab exchange randomly attacking authors. Most authors just ignore them, but I fear some may be turned away from the site because of them. So thank you for sticking up for yourself and demonstrating how shrill, petty, and useless their comments usually are. These people seem to want published, polished, ready-for-sale code. The truth is matlab isn't for hardcore programmers. It is for people who want a TOOL to solve a particular PROBLEM THEY need solved, not a product development suite. So most programs found on the file exchange will be imperfect and have a very limited scope.

I'm amazed by the sense of entitlement these commenters display. Contributors have, out of the goodness of their hearts, given free code to the world. It is FREE. If they don't like it, they need not use it.

Thanks again Mr. Godfrey.
I will give Mr. D'Errico credit for engaging you in the discussion, regardless of the tone. I am also very impressed that he changed is rating after some discussion.
Jos x@y.z however seems to have a habit of being the first and least knowledgable commenter on many programs. He also rarely returns to apologize or retract his statements when he is shown to be very wrong. His rediculous snap judgments are exactly the type of thing my whole comment is about.

[Tim Davis]

The determinant of a matrix can easily be extremely large or extremely small for a modest-sized matrix, leading to numeric underflow or overflow. I would suggest an option that sums the log of the diagonal of R, instead of prod(diag(R)). I do the same trick in umfpack with [d e] = umfpack(A,'det') which returns the determinant as d*10^e.

The function det in MATLAB has the same underflow/overflow problem, but det(Q) I think is +1 or -1 so perhaps it's not a issue.

Seems like there should be a better way of computing det(Q) than with det(Q), since Q is an orthogonal matrix and has such nice properties.

Just a suggestion.

[Ben Petschel]

Interesting idea, it would be good to see references to the literature; for example does it relate to the one in
 http://www.emis.de/journals/BAG/vol.46/no.2/b46h2rad.pdf

which can be used to calculate the areas of polygons, for example.

Any definition for rectangular matrices, including the above, would be a major rewrite of the theory for square determinants which are defined as alternating multilinear forms; this would imply that det(A)=0 whenever m<n.