# Speeding up the calculation of a binomial sum

Is it possible to speed up this calculation?

A[j_, p_, n_] := Sum[Binomial[n, k] p^k (1 - p)^(n - k), {k, 0, j}]
Plot[{A[j, 0.5, 12000]}, {j, 0, 12000}]


Thanks!

Mathematica complains due to the occurence of very small numbers that cannot be represented in machine precision. However, the built-in cummulative distribution function (CDF) of the BinomialDistribution is sufficiently robust against that (replacing very small numbers simply by 0.).

A[j_, p_, n_] := CDF[BinomialDistribution[n, p], j]
nn = 120;
Plot[{A[j, 0.5, nn]}, {j, 0, nn}]


nn = 12000;
Plot[{A[j, 0.5, nn]}, {j, 0, nn}]


The summation can be evaluated to a closed form, which can speed up the further usage:

Sum[Binomial[n, k] p^k (1 - p)^(n - k), {k, 0, j}] // FullSimplify

(1 - p)^n ((1/(1 - p))^n - (1 - p)^(-1 - j) p^(1 + j)
Binomial[n, 1 + j] Hypergeometric2F1[1, 1 + j - n, 2 + j,
p/(-1 + p)])


The expression inside the Sum[] is actually built into Mathematica as BernsteinBasis[], which gets primarily used in the construction of Bézier curves. Applied to this situation, we have the formula

A[j_, p_, n_] := Sum[BernsteinBasis[n, k, p], {k, 0, j}, Method -> "Procedural"]


which ought to be very stable to evaluate. If needed, however, one can also implement the de Casteljau formula so that the Bernstein polynomials are recursively generated instead.

(back to gedanken Mathematica, so untested)