Skip to content

Change PoissonADSampler threshold to 10 - #2060

Merged
devmotion merged 3 commits into
JuliaStats:masterfrom
k88hudson:k88/poisson_threshold
Jun 19, 2026
Merged

devmotion merged 3 commits into
JuliaStats:masterfrom
k88hudson:k88/poisson_threshold

Conversation

@k88hudson

Copy link
Copy Markdown
Contributor

Another one that came up when I was going through some of the older papers for sampling algorithms:

In the Poisson sampler, the threshold at which PoissonADSampler kicks in is currently μ >= 6, but the Ahrens-Dieter PD algorithm says it is designed for μ >= 10. I am not sure why 6 was chosen (it looks like it was originally from PoissonRandom.jl, see #1021 but I don't see anything beyond that). To verify that this is actually a problem: if you plot the hat function (equation 24) against the residual between the discretized normal and the poisson PMF, you will see that at 6 it definitely goes outside the bounds. With the following:

pmf_diff <- function(x, mu) {
  s <- sqrt(mu)
  dpois(floor(x), lambda = mu) - 
    pnorm(floor(x+1), mean = mu, sd = s) +
    pnorm(floor(x), mean = mu, sd = s)
}

hat <- function(t, mu) (0.1069 / mu) * exp(-abs(t - 1.8))
plot of the hat function vs the pmf

In terms of demonstrating the bias, right now the sampler tests only check down to 10.0, not 6.0:

(PoissonADSampler, [10.0, 15.0, 20.0, 30.0])]

You can detect it at 6, but I found that you need to use n=10^8 (I included a regression test in the PR which fails on master). At μ=7 it requires n=10^9, I didn't go further than that.

Note: I think you would expect performance regressions for this change for values of μ from 6 to 9, since PoissonADSampler should be faster down to 6.

Ahrens-Dieter PD (`PoissonADSampler`) is designed for μ ≥ 10; the hat
function bounds the residual between the discretized normal and the
poisson pmf correctly for that threshold, below it the rejection step
is no longer exact. Threshold now requires μ ≥ 10 for AD;
μ < 10 uses PoissonCountSampler.
@codecov-commenter

codecov-commenter commented Apr 25, 2026 •

Copy link
Copy Markdown

Codecov Report

✅ All modified and coverable lines are covered by tests.
✅ Project coverage is 86.63%. Comparing base (a7ed7b6) to head (b022b0e).
⚠️ Report is 5 commits behind head on master.

Additional details and impacted files
@@            Coverage Diff             @@
##           master    #2060      +/-   ##
==========================================
- Coverage   86.65%   86.63%   -0.02%     
==========================================
  Files         148      149       +1     
  Lines        8919     8921       +2     
==========================================
  Hits         7729     7729              
- Misses       1190     1192       +2     

☔ View full report in Codecov by Harness.
📢 Have feedback on the report? Share it here.

🚀 New features to boost your workflow:
  • ❄️ Test Analytics: Detect flaky tests, report on failures, and find test suite problems.

Comment thread test/univariate/discrete/poisson.jl Outdated
Comment on lines +10 to +12
let n = 100_000_000, μ = 6.0, tol = 5 * sqrt(μ / n)
@test abs(sum(rand(Poisson(μ)) for _ in 1:n) / n - μ) < tol
end

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Minor nitpick: Could you change this to a @testset? It would also be a bit more natural to use @test ... ≈ ... atol = ....

Comment thread src/univariate/discrete/poisson.jl
Co-authored-by: David Müller-Widmann <devmotion@users.noreply.github.com>

@devmotion devmotion left a comment

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Thank you!

@devmotion
devmotion merged commit c6c9d82 into JuliaStats:master Jun 19, 2026
11 checks passed
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Labels

None yet

Projects

None yet

Development

Successfully merging this pull request may close these issues.

3 participants