binomial_distribution overflows at 171 trials even when the result is representable #15656

Repository commit

6e03a48579ca212a08b36be00378f887dc41a003

Python version (python --version)

Python 3.11.0, macOS

Dependencies version (pip freeze)

No third-party dependencies are needed for this reproduction; both the implementation and reference use the Python standard library.

Expected behavior

The probability mass function should return a finite probability for valid inputs such as 85 successes in 171 fair trials. The expected value for this case is approximately 0.060749764626165485.

Actual behavior

Run from the repository root:

from math import comb
from maths.binomial_distribution import binomial_distribution

k, n, p = 85, 171, 0.5
print('reference:', comb(n, k) * p**k * (1-p)**(n-k))
print('actual:', binomial_distribution(k, n, p))

Output:

reference: 0.060749764626165485
Traceback (most recent call last):
  ...
OverflowError: int too large to convert to float

In maths/binomial_distribution.py, float(factorial(trials)) overflows for trials >= 171, even though the final probability for this input is representable. Using math.comb(trials, successes) avoids converting the much larger intermediate factorial to a float. For substantially larger inputs, an implementation should also consider underflow when multiplying the probability terms; the reference expression above is intended only for the concrete reproduction.

Searched existing issues and pull requests for binomial_distribution and binomial overflow; I did not find an existing report of this failure.