Turing.jl: logsum -> logsumexp

[mohamed82008]

https://github.com/TuringLang/Turing.jl/blob/70640d35805bb46a54c04d9ff25e5a90f22eee2b/src/core/util.jl#L23

I believe this is the famous log-sum-exp function.

[mohamed82008]

Can also use a mapreduce there to avoid allocating ys.

[willtebbutt]

Huh. We should probably be using a standard implementation for this stuff eg. that in StatsFuns.jl
This probably also applies to logit and invlogit / logistic.
See here for their implementations.

[mohamed82008]

Apparently, StatsFuns is not registered so I am just renaming logsum to logsumexp in Turing.

[trappmartin]

Maybe we could do the performance improvement you proposed while renaming? Log-sum-exp can be quite a performance bottleneck for certain types of models, e.g. deep mixtures.

Or maybe something like this? But it might be a bit overly complicated. =)

function logsumexp(x::AbstractVector{T}) where {T <: Real}
  alpha = map(T, -Inf)
  r = zero(T)

  @inbounds for i in 1:length(x)
      if isinf(x[i])
        continue
      elseif x[i] <= alpha
        r += exp(x[i] - alpha)
      else
        r *= exp(alpha - x[i])
        r += one(T)
        alpha = x[i]
      end
    end
    return log(r) + alpha
end

[yebai]

Apparently, StatsFuns is not registered so I am just renaming logsum to logsumexp in Turing.

It seems the following command works on my computer:

(v1.0) pkg> add StatsFuns
  Updating registry at `~/.julia/registries/General`
  Updating git-repo `https://github.com/JuliaRegistries/General.git`
 Resolving package versions...
 Installed Graphics ─ v0.4.0
  Updating `~/.julia/environments/v1.0/Project.toml`
  [4c63d2b9] + StatsFuns v0.7.0
  Updating `~/.julia/environments/v1.0/Manifest.toml`
  [a2bd30eb] ↑ Graphics v0.3.0 ⇒ v0.4.0

[mohamed82008]

hmm that's weird, it doesn't work on mine. I will upgrade everything and try again.

[mohamed82008]

Actually now it works no idea what changed. Anyways, back to plan A then.

[nalimilan]

Note that with https://github.com/JuliaStats/StatsFuns.jl/pull/97, StatsFuns is going to implement a faster version of logsumexp, so it would be worth using it instead of an internal implementation.

[devmotion]

I replaced the internal implementation with StatsFuns.logsumexp last year (https://github.com/TuringLang/Turing.jl/commit/6e5843517d71cd55429afdeb6abcb97a0bbe7361).