Fast Poisson Random Numbers in pure Julia
Switch branches/tags
Clone or download
Fetching latest commit…
Cannot retrieve the latest commit at this time.
Permalink
Failed to load latest commit information.
src
test
.codecov.yml
.gitignore
.travis.yml
LICENSE.md
README.md
REQUIRE
appveyor.yml

README.md

PoissonRandom

Build Status

Usage

Pkg.add("PoissonRandom")

# Simple Poisson random
pois_rand(λ)

# Using another RNG
using RandomNumbers
rng = Xorshifts.Xoroshiro128Plus()
pois_rand(λ,rng)

Implementation

It mixes two methods. A simple count method and a method from a normal approximation. See this blog post for details.

Benchmark

using RandomNumbers, Distributions, BenchmarkTools, StaticArrays,
      RecursiveArrayTools, Plots, PoissonRandom
labels = ["count_rand","ad_rand","pois_rand","Distributions.jl"]
rng = Xorshifts.Xoroshiro128Plus()

function n_count(λ,rng,n)
  tmp = 0
  for i in 1:n
    tmp += PoissonRandom.count_rand(λ,rng)
  end
end

function n_pois(λ,rng,n)
  tmp = 0
  for i in 1:n
    tmp += pois_rand(λ,rng)
  end
end

function n_ad(λ,rng,n)
  tmp = 0
  for i in 1:n
    tmp += PoissonRandom.ad_rand(λ,rng)
  end
end

function n_dist(λ,rng,n)
  tmp = 0
  for i in 1:n
    tmp += rand(Poisson(λ))
  end
end

function time_λ(λ,rng,n)
  t1 = @elapsed n_count(λ,rng,n)
  t2 = @elapsed n_ad(λ,rng,n)
  t3 = @elapsed n_pois(λ,rng,n)
  t4 = @elapsed n_dist(λ,rng,n)
  @SArray [t1,t2,t3,t4]
end

# Compile
time_λ(5,rng,5000000)
# Run with a bunch of λ
times = VectorOfArray([time_λ(n,rng,5000000) for n in 1:20])'
plot(times,labels = labels, lw = 3)

benchmark result

So this package ends up about 30% or so faster than Distributions.jl (the method at the far edge is λ-independent so that goes on forever).