The Law of Small Numbers
The consecutive odds ratios of the binomial (n,p) distribution help us derive an approximation for the distribution when n is large and p is small. The approximation is sometimes called “the law of small numbers” because it approximates the distribution of the number of successes when the chance of success is small: you only expect a small number of successes.
As an example, here is the binomial (1000,2/1000) distribution. Note that 1000 is large, 2/1000 is pretty small, and 1000×(2/1000)=2 is the natural number of successes to be thinking about.
n = 1000
p = 2/1000
k = np.arange(16)
binom_probs = stats.binom.pmf(k, n, p)
binom_dist = Table().values(k).probability(binom_probs)
Plot(binom_dist)
Though the possible values of the number of successes in 1000 trials can be anywhere between 0 and 1000, the probable values are all rather small because p is small. That is why we didn’t even bother computing the probabilities beyond k=15.
Since the histogram is all scrunched up near 0, only very few bars have noticeable probability. It really should be possible to find or approximate the chances of the corresponding values by a simpler calculation than the binomial formula.
To see how to do this, we will start with P(0).
Approximation to P(0)
Let n→∞ and pn→0 in such a way that npn→μ>0. It’s important to ensure that pn doesn’t go to 0 so fast that npn→0 as well, because in that case all the probability just gets concentrated at the value 0 when n is large.
Let Pn(k) be the binomial (n,pn) probability of k successes.
Then
Pn(0)=(1−pn)n=(1−npnn)n→e−μ as n→∞If you can’t see the limit directly, appeal to our familiar exponential approxmation:
log(Pn(0))=nlog(1−npnn)=n⋅log(1−pn)∼n(−pn)=−npn∼−μwhen n is large, because pn∼0 and npn∼μ.
Approximation to P(k)
In general, for fixed k>1,
Pn(k)=Pn(k−1)Rn(k)=Pn(k−1)n−k+1k⋅pn1−pn=Pn(k−1)(npnk−(k−1)pnk)11−pn∼Pn(k−1)⋅μkwhen n is large, because k is constant, npn→μ, pn→0, and 1−pn→1. By induction, this implies the following approximation for each fixed k.
Pn(k) ∼ e−μ⋅μ1⋅μ2⋯μk = e−μμkk!if n is large, under all the additional conditions we have assumed. Here is a formal statement.
Poisson Approximation to the Binomial
Let n→∞ and pn→0 in such a way that npn→μ>0. Let Pn(k) be the binomial (n,pn) probability of k successes. Then for each k such that 0≤k≤n,
Pn(k)∼e−μμkk! for large nThis is called the Poisson approximation to the binomial. The parameter of the Poisson distribution is μ∼npn for large n.
The distribution is named after its originator, the French mathematician Siméon Denis Poisson (1781-1840).
The terms in the approximation are proportional to the terms in the series expansion of eμ:
μkk!, k≥0The expansion is infinite, but we are only going up to a finite (though large) number of terms n. You now start to see the value of being able to work with probability spaces that have an infinite number of possible outcomes.
We’ll get to that in a later section. For now, let’s see if the approximation we derived is any good.
Poisson Probabilities in Python
Use stats.poisson.pmf
just as you would use stats.binomial.pmf
, but keep in mind that the Poisson has only one parameter.
Suppose n=1000 and p=2/1000. Then the exact binomial chance of 3 successes is
stats.binom.pmf(3, 1000, 2/1000)
0.18062773231746918
The approximating Poisson distribution has parameter 1000×(2/1000)=2, and so the Poisson approximation to the probability above is
stats.poisson.pmf(3, 2)
0.18044704431548356
Not bad. To compare the entire distributions, first create the two distribution objects:
k = range(16)
bin_probs = stats.binom.pmf(k, 1000, 2/1000)
bin_dist = Table().values(k).probability(bin_probs)
poi_probs = stats.poisson.pmf(k, 2)
poi_dist = Table().values(k).probability(poi_probs)
The prob140
function that draws overlaid histograms is called Plots
(note the plural). The syntax has alternating arguments: a string label you provide for a distribution, followed by that distribution, then a string label for the second distribution, then that distribution.
Plots('Binomial (1000, 2/1000)', bin_dist, 'Poisson (2)', poi_dist)
Does it look as though there is only one histogram? That’s because the approximation is great! Here are the two histograms individually.
Plot(bin_dist)
Plot(poi_dist)
In lab, you will use total variation distance to get a bound on the error in the approximation.
A reasonable question to ask at this stage is, “Well that’s all very nice, but why should I bother with approximations when I can just use Python to compute the exact binomial probabilities using stats.binom.pmf
?”
Part of the answer is that if a function involves parameters, you can’t understand how it behaves by just computing its values for some particular choices of the parameters. In the case of Poisson probabilities, we will also see shortly that they form a powerful distribution in their own right, on an infinite set of values.