# Fast logsumexp

**URL:** <https://discourse.julialang.org/t/fast-logsumexp/22827>\
**Category:** Performance\
**Tags:** benchmark\
**Created:** [April 6, 2019, 5:55am UTC](https://discourse.julialang.org/t/fast-logsumexp/22827 "2019-04-06T05:55:42Z")\
**Posts on this page:** 1\
**Showing post:** 7

<div class="post-metadata">

**Author:** ![stevengj](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/stevengj/32/71_2.png) [@stevengj](https://discourse.julialang.org/u/stevengj)\
**Post date:** [April 22, 2019, 6:32am UTC](https://discourse.julialang.org/t/fast-logsumexp/22827/7 "2019-04-22T06:32:56Z")

</div>

I should comment that the offset algorithm that you are using for logsumexp is potentially problematic — you can easily contrive a case where it is off by a factor of 2. In particular, try `[1e-20, log(1e-20)]` with your logsumexp algorithm

```julia
function f(x)
    X = maximum(x)
    return X + log(sum(exp.(x .- X)))
end

```

You get:

```julia
julia> x = [1e-20, log(1e-20)];

julia> f(x) # inaccurate!
1.0e-20

julia> Float64(f(big.(x))) # accurate
1.9999999999999993e-20

```

A possible fix is to pull the maximum `x` term out of the sum and use the `log1p` function. I actually assigned this as an exam question recently, so you can see the explanation in [my problem 3 solutions](https://github.com/mitmath/18335/blob/spring19/psets/midtermsol.pdf).

Not sure if this matters for the machine-learning application, however, since there you are adding the logsumexp to a posterior probability and so errors in tiny values like this may get rounded away in your final result.

---

_[View the full topic](https://discourse.julialang.org/t/fast-logsumexp/22827)._