On 20 October 2015 at 16:44, Laurent <[email protected]> wrote:
> Thanks Paolo, yes, I saw that indeed and it is most likely relevant. I was
> hoping my problem would be sufficiently more specific; in particular, I'd
> like to detect whether the value being above 0 is likely a rounding error
> (in which case I can just round it down to 0 again) or likely a mistake on
> my part in my equations or in their implementation.

I'm sorry, but your dichotomy seems a fallacy. Rounding errors are the
rule in floating point, and you need to account for them, otherwise
you have a bug, as every introduction to numerical analysis should
teach.

For instance, with infinite precision reals, (1 - x - y - z) = 1 - (x
+ y + z), but you can't use that in floating point. I'll focus on the
case of log-space probabilities: (lg- (lg- (lg- (fllog 1.0) x) y) z)
and (lg- (fllog 1.0) (lgsum x y z)) can have different behavior.
Log probabilities are good for probabilities so small they'd
underflow, whose logs are never close to 0.
But (lg- (fllog 1.0) x) *is* close to zero, where lg+ and lg- behave
badly. If you instead first compute (lgsum x y z), you get a bigger
log-space probability, so you're somewhat safer; plus (lgsum x y z) is
better than using lg+.

If x is in fact close to 0 itself, the reverse happens, but there's
high danger that x has little precision, so you should try hard to
restructure your computation to avoid such intermediate results.

In general, without an expert inspecting the code, and assuming your
equations themselves (and their implementation) are bug-free on
infinite precision reals, they're still likely to be buggy on floating
point numbers, just because it's floating point.
if you are not studying error propagation in your code, you might well
get nonsense in tons of other ways.

I understand experts would actually compute the error bounds on the
result based on the values and errors in input, which seems very
time-consuming. Luckily, Racket's library has support for studying
error propagation without pen and paper:
http://docs.racket-lang.org/math/flonum.html#%28part._.Debugging_.Flonum_.Functions%29

I'll admit that some mortals can get this right, and some seem to hang
around https://math.stackexchange.com/questions/tagged/numerical-methods,
but I don't think I'm one of them, so I stay away from floating point
if I can :-).

Cheers,
Paolo

> On Tue, Oct 20, 2015 at 2:24 PM, Paolo Giarrusso <[email protected]>
> wrote:
>>
>> On Tuesday, October 20, 2015 at 12:15:21 PM UTC+2, Laurent Orseau wrote:
>> > The built-in log-space arithmetic operations are wonderful. (Thanks so
>> > much Neil!)
>> >
>> >
>> > It sometimes happens that after some log-operations where the result is
>> > very close to 0 that it actually goes slightly above 0 as a result of
>> > floating point errors, and thus is not a log-probability anymore, which 
>> > gave
>> > me a few headaches.
>> >
>> >
>> >
>> > One obvious approach is to surround the result with a check, and if it
>> > is positive but not, say, above 1e-8 then round it down to 0, otherwise
>> > raise an exception if positive.
>> >
>> >
>> > But this feels like a hack and I was wondering if there was any "good"
>> > way to ensure that any correct log-space operation remains a 
>> > log-probability
>> > (i.e., never goes above 0).
>> >
>> >
>> > Thanks,
>> > Laurent
>>
>> If lg+ is involved, this discussion in the docs seems relevant — apologies
>> if you've seen that already:
>>
>> http://docs.racket-lang.org/math/flonum.html#%28def._%28%28lib._math%2Fflonum..rkt%29._lg-%29%29
>> ?
>>
>> After explaining that a good solution is an open research problem, docs
>> state:
>> >  Further, catastrophic cancellation is unavoidable when logx and logy
>> > themselves have error, which is by far the common case. [...] There are
>> > currently no reasonably fast algorithms for computing lg+ and lg- with low
>> > relative error. For now, if you need that kind of accuracy, use
>> > math/bigfloat.
>>
>> Disclaimer: not an expert here, might well be missing something.
>>
>> --
>> You received this message because you are subscribed to the Google Groups
>> "Racket Users" group.
>> To unsubscribe from this group and stop receiving emails from it, send an
>> email to [email protected].
>> For more options, visit https://groups.google.com/d/optout.
>
>



-- 
Paolo G. Giarrusso - Ph.D. Student, Tübingen University
http://ps.informatik.uni-tuebingen.de/team/giarrusso/

-- 
You received this message because you are subscribed to the Google Groups 
"Racket Users" group.
To unsubscribe from this group and stop receiving emails from it, send an email 
to [email protected].
For more options, visit https://groups.google.com/d/optout.

Reply via email to