Closed bsxfan closed 9 years ago
There is a bug in cholfact:
julia> logdet(rand(3,3))
ERROR: PosDefException not defined
in cholfact! at linalg/factorization.jl:13
in logdet at linalg/dense.jl:422
Also, is there or can there be some interface for getting log(abs(det)) and the sign? Or det of non-positive-definite matrices?
@JeffBezanson A det method for the LU type already exists. So logdet for a square matrix should be defined similarly although it would need to return log(abs(det)) and the sign separately.
det()
and LU
are already equipped to return a complex result for a complex matrix. So maybe logdet(::LU)
could always return a complex value too, which would accomodate the case of both negative and complex determinants.
The inconsistency is that at present logdet(::Matrix)
returns a real result. Maybe one should rather define complex_logdet()
which always returns a complex result.
Having logdet always return a complex result occurred to me too and feels more appealing than returning both a value and a sign, bu it does seem inconsistent with our other functions, which return NaN when applied to a real value in such a way that the result would be non-real.
Cc: @dmbates @andreasnoack
We now throw domain errors instead of returning NaN.
Oh, right. But same idea.
The thinking among the R developers was that returning the logarithm of the absolute value of the determinant along with its sign was done to provide a more stable calculation and to avoid overflow/underflow. However I will admit that almost all the time I use logdet it is for a positive definite matrix, usually some form of a variance-covariance matrix.
While we're on logdet, it might be nice to supply Woodbury with a logdet too. Woodbury is e.g. useful in probability theory when doing factor analysis, in which case the determinants also end up being positive. Since the woodbury matrices can be huge, the determinant may well under/overflow, but logdet won´t.
Restricting logdet
to Hermitian positive-definite matrices seems wrong to me.
For example, in some problems arising in integral-equation methods for quantum field theory, we need the logdet of a real matrix A
that is nonsymmetric but its symmetric part is positive definite, i.e. A + A'
is symmetric positive-definite. Such a matrix always has a positive determinant and so its logdet
is real.
Also, currently the logdet
of a Hermitian positive-definite matrix returns a complex result (with zero imaginary part), which doesn't seem right.
I agree with @stevengj that returning a complex result from logdet
of a Hermitian positive-definite matrix doesn't seem right.
So returning a complex number from logdet
doesn't seem right and restricting logdet
to positive-definite matrices doesn't seem right either. I can see only two other options:
logdet
.logdet
to some larger set of matrix inputs – possibly any matrix for which the determinant is positive?Are there other options I'm missing? If not, which of these is the better choice?
It could work like eig
and return a complex result when the determinant is negative and real otherwise.
I'm not a fan of eig
behaving that way, honestly.
Then I'd prefer DomainError
over tuple output.
Another option is to have two functions, of which one could be called logabsdet
for example.
@StefanKarpinski, I would tend to (2) restrict logdet
to matrices with positive real determinants, in which case it always returns a real result, and throw a DomainError
otherwise. To begin with, there are three main cases:
In theory, I suppose that you could also have an unsymmetrical complex matrix that nevertheless has a positive real determinant, but I'm not sure that this case arises enough that it is worth the trouble to figure out how to compute logdet
efficiently for it.
How about having logdet
as proposed by @stevengj and logdetabs
, a slight modification of the name proposed by @timholy, to make the name more discoverable?
@ViralBShah, @timholy Tim proposed logabsdet
not logdetabs
. I think the logabsdet
name more clearly reflects the operation being performed.
I think Viral was acknowledging that, but logdetabs
is tab-completable from logdet
.
Another option would be a keyword argument indicating whether to take the absolute value or not.
Of course - let's go with a keyword argument.
Could someone enlighten me on the advantage of logdet(A, abs=true)
over logdetabs(A)
?
I was proposing logdet(A; abs=true)
. The real question is whether the two routines are similar enough that they should share a name.
I think @stevengj 's point was that logdetabs
is shorter to type than the keyword/optional input version.
But if we do give it a separate name, I think logabsdet
is really much better than logdetabs
, despite the tab-completion issue. logdetabs(A) = log(det(abs(A)))
, but logabsdet(A) = log(abs(det(A)))
. These are very different, and it would be dangerous to confuse them.
Yeah, I'm not sold on using a keyword, just an idea. Calling it logdetabs seems outright obtuse.
It sounds like logabsdet
and logdet
return different things. Is it good for keyword arguments to change the function output type?
Yeah, good point. I was imagining the keyword as only affecting whether a negative determinant causes an error or a sign flip (I realize it may not be implemented that way).
Ok. The arguments for logabsdet are quite persuasive.
See my comment in #3295, where I give a tyre-kicking implementation of a generalized, but type-stable logdet for the LU factorization. There is a logdet2
which returns log(abs(det)) and sign(det) as well as a logdet
which returns complex log for complex arguments and real log for real arguments. The latter gives Domainerror when the det<0, but works for non-positive definite matrices with positive determinants, e.g. diagm(-ones(2))
.
The current logdet implementation in linalg/dense.jl,
would probably be cleaner if implemented as: