log(Dual(0.0,0.0)) is Dual(-Inf,NaN)

[tpapp]

I have been tracking an issue of NaN propagation in my code, and narrowed it down to the following MWE:

julia>  log(Dual(0.0,0.0))
Dual(-Inf,NaN)

Would Dual(-Inf,0.0) make sense? After all, the value is not sensitive to that variable.

The context is a Hidden Markov Chain problem for which I am fitting a model using derivative-based Bayesian methods. It makes sense to keep track of the log likelihood because of underflow/overflow.

A quick workaround (that would overwrite the method) would be appreciated so that I could proceed with work (I haven't had time to dig into ForwardDiff internals yet), but given some advice on what to do I would also submit a PR.

[KristofferC]

Dual(-Inf,0.0) does not make sense in my opinion. The derivative of log at the point 0.0 is not 0.0 which this would give. Perhaps Dual(-Inf, Inf).

[jrevels]

Would Dual(-Inf,0.0) make sense? After all, the value is not sensitive to that variable.

Hmm, this does make sense to me. Dual(0.0, 0.0) should behave as close to 0.0 as possible since the derivative isn't actually being tracked. Throwing in NaNs in this case just "poisons the well" unnecessarily. I'll do some thinking about how to implement this with as little performance loss as possible.

The derivative of log at the point 0.0 is not 0.0 which this would give. Perhaps Dual(-Inf, Inf).

Note that this is already the current behavior when the perturbation is actually present, and it "correctly" picks negative infinity or positive infinity depending on the sign bit ("correctly" in quotes because it really should be complex infinity, but this pick-a-side convention is what we already use for abs):

julia> log(Dual(0.0, 1.0))
Dual(-Inf,Inf)

julia> log(Dual(-0.0, 1.0))
Dual(-Inf,-Inf)

[tpapp]

To give a bit of context: in my particular problem (hidden Markov model likelihood calculation) all I am using is sums and products. The objects are probabilities, some of which are constant 0.0. The calculation would go through this way, and these objects would just give a Dual(0.0,0.0).

However, because of underflow, we work with logarithms, and 0.0 becomes -Inf. Products are now sums, and sums are logsumexps.

Intuitively, I think that if a partial is 0.0 in some x, it should remain 0.0 for f(x). But I am not experienced with AD so I don't know how people handle this.

@jrevels: do you think that overwriting log(::Dual) to special-case this would be reasonable workaround for me until the issue is fixed?

[jrevels]

@jrevels: do you think that overwriting log(::Dual) to special-case this would be reasonable workaround for me until the issue is fixed?

Yeah, you should be able to just copy-paste this definition:

function Base.log(n::ForwardDiff.Dual)
    x = log(ForwardDiff.value(n))
    d = 1 / ForwardDiff.value(n)
    ps = ifelse(ForwardDiff.isconstant(n), ForwardDiff.partials(n), d * ForwardDiff.partials(n))
    return ForwardDiff.Dual(x, ps)
end

EDIT: Left out the derivative, oops!

Intuitively, I think that if a partial is 0.0 in some x, it should remain 0.0 for f(x).

This is basically the ideal situation, but some weirdness can pop up if you start propagating NaNs and Infs (as we've seen here).

[tpapp]

@jrevels : thank you so much! now it works. after you have thought about this, please consider adding this to the library (or I can make a PR, with tests).

[mlubin]

I don't claim to know what the correct solution is, but we have to work around what I believe is this same issue in JuMP:
https://github.com/mlubin/ReverseDiffSparse.jl/blob/72b2a3002a75c3668b53336cb4d2988d9e07eb26/src/forward.jl#L239
https://github.com/mlubin/ReverseDiffSparse.jl/blob/72b2a3002a75c3668b53336cb4d2988d9e07eb26/src/reverse.jl#L27