# Defining the Poisson distribution

The Poisson distribution is a family of discrete distributions with positive probabilities on the non-negative numbers $0,1,2,\cdots$. Each distribution in this family is indexed by a positive number $\lambda>0$. One way to define this distribution is to give its probability function given the parameter $\lambda$ and then derive various distributional quantities such as mean and variance. Along with other mathematical facts, it can be shown that both the mean and the variance are $\lambda$. In this post, we take a different tack. We look at two view points that give rise to the Poisson distribution. Taking this approach will make it easier to appreciate some of the possible applications of the Poisson distribution. The first view point is that the Poisson distribution is the limiting case of the binomial distribution. The second view point is through the Poisson process, a stochastic process that, under some conditions, counts the number of events and the time points at which these events occur in a given time (or physical) interval.

________________________________________________________________________

Poisson as a limiting case of binomial

A binomial distribution where the number of trials $n$ is large and the probability of success $p$ is small such that $np$ is moderate in size can be approximated using the Poisson distribution with mean $\lambda=np$. This fact follows from Theorem 1, which indicates that the Poisson distribution is the limiting case of the binomial distribution.

Theorem 1
Let $\lambda$ be a fixed positive constant. Then for each integer $x=0,1,2,\cdots$, the following is true:

$\displaystyle \lim_{n \rightarrow \infty} \binom{n}{x} \ p^x \ (1-p)^{n-x}=\lim_{n \rightarrow \infty} \frac{n!}{x! \ (n-x)!} \ p^x \ (1-p)^{n-x}=\frac{e^{-\lambda} \ \lambda^x}{x!}$

where $\displaystyle p=\frac{\lambda}{n}$.

Proof of Theorem 1
We start with a binomial distribution with $n$ trials and with $\displaystyle p=\frac{\lambda}{n}$ being the probability of success, where $n>\lambda$. Let $X_n$ be the count of the number of successes in these $n$ Bernoulli trials. The following is the probability that $X_n=k$.

\displaystyle \begin{aligned} P(X_n=k)&=\binom{n}{k} \biggl(\frac{\lambda}{n}\biggr)^k \biggr(1-\frac{\lambda}{n}\biggr)^{n-k} \\&=\frac{n!}{k! (n-k)!} \biggl(\frac{\lambda}{n}\biggr)^k \biggr(1-\frac{\lambda}{n}\biggr)^{n-k} \\&=\frac{n(n-1)(n-2) \cdots (n-k+1)}{n^k} \biggl(\frac{\lambda^k}{k!}\biggr) \biggr(1-\frac{\lambda}{n}\biggr)^{n} \biggr(1-\frac{\lambda}{n}\biggr)^{-k} \\&=\biggl(\frac{\lambda^k}{k!}\biggr) \ \biggl[ \frac{n(n-1)(n-2) \cdots (n-k+1)}{n^k} \ \biggr(1-\frac{\lambda}{n}\biggr)^{n} \ \biggr(1-\frac{\lambda}{n}\biggr)^{-k} \biggr] \end{aligned}

In the last step, the terms that contain $n$ are inside the square brackets. Let’s see what they are when $n$ approaches infinity.

$\displaystyle \lim \limits_{n \rightarrow \infty} \ \frac{n(n-1)(n-2) \cdots (n-k+1)}{n^k}=1$

$\displaystyle \lim \limits_{n \rightarrow \infty} \biggr(1-\frac{\lambda}{n}\biggr)^{n}=e^{-\lambda}$

$\displaystyle \lim \limits_{n \rightarrow \infty} \biggr(1-\frac{\lambda}{n}\biggr)^{-k}=1$

The reason that the first result is true is that the numerator is a polynomial where the leading term is $n^k$. Upon dividing by $n^k$ and taking the limit, we get 1. The second result is true since the following limit is one of the definitions of the exponential function $e^x$.

$\displaystyle \lim \limits_{n \rightarrow \infty} \biggr(1+\frac{x}{n}\biggr)^{n}=e^{x}$

The third result is true since the exponent $-k$ is a constant. Thus the following is the limit of the probability $P(X_n=k)$ as $n \rightarrow \infty$.

\displaystyle \begin{aligned} \lim \limits_{n \rightarrow \infty} P(X_n=k)&= \biggl(\frac{\lambda^k}{k!}\biggr) \ \lim \limits_{n \rightarrow \infty} \biggl[ \frac{n(n-1)(n-2) \cdots (n-k+1)}{n^k} \ \biggr(1-\frac{\lambda}{n}\biggr)^{n} \ \biggr(1-\frac{\lambda}{n}\biggr)^{-k} \biggr] \\&=\biggl(\frac{\lambda^k}{k!}\biggr) \cdot 1 \cdot e^{-\lambda} \cdot 1 \\&=\frac{e^{-\lambda} \lambda^k}{k!} \end{aligned}

This above derivation completes the proof. $\blacksquare$

In a given binomial distribution, whenever the number of trials $n$ is large and the probability $p$ of success in each trial is small (i.e. each of the Bernoulli trial rarely results in a success), Theorem 1 tells us that we can use the Poisson distribution with parameter $\lambda=np$ to estimate the binomial distribution.

Example 1
The probability of being dealt a full house in a hand of poker is approximately 0.001441. Out of 5000 hands of poker that are dealt at a certain casino, what is the probability that there will be at most 4 full houses?

Let $X$ be the number of full houses in these 5000 poker hands. The exact distribution for $X$ is the binomial distribution with $n=$ 5000 and $p=$ 0.001441. Thus example deals with a large number of trials where each trial is a rare event. So the Poisson estimation is applicable. Let $\lambda=$ 5000(0.001441) = 7.205. Then $P(X \le 4)$ can be approximated by the Poisson random variable $Y$ with parameter $\lambda$. The following is the probability function of $Y$:

$\displaystyle P(Y=y)=e^{-7.205} \ \frac{7.205^y}{y!}$

The following is the approximation of $P(X \le 4)$:

\displaystyle \begin{aligned} P(X \le 4)&\approx P(Y \le 4) \\&=P(Y=0)+P(Y=1)+P(Y=2)+P(Y=3)+P(Y=4) \\&= e^{-7.205} \biggl[ 1+7.205+\frac{7.205^2}{2!}+\frac{7.205^3}{3!}+\frac{7.205^4}{4!}\biggr] \\&=0.155098087 \end{aligned}

The following is a side by side comparison between the binomial distribution and its Poisson approximation. For all practical purposes, they are indistingusihable from one another.

$\displaystyle \begin{bmatrix} \text{Count of}&\text{ }&\text{ }&\text{Binomial } &\text{ }&\text{ }&\text{Poisson } \\\text{Full Houses}&\text{ }&\text{ }&P(X \le x) &\text{ }&\text{ }&P(Y \le x) \\\text{ }&\text{ }&\text{ }&n=5000 &\text{ }&\text{ }&\lambda=7.205 \\\text{ }&\text{ }&\text{ }&p=0.001441 &\text{ }&\text{ }&\text{ } \\\text{ }&\text{ }&\text{ } &\text{ }&\text{ } \\ 0&\text{ }&\text{ }&0.000739012&\text{ }&\text{ }&0.000742862 \\ 1&\text{ }&\text{ }&0.006071278&\text{ }&\text{ }&0.006095184 \\ 2&\text{ }&\text{ }&0.025304641&\text{ }&\text{ }&0.025376925 \\ 3&\text{ }&\text{ }&0.071544923&\text{ }&\text{ }&0.071685238 \\ 4&\text{ }&\text{ }&0.154905379&\text{ }&\text{ }&0.155098087 \\ 5&\text{ }&\text{ }&0.275104906&\text{ }&\text{ }&0.275296003 \\ 6&\text{ }&\text{ }&0.419508250&\text{ }&\text{ }&0.419633667 \\ 7&\text{ }&\text{ }&0.568176421 &\text{ }&\text{ }&0.568198363 \\ 8&\text{ }&\text{ }&0.702076190 &\text{ }&\text{ }&0.701999442 \\ 9&\text{ }&\text{ }&0.809253326&\text{ }&\text{ }&0.809114639 \\ 10&\text{ }&\text{ }&0.886446690&\text{ }&\text{ }&0.886291139 \\ 11&\text{ }&\text{ }&0.936980038&\text{ }&\text{ }&0.936841746 \\ 12&\text{ }&\text{ }&0.967298041&\text{ }&\text{ }&0.967193173 \\ 13&\text{ }&\text{ }&0.984085073&\text{ }&\text{ }&0.984014868 \\ 14&\text{ }&\text{ }&0.992714372&\text{ }&\text{ }&0.992672033 \\ 15&\text{ }&\text{ }&0.996853671&\text{ }&\text{ }&0.996830358 \end{bmatrix}$

The above table is calculated using the functions BINOM.DIST and POISSON.DIST in Excel. The following shows how it is done. The parameter TRUE indicates that the result is a cumulative distribution. When it is set to FALSE, the formula gives the probability function.

$P(X \le x)=\text{BINOM.DIST(x, 5000, 0.001441, TRUE)}$

$P(Y \le x)=\text{POISSON.DIST(x, 7.205, TRUE)}$

________________________________________________________________________

The Poisson distribution

The limit in Theorem 1 is a probability function and the resulting distribution is called the Poisson distribution. We now gives the formal definition. A random variable $X$ that takes on one of the numbers $0,1,2,\cdots$ is said to be a Poisson random variable with parameter $\lambda>0$ if

$\displaystyle P(X=x)=\frac{e^{-\lambda} \ \lambda^x}{x!} \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ x=0,1,2,\cdots$

It can be shown that the above function is indeed a probability function, i.e., the probabilities sum to 1. Any random variable that has a probability function of the above form is said to follow (or to have) a Poisson distribution. Furthermore, it can be shown that $E(X)=var(X)=\lambda$, i.e., the Poisson parameter is both the mean and variance. Thus the Poisson distribution may be a good fit if the observed data indicate that the sample mean and the sample variance are nearly identical.

The following is the moment generating function of the Poisson distribution with parameter $\lambda$.

$\displaystyle M(t)=E(e^{tX})=e^{\lambda \ (e^t-1)}$

One consequence of the Poisson moment generating function is that any independent sum of Poisson distributions is again a Poisson distribution.

________________________________________________________________________

The Poisson process

Another way, the more important way, to look at the Poisson distribution is the view point of the Poisson process. Consider an experiment in which events that are of interest occur at random in a time interval. The goal here is to record the time of the occurrence of each random event and for the purpose at hand, count the number of random events occurring in a fixed time interval. Starting at time 0, note the time of the occurrence of the first event. Then the time at which the second random event occurs and so on. Out of these measurements, we can derive the length of time between the occurrences of any two consecutive random events. Such measurements belong to a continuous random variable. In this post, we focus on the discrete random variable of the count of the random events in a fixed time interval.

A good example of a Poisson process is the well known experiment in radioactivity conducted by Rutherford and Geiger in 1910. In this experiment, $\alpha$-particles were emitted from a polonium source and the number of $\alpha$-particles were counted during an interval of 7.5 seconds (2608 many such time intervals were observed). A Poisson process is a random process in which several criteria are satisfied. We will show that in a Poisson process, the number of these random occurrences in the fixed time interval will follow a Poisson distribution. First, we discuss the criteria to which a Poisson process must conform.

One of the criteria is that in a very short time interval, the chance of having more than one random event is essentially zero. So either one random event will occur or none will occur in a very short time interval. Considering the occurrence of a random event as a success, there is either a success or a failure in a very short time interval. So a very short time interval in a Poisson process can be regarded as a Bernoulli trial.

The second criterion is that the experiment remains constant over time. Specifically this means that the probability of a random event occurring in a given subinterval is proportional to the length of that subinterval and not on where the subinterval is in the original interval. For example, in the 1910 radioactivity study, $\alpha$-particles were emitted at the rate of $\lambda=$ 3.87 per 7.5 seconds. So the probability of one $\alpha$-particle emitted from the radioactive source in a one-second interval is 3.87/7.5 = 0.516. Then the probability of observing one $\alpha$-particle in a half-second interval is 0.516/2 = 0.258. For a quarter-second interval, the probability is 0.258/2 = 0.129. So if we observe half as long, it will be half as likely to observe the occurrence of a random event. On the other hand, it does not matter when the quarter-second subinterval is, whether at the beginning or toward the end of the original interval of 7.5 seconds.

The third criterion is that non-overlapping subintervals are mutually independent in the sense that what happens in one subinterval (i.e. the occurrence or non-occurrence of a random event) will have no influence on the occurrence of a random event in another subinterval. To summarize, the following are the three criteria of a Poisson process:

Suppose that on average $\lambda$ random events occur in a time interval of length 1.

1. The probability of having more than one random event occurring in a very short time interval is essentially zero.
2. For a very short subinterval of length $\frac{1}{n}$ where $n$ is a sufficiently large integer, the probability of a random event occurring in this subinterval is $\frac{\lambda}{n}$.
3. The numbers of random events occurring in non-overlapping time intervals are independent.

Consider a Poisson process in which the average rate is $\lambda$ random events per unit time interval. Let $Y$ be the number of random events occurring in the unit time interval. In the 1910 radioactivity study, the unit time interval is 7.5 seconds and $Y$ is the count of the number of $\alpha$-particles emitted in 7.5 seconds. It follows that $Y$ has a Poisson distribution with parameter $\lambda$. To see this, subdivide the unit interval into $n$ non-overlapping subintervals of equal length where $n$ is a sufficiently large integer. Let $X_{n,j}$ be the number of random events in the the $j$th subinterval ($1 \le j \le n$). Based on the three assumptions, $X_{n,1},X_{n,2},\cdots,X_{n,n}$ are independent Bernoulli random variables, where the probability of success for each $X_{n,j}$ is $\frac{\lambda}{n}$. Then $X_n=X_{n,1}+X_{n,2}+\cdots+X_{n,n}$ has a binomial distribution with parameters $n$ and $p=\frac{\lambda}{n}$. Theorem 1 tells us that the limiting case of the binomial distributions for $X_n$ is the Poisson distribution with parameter $\lambda$. This Poisson distribution should agree with the distribution for $Y$. The Poisson is also discussed in quite a lot of details in the previous post called Poisson as a Limiting Case of Binomial Distribution.

We now examine the 1910 radioactivity study a little more closely.

Example 2
The basic idea of the 1910 radioactivity study conducted by Rutherford and Geiger is that a polonium source was placed a short distance from an observation point. The number of $\alpha$-particles emitted from the source were counted in 7.5-second intervals for 2608 times. The following is the tabulated results.

$\displaystyle \begin{bmatrix} \text{Number of alpha particles}&\text{ }&\text{Observed} \\ \text{recorded per 7.5 seconds }&\text{ }&\text{counts} \\ \text{ }&\text{ }&\text{ } \\ 0&\text{ }&57 \\ 1&\text{ }&203 \\ 2&\text{ }&383 \\ 3&\text{ }&525 \\ 4&\text{ }&532 \\ 5&\text{ }&408 \\ 6&\text{ }&273 \\ 7&\text{ }&139 \\ 8&\text{ }&45 \\ 9&\text{ }&27 \\ 10&\text{ }&10 \\ 11+&\text{ }&6 \\ \text{ }&\text{ }&\text{ } \\ \text{Total }&\text{ }&2608 \end{bmatrix}$

What is the average number of particles observed per 7.5 seconds? The total number of $\alpha$-particles in these 2608 periods is

$0 \times 57+1 \times 203+2 \times 383+ 3 \times 525 + \cdots=10097$.

The mean count per period is $\lambda=\frac{10097}{2608}=3.87$. Consider the Poisson distribution with parameter 3.87. The following is its probability function.

$\displaystyle P(X=x)=\frac{e^{-3.87} \ 3.87^x}{x!} \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ x=0,1,2,\cdots$

Out of 2608 periods, the expected number of periods with $x$ particles in emission is $2608P(X=x)$. The following is a side by side comparison in the observed counts and the expected counts.

$\displaystyle \begin{bmatrix} \text{Number of alpha particles}&\text{ }&\text{Observed}&\text{ }&\text{Expected} \\ \text{recorded per 7.5 seconds }&\text{ }&\text{counts}&\text{ }&\text{counts} \\ \text{ }&\text{ }&\text{ }&\text{ }&2608 \times P(X=x) \\ \text{ }&\text{ }&\text{ }&\text{ }&\text{ } \\ 0&\text{ }&57&\text{ }&54.40 \\ 1&\text{ }&203&\text{ }&210.52 \\ 2&\text{ }&383&\text{ }&407.36 \\ 3&\text{ }&525&\text{ }&525.50 \\ 4&\text{ }&532&\text{ }&508.42 \\ 5&\text{ }&408&\text{ }&393.52 \\ 6&\text{ }&273&\text{ }&253.82 \\ 7&\text{ }&139&\text{ }&140.32 \\ 8&\text{ }&45&\text{ }&67.88 \\ 9&\text{ }&27&\text{ }&29.19 \\ 10&\text{ }&10&\text{ }&11.30 \\ 11+&\text{ }&6&\text{ }&5.78 \\ \text{ }&\text{ }&\text{ }&\text{ }&\text{ } \\ \text{Total }&\text{ }&2608&\text{ }&2608 \end{bmatrix}$

The expected counts are quite close to the observed counts, showing that the Poisson distribution is a very good fit to the observed data from the 1910 study.

________________________________________________________________________

We have described the Poisson process as the distribution of random events in a time interval. The same idea can be used to describe random events occurring along a spatial interval, i.e. intervals in terms of distance or volume or other spatial measurements (see Examples 5 and 6 below).

Another point to make is that sometimes it may be necessary to consider an interval other than the unit length. Instead of counting the random events occurring in an interval of length 1, we may want to count the random events in an interval of length $t$. As before, let $\lambda$ be the rate of occurrences in a unit interval. Then the rate of occurrences of the random events is over the interval of length $t$ is $\lambda t$. The same idea will derive that fact that the number of occurrences of the random events of interest in the interval of length $t$ is a Poisson distribution with parameter $\lambda t$. The following is its probability function.

$\displaystyle P(X_t=x)=\frac{e^{-\lambda t} \ (\lambda t)^x}{x!} \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ x=0,1,2,\cdots$

where $X_t$ is the count of the random events in an interval of length $t$.

For example, in the 1910 radioactive study, the unit length is 7.5 seconds. The number of $\alpha$-particles observed in a half unit interval (3.75 seconds) will follow a Poisson distribution with parameter $0.5 \lambda=$ 0.5(3.87) = 1.935 with the following probability function:

$\displaystyle P(X_{0.5}=x)=\frac{e^{-1.935} \ (1.935)^x}{x!} \ \ \ \ \ \ \ \ \ \ \ \ \ x=0,1,2,\cdots$

________________________________________________________________________

More examples

Example 3
A radioactive source is metered for 5 hours. During this period, 9638 $\alpha$-particles are counted. What is the probability that during the next minute, between 30 and 34 particles (both inclusive ) will be counted?

The average number of $\alpha$-particles counted per minute is $\lambda=\frac{9638}{300}=32.12$. Let $X$ be the number of $\alpha$-particles counted per minute. Then $X$ has a Poisson distribution with parameter $\lambda=32.12$. The following calculates $P(30 \le X \le 34)$.

\displaystyle \begin{aligned} P(30 \le X \le 34)&=e^{-32.12} \biggl[ \frac{32.12^{30}}{30!}+\frac{32.12^{31}}{31!}+\frac{32.12^{32}}{32!}+\frac{32.12^{33}}{33!}+\frac{32.12^{34}}{34!} \biggr] \\&=0.341118569 \end{aligned}

Alternatively, the POISSON.DIST function in Excel can be used as follows:

\displaystyle \begin{aligned} P(30 \le X \le 34)&=P(X \le 34)-P(X \le 29) \\&=\text{POISSON.DIST(34,32.12,TRUE)} \\& \ \ \ \ \ \ -\text{POISSON.DIST(29,32.12,TRUE)} \\&=0.671501917-0.330383348 \\&=0.341118569 \end{aligned}

Example 4
The side effect of dry mouth is known to be experienced, on the average, by 5 out of 10,000 individuals taking a certain medication. About 20,000 patients are expected to take this medication next year. What is the probability that between 12 and 16 (both inclusive) patients will experience the side effect of dry mouth? What is the exact probability model that can also be used to work this problem?

The exact model is a binomial distribution. The number of trials $n=$ 20000 and the probability of success in each trial is $p=$ 0.0005 (experiencing the side effect). Here, we use Poisson to estimate the binomial. The average number of patients experiencing side effect is $\lambda=20000(0.0005)=10$. Let $X$ be the number of patients experiencing the side effect. The following calculates the Poisson probability for $P(12 \le X \le 16)$ in two different ways.

\displaystyle \begin{aligned} P(12 \le X \le 16)&=e^{-10} \biggl[ \frac{10^{12}}{12!}+\frac{10^{13}}{13!}+\frac{10^{14}}{14!}+\frac{10^{15}}{15!}+\frac{10^{16}}{16!} \biggr] \\&=0.276182244 \end{aligned}
\displaystyle \begin{aligned} P(12 \le X \le 16)&=P(X \le 11)-P(X \le 16) \\&=\text{POISSON.DIST(16,10,TRUE)} \\& \ \ \ \ \ \ -\text{POISSON.DIST(11,10,TRUE)} \\&=0.97295839-0.696776146 \\&=0.276182244 \end{aligned}

Example 5
In a 10-mile stretch of a highway, car troubles (e.g. tire punctures, dead batteries, and mechanical breakdown) occur at a rate of 1.5 per hour. A tow truck driver can respond to such car troubles and offer roadside assistance, which can include towing and minor repair. Assume that the number of such incidences per hour follows a Poisson distribution. At the beginning of the hour, three tow trucks (and their drivers) are available to respond to any car troubles in this stretch of highway. What is the probability that in the next hour all three tow trick drivers will be busy helping motorists with car troubles in this stretch of highway?

Let $X$ be the number of car troubles that occur in this 10-mile stretch of highway in the one-hour period in question. If in this one hour there are 3 or more car troubles ($X \ge 3$), then all three tow truck drivers will be busy.

\displaystyle \begin{aligned} P(X \ge 3)&=1-P(X \le 2) \\&=1-e^{-1.5} \biggl[ 1+1.5+\frac{1.5^{2}}{2!} \biggr] \\&=1-0.808846831\\&=0.191153169 \end{aligned}

Example 6
Continuing Example 5. Considering that there is only 19% chance that all 3 tow truck drivers will be busy, there is a good chance that the resources are under utilized. What if one of the drivers is assigned to another stretch of highway?

With only two tow trucks available for this 10-mile stretch of highway, the following is the probability that all two tow truck drivers will be busy:

\displaystyle \begin{aligned} P(X \ge 2)&=1-P(X \le 1) \\&=1-e^{-1.5} \biggl[ 1+1.5 \biggr] \\&=1-0.5578254\\&=0.4421746 \end{aligned}

Assigning one driver to another area seems to be a better way to make good use of the available resources. With only two tow truck drivers available, there is much reduced chance (56%) that one of the drivers will be idle, and there is a much increased chance (44%) that all available drivers will be busy.

________________________________________________________________________

Remarks

The Poisson distribution is one of the most important of all probability models and has shown to be an excellent model for a wide array of phenomena such as

• the number of $\alpha$-particles emitted from radioactive source in a given amount of time,
• the number of vehicles passing a particular location on a busy highway,
• the number of traffic accidents in a stretch of highway in a given period of time,
• the number of phone calls arriving at a particular point in a telephone network in a fixed time period,
• the number of insurance losses/claims in a given period of time,
• the number of customers arriving at a ticket window,
• the number of earthquakes occurring in a fixed period of time,
• the number of mutations on a strand of DNA.
• the number of hurricanes in a year that originate in the Atlantic ocean.

What is the Poisson distribution so widely applicable in these and many other seemingly different and diverse phenomena? What is the commonality that ties all these different and diverse phenomena? The commonality is that all these phenomena are basically a series of independent Bernoulli trials. If a phenomenon is a Binomial model where $n$ is large and $p$ is small, then it has a strong connection to Poisson model mathematically through Theorem 1 above (i.e. it has a Poisson approximation). On the other hand, if the random phenomenon follows the criteria in a Poisson process, then the phenomenon is also approximately a Binomial model, which means that in the limiting case it is Poisson.

In both view points discussed in this post, the Poisson distribution can be regarded as a binomial distribution taken at a very granular level. This connection with the binomial distribution points to a vast arrays of problems that can be solved using the Poisson distribution.

________________________________________________________________________

Exercises

Practice problems for the Poisson concepts discussed above can be found in the companion blog (go there via the following link). Working on these exercises is strongly encouraged (you don’t know it until you can do it).

________________________________________________________________________
$\copyright \ \text{2015 by Dan Ma}$

# Poisson as a Limiting Case of Binomial Distribution

In many binomial problems, the number of Bernoulli trials $n$ is large, relatively speaking, and the probability of success $p$ is small such that $n p$ is of moderate magnitude. For example, consider problems that deal with rare events where the probability of occurrence is small (as a concrete example, counting the number of people with July 1 as birthday out of a random sample of 1000 people). It is often convenient to approximate such binomial problems using the Poisson distribution. The justification for using the Poisson approximation is that the Poisson distribution is a limiting case of the binomial distribution. Now that cheap computing power is widely available, it is quite easy to use computer or other computing devices to obtain exact binomial probabiities for experiments up to 1000 trials or more. Though the Poisson approximation may no longer be necessary for such problems, knowing how to get from binomial to Poisson is important for understanding the Poisson distribution itself.

Consider a counting process that describes the occurrences of a certain type of events of interest in a unit time interval subject to three simplifying assumptions (discussed below). We are interested in counting the number of occurrences of the event of interest in a unit time interval. As a concrete example, consider the number of cars arriving at an observation point in a certain highway in a period of time, say one hour. We wish to model the probability distribution of how many cars that will arrive at the observation point in this particular highway in one hour. Let $X$ be the random variable described by this probability distribution. We wish to konw the probability that there are $k$ cars arriving in one hour. We start with using a binomial distribution as an approximation to the probability $P(X=k)$. We will see that upon letting $n \rightarrow \infty$, the $P(X=k)$ is a Poisson probability.

Suppose that we know $E(X)=\alpha$, perhaps an average obtained after observing cars at the observation points for many hours. The simplifying assumptions alluded to earlier are the following:

1. The numbers of cars arriving in nonoverlapping time intervals are independent.
2. The probability of one car arriving in a very short time interval of length $h$ is $\alpha h$.
3. The probability of having more than one car arriving in a very short time interval is esstentially zero.

Assumption 1 means that a large number of cars arriving in one period does not imply fewer cars will arrival in the next period and vice versa. In other words, the number of cars that arrive in any one given moment does not affect the number of cars that will arrive subsequently. Knowing how many cars arriving in one minute will not help predict the number of cars arriving at the next minute.

Assumption 2 means that the rate of cars arriving is dependent only on the length of the time interval and not on when the time interval occurs (e.g. not on whether it is at the beginning of the hour or toward the end of the hour).

Assumptions 3 allows us to think of a very short period of time as a Bernoulli trial. Thinking of the arrival of a car as a success, each short time interval will result in only one success or one failure.

Assumption 2 tells us that non-overlapping short time intervals of the same length have identical probability of success. Overall, all three assumptions imply that non-overlapping short intervals of the same length are inpdendent Bernoulli trials with identical probability of success, giving us the basis for applying binomial distribution.

To start, we can break up the hour into 60 minutes (into 60 Bernoulli trials). We then consider the binomial distribution with $n=60$ and $p=\frac{\alpha}{60}$. So the following is an approximation to our desired probability distribution.

$\displaystyle (1) \ \ \ \ \ P(X=k) \approx \binom{60}{k} \biggl(\frac{\alpha}{60}\biggr)^k \biggr(1-\frac{\alpha}{60}\biggr)^{60-k} \ \ \ \ \ k=0,1,2,\cdots, 60$

Conceivably, there can be more than 1 car arriving in a minute and observing cars in a one-minute interval may not be a Bernoulli trial. For a one-minute interval to qualify as a Bernoulli trial, there is either no car arriving or 1 car arriving in that one minute. So we can break up an hour into 3,600 seconds (into 3,600 Bernoulli trials). We now consider the binomial distribution with $n=3600$ and $p=\frac{\alpha}{3600}$.

$\displaystyle (2) \ \ \ \ \ P(X=k) \approx \binom{3600}{k} \biggl(\frac{\alpha}{3600}\biggr)^k \biggr(1-\frac{\alpha}{3600}\biggr)^{3600-k} \ \ \ \ \ k=0,1,2,\cdots, 3600$

It is also conceivable that more than 1 car can arrive in one second and observing cars in one-second interval may still not qualify as a Bernoulli trial. So we need to get more granular. We can divide up the hour into $n$ equal subintervals, each of length $\frac{1}{n}$. Assumption 3 ensures that each subinterval is a Bernoulli trial (either it is a success or a failure; one car arriving or no car arriving). Assumption 2 ensures that the non-overlapping subintervals all have the same probability of success. Assumption 1 tells us that all the $n$ subintervals are independent. So breaking up the hour into $n$ moments and counting the number of moments that are successes will result in a binomial distribution with parameters $n$ and $p=\frac{\alpha}{n}$. So we are ready to proceed with the following approximation to our probability distribution $P(X=k)$.

$\displaystyle (3) \ \ \ \ \ P(X=k) \approx \binom{n}{k} \biggl(\frac{\alpha}{n}\biggr)^k \biggr(1-\frac{\alpha}{n}\biggr)^{n-k} \ \ \ \ \ k=0,1,2,\cdots, n$

As we get more granular, $n \rightarrow \infty$. We show that the limit of the binomial probability in $(3)$ is the Poisson distribution with parameter $\alpha$. We show the following.

$\displaystyle (4) \ \ \ \ \ P(X=k) = \lim \limits_{n \rightarrow \infty} \binom{n}{k} \biggl(\frac{\alpha}{n}\biggr)^k \biggr(1-\frac{\alpha}{n}\biggr)^{n-k}=\frac{e^{-\alpha} \alpha^k}{k!} \ \ \ \ \ \ k=0,1,2,\cdots$

In the derivation of $(4)$, we need the following two mathematical tools. The statement $(5)$ is one of the definitions of the mathematical constant $e$. In the statement $(6)$, the integer $n$ in the numerator is greater than the integer $k$ in the denominator. It says that whenever we work with such a ratio of two factorials, the result is the product of $n$ with the smaller integers down to $(n-(k-1))$. There are exactly $k$ terms in the product.

$\displaystyle (5) \ \ \ \ \ \lim \limits_{n \rightarrow \infty} \biggl(1+\frac{x}{n}\biggr)^n=e^x$

$\displaystyle (6) \ \ \ \ \ \frac{n!}{(n-k)!}=n(n-1)(n-2) \cdots (n-k+1) \ \ \ \ \ \ \ \ k

The following is the derivation of $(4)$.

\displaystyle \begin{aligned}(7) \ \ \ \ \ P(X=k)&=\lim \limits_{n \rightarrow \infty} \binom{n}{k} \biggl(\frac{\alpha}{n}\biggr)^k \biggr(1-\frac{\alpha}{n}\biggr)^{n-k} \\&=\lim \limits_{n \rightarrow \infty} \ \frac{n!}{k! (n-k)!} \biggl(\frac{\alpha}{n}\biggr)^k \biggr(1-\frac{\alpha}{n}\biggr)^{n-k} \\&=\lim \limits_{n \rightarrow \infty} \ \frac{n(n-1)(n-2) \cdots (n-k+1)}{n^k} \biggl(\frac{\alpha^k}{k!}\biggr) \biggr(1-\frac{\alpha}{n}\biggr)^{n} \biggr(1-\frac{\alpha}{n}\biggr)^{-k} \\&=\biggl(\frac{\alpha^k}{k!}\biggr) \lim \limits_{n \rightarrow \infty} \ \frac{n(n-1)(n-2) \cdots (n-k+1)}{n^k} \\&\times \ \ \ \lim \limits_{n \rightarrow \infty} \biggr(1-\frac{\alpha}{n}\biggr)^{n} \ \lim \limits_{n \rightarrow \infty} \biggr(1-\frac{\alpha}{n}\biggr)^{-k} \\&=\frac{e^{-\alpha} \alpha^k}{k!} \end{aligned}

In $(7)$, we have $\displaystyle \lim \limits_{n \rightarrow \infty} \ \frac{n(n-1)(n-2) \cdots (n-k+1)}{n^k}=1$. The reason being that the numerator is a polynomial where the leading term is $n^k$. Upon dividing by $n^k$ and taking the limit, we get 1. Based on $(5)$, we have $\displaystyle \lim \limits_{n \rightarrow \infty} \biggr(1-\frac{\alpha}{n}\biggr)^{n}=e^{-\alpha}$. For the last limit in the derivation we have $\displaystyle \lim \limits_{n \rightarrow \infty} \biggr(1-\frac{\alpha}{n}\biggr)^{-k}=1$.

We conclude with some comments. As the above derivation shows, the Poisson distribution is at heart a binomial distribution. When we divide the unit time interval into more and more subintervals (as the subintervals get more and more granular), the resulting binomial distribution behaves more and more like the Poisson distribution.

The three assumtions used in the derivation are called the Poisson postulates, which are the underlying assumptions that govern a Poisson process. Such a random process describes the occurrences of some type of events that are of interest (e.g. the arrivals of cars in our example) in a fixed period of time. The positive constant $\alpha$ indicated in Assumption 2 is the parameter of the Poisson process, which can be interpreted as the rate of occurrences of the event of interest (or rate of changes, or rate of arrivals) in a unit time interval, meaning that the positive constant $\alpha$ is the mean number of occurrences in the unit time interval. The derivation in $(7)$ shows that whenever a certain type of events occurs according to a Poisson process with parameter $\alpha$, the counting variable of the number of occurrences in the unit time interval is distributed according to the Poisson distribution as indicated in $(4)$.

If we observe the occurrences of events over intervals of length other than unit length, say, in an interval of length $t$, the counting process is governed by the same three postulates, with the modification to Assumption 2 that the rate of changes of the process is now $\alpha t$. The mean number of occurrences in the time interval of length $t$ is now $\alpha t$. The Assumption 2 now states that for any very short time interval of length $h$ (and that is also a subinterval of the interval of length $t$ under observation), the probability of having one occurrence of event in this short interval is $(\alpha t)h$. Applyng the same derivation, it can be shown that the number of occurrences ($X_t$) in a time interval of length $t$ has the Poisson distribution with the following probability mass function.

$\displaystyle (8) \ \ \ \ \ P(X_t=k)=\frac{e^{-\alpha t} \ (\alpha t)^k}{k!} \ \ \ \ \ \ \ \ k=0,1,2,\cdots$

# Relating Binomial and Negative Binomial

The negative binomial distribution has a natural intepretation as a waiting time until the arrival of the rth success (when the parameter r is a positive integer). The waiting time refers to the number of independent Bernoulli trials needed to reach the rth success. This interpretation of the negative binomial distribution gives us a good way of relating it to the binomial distribution. For example, if the rth success takes place after $k$ failed Bernoulli trials (for a total of $k+r$ trials), then there can be at most $r-1$ successes in the first $k+r$ trials. This tells us that the survival function of the negative binomial distribution is the cumulative distribution function (cdf) of a binomial distribution. In this post, we gives the details behind this observation. A previous post on the negative binomial distribution is found here.

A random experiment resulting in two distinct outcomes (success or failure) is called a Bernoulli trial (e.g. head or tail in a coin toss, whether or not the birthday of a customer is the first of January, whether an insurance claim is above or below a given threshold etc). Suppose a series of independent Bernoulli trials are performed until reaching the rth success where the probability of success in each trial is $p$. Let $X_r$ be the number of failures before the occurrence of the rth success. The following is the probablity mass function of $X_r$.

$\displaystyle (1) \ \ \ \ P(X_r=k)=\binom{k+r-1}{k} p^r (1-p)^k \ \ \ \ \ \ k=0,1,2,3,\cdots$

Be definition, the survival function and cdf of $X_r$ are:

$\displaystyle (2) \ \ \ \ P(X_r > k)=\sum \limits_{j=k+1}^\infty \binom{j+r-1}{j} p^r (1-p)^j \ \ \ \ \ \ k=0,1,2,3,\cdots$

$\displaystyle (3) \ \ \ \ P(X_r \le k)=\sum \limits_{j=0}^k \binom{j+r-1}{j} p^r (1-p)^j \ \ \ \ \ \ k=0,1,2,3,\cdots$

For each positive integer $k$, let $Y_{r+k}$ be the number of successes in performing a sequence of $r+k$ independent Bernoulli trials where $p$ is the probability of success. In other words, $Y_{r+k}$ has a binomial distribution with parameters $r+k$ and $p$.

If the random experiment requires more than $k$ failures to reach the rth success, there are at most $r-1$ successes in the first $k+r$ trails. Thus the survival function of $X_r$ is the same as the cdf of a binomial distribution. Equivalently, the cdf of $X_r$ is the same as the survival function of a binomial distribution. We have the following:

\displaystyle \begin{aligned}(4) \ \ \ \ P(X_r > k)&=P(Y_{k+r} \le r-1) \\&=\sum \limits_{j=0}^{r-1} \binom{k+r}{j} p^j (1-p)^{k+r-j} \ \ \ \ \ \ k=0,1,2,3,\cdots \end{aligned}

\displaystyle \begin{aligned}(5) \ \ \ \ P(X_r \le k)&=P(Y_{k+r} > r-1) \ \ \ \ \ \ k=0,1,2,3,\cdots \end{aligned}

Remark
The relation $(4)$ is analogous to the relationship between the Gamma distribution and the Poisson distribution. Recall that a Gamma distribution with shape parameter $\alpha$ and scale parameter $n$, where $n$ is a positive integer, can be interpreted as the waiting time until the nth change in a Poisson process. Thus, if the nth change takes place after time $t$, there can be at most $n-1$ arrivals in the time interval $[0,t]$. Thus the survival function of this Gamma distribution is the same as the cdf of a Poisson distribution. The relation $(4)$ is analogous to the following relation.

$\displaystyle (5) \ \ \ \ \int_t^\infty \frac{\alpha^n}{(n-1)!} \ x^{n-1} \ e^{-\alpha x} \ dx=\sum \limits_{j=0}^{n-1} \frac{e^{-\alpha t} \ (\alpha t)^j}{j!}$

A previous post on the negative binomial distribution is found here.

# The Negative Binomial Distribution

A counting distribution is a discrete distribution with probabilities only on the nonnegative integers. Such distributions are important in insurance applications since they can be used to model the number of events such as losses to the insured or claims to the insurer. Though playing a prominent role in statistical theory, the Poisson distribution is not appropriate in all situations, since it requires that the mean and the variance are equaled. Thus the negative binomial distribution is an excellent alternative to the Poisson distribution, especially in the cases where the observed variance is greater than the observed mean.

The negative binomial distribution arises naturally from a probability experiment of performing a series of independent Bernoulli trials until the occurrence of the rth success where r is a positive integer. From this starting point, we discuss three ways to define the distribution. We then discuss several basic properties of the negative binomial distribution. Emphasis is placed on the close connection between the Poisson distribution and the negative binomial distribution.

________________________________________________________________________

Definitions
We define three versions of the negative binomial distribution. The first two versions arise from the view point of performing a series of independent Bernoulli trials until the rth success where r is a positive integer. A Bernoulli trial is a probability experiment whose outcome is random such that there are two possible outcomes (success or failure).

Let $X_1$ be the number of Bernoulli trials required for the rth success to occur where r is a positive integer. Let $p$ is the probability of success in each trial. The following is the probability function of $X_1$:

$\displaystyle (1) \ \ \ \ \ P(X_1=x)= \binom{x-1}{r-1} p^r (1-p)^{x-r} \ \ \ \ \ \ \ x=r,r+1,r+2,\cdots$

The idea for $(1)$ is that for $X_1=x$ to happen, there must be $r-1$ successes in the first $x-1$ trials and one additional success occurring in the last trial (the $x$th trial).

A more common version of the negative binomial distribution is the number of Bernoulli trials in excess of r in order to produce the rth success. In other words, we consider the number of failures before the occurrence of the rth success. Let $X_2$ be this random variable. The following is the probability function of $X_2$:

$\displaystyle (2) \ \ \ \ \ P(X_2=x)=\binom{x+r-1}{x} p^r (1-p)^x \ \ \ \ \ \ \ x=0,1,2,\cdots$

The idea for $(2)$ is that there are $x+r$ trials and in the first $x+r-1$ trials, there are $x$ failures (or equivalently $r-1$ successes).

In both $(1)$ and $(2)$, the binomial coefficient is defined by

$\displaystyle (3) \ \ \ \ \ \binom{y}{k}=\frac{y!}{k! \ (y-k)!}=\frac{y(y-1) \cdots (y-(k-1))}{k!}$

where $y$ is a positive integer and $k$ is a nonnegative integer. However, the right-hand-side of $(3)$ can be calculated even if $y$ is not a positive integer. Thus the binomial coefficient $\displaystyle \binom{y}{k}$ can be expanded to work for all real number $y$. However $k$ must still be nonnegative integer.

$\displaystyle (4) \ \ \ \ \ \binom{y}{k}=\frac{y(y-1) \cdots (y-(k-1))}{k!}$

For convenience, we let $\displaystyle \binom{y}{0}=1$. When the real number $y>k-1$, the binomial coefficient in $(4)$ can be expressed as:

$\displaystyle (5) \ \ \ \ \ \binom{y}{k}=\frac{\Gamma(y+1)}{\Gamma(k+1) \Gamma(y-k+1)}$

where $\Gamma(\cdot)$ is the gamma function.

With the more relaxed notion of binomial coefficient, the probability function in $(2)$ above can be defined for all real number r. Thus the general version of the negative binomial distribution has two parameters r and $p$, both real numbers, such that $0. The following is its probability function.

$\displaystyle (6) \ \ \ \ \ P(X=x)=\binom{x+r-1}{x} p^r (1-p)^x \ \ \ \ \ \ \ x=0,1,2,\cdots$

Whenever r in $(6)$ is a real number that is not a positive integer, the interpretation of counting the number of failures until the occurrence of the rth success is no longer important. Instead we can think of it simply as a count distribution.

The following alternative parametrization of the negative binomial distribution is also useful.

$\displaystyle (6a) \ \ \ \ \ P(X=x)=\binom{x+r-1}{x} \biggl(\frac{\alpha}{\alpha+1}\biggr)^r \biggl(\frac{1}{\alpha+1}\biggr)^x \ \ \ \ \ \ \ x=0,1,2,\cdots$

The parameters in this alternative parametrization are r and $\alpha>0$. Clearly, the ratio $\frac{\alpha}{\alpha+1}$ takes the place of $p$ in $(6)$. Unless stated otherwise, we use the parametrization of $(6)$.
________________________________________________________________________

What is negative about the negative binomial distribution?
What is negative about this distribution? What is binomial about this distribution? The name is suggested by the fact that the binomial coefficient in $(6)$ can be rearranged as follows:

\displaystyle \begin{aligned}(7) \ \ \ \ \ \binom{x+r-1}{x}&=\frac{(x+r-1)(x+r-2) \cdots r}{x!} \\&=(-1)^x \frac{(-r-(x-1))(-r-(x-2)) \cdots (-r)}{x!} \\&=(-1)^x \frac{(-r)(-r-1) \cdots (-r-(x-1))}{x!} \\&=(-1)^x \binom{-r}{x} \end{aligned}

The calculation in $(7)$ can be used to verify that $(6)$ is indeed a probability function, that is, all the probabilities sum to 1.

\displaystyle \begin{aligned}(8) \ \ \ \ \ 1&=p^r p^{-r}\\&=p^r (1-q)^{-r} \\&=p^r \sum \limits_{x=0}^\infty \binom{-r}{x} (-q)^x \ \ \ \ \ \ \ \ (8.1) \\&=p^r \sum \limits_{x=0}^\infty (-1)^x \binom{-r}{x} q^x \\&=\sum \limits_{x=0}^\infty \binom{x+r-1}{x} p^r q^x \end{aligned}

In $(8)$, we take $q=1-p$. The step $(8.1)$ above uses the following formula known as the Newton’s binomial formula.

$\displaystyle (9) \ \ \ \ \ (1+t)^w=\sum \limits_{k=0}^\infty \binom{w}{k} t^k$

For a detailed discussion of (8) with all the details worked out, see the post called Deriving some facts of the negative binomial distribution.

________________________________________________________________________

The Generating Function
By definition, the following is the generating function of the negative binomial distribution, using :

$\displaystyle (10) \ \ \ \ \ g(z)=\sum \limits_{x=0}^\infty \binom{r+x-1}{x} p^r q^x z^x$

where $q=1-p$. Using a similar calculation as in $(8)$, the generating function can be simplified as:

$\displaystyle (11) \ \ \ \ \ g(z)=p^r (1-q z)^{-r}=\frac{p^r}{(1-q z)^r}=\frac{p^r}{(1-(1-p) z)^r}; \ \ \ \ \ z<\frac{1}{1-p}$

As a result, the moment generating function of the negative binomial distribution is:

$\displaystyle (12) \ \ \ \ \ M(t)=\frac{p^r}{(1-(1-p) e^t)^r}; \ \ \ \ \ \ \ t<-ln(1-p)$

For a detailed discussion of (12) with all the details worked out, see the post called Deriving some facts of the negative binomial distribution.

________________________________________________________________________

Independent Sum

One useful property of the negative binomial distribution is that the independent sum of negative binomial random variables, all with the same parameter $p$, also has a negative binomial distribution. Let $Y=Y_1+Y_2+\cdots+Y_n$ be an independent sum such that each $X_i$ has a negative binomial distribution with parameters $r_i$ and $p$. Then the sum $Y=Y_1+Y_2+\cdots+Y_n$ has a negative binomial distribution with parameters $r=r_1+\cdots+r_n$ and $p$.

Note that the generating function of an independent sum is the product of the individual generating functions. The following shows that the product of the individual generating functions is of the same form as $(11)$, thus proving the above assertion.

$\displaystyle (13) \ \ \ \ \ h(z)=\frac{p^{\sum \limits_{i=1}^n r_i}}{(1-(1-p) z)^{\sum \limits_{i=1}^n r_i}}$
________________________________________________________________________

Mean and Variance
The mean and variance can be obtained from the generating function. From $E(X)=g'(1)$ and $E(X^2)=g'(1)+g^{(2)}(1)$, we have:

$\displaystyle (14) \ \ \ \ \ E(X)=\frac{r(1-p)}{p} \ \ \ \ \ \ \ \ \ \ \ \ \ Var(X)=\frac{r(1-p)}{p^2}$

Note that $Var(X)=\frac{1}{p} E(X)>E(X)$. Thus when the sample data suggest that the variance is greater than the mean, the negative binomial distribution is an excellent alternative to the Poisson distribution. For example, suppose that the sample mean and the sample variance are 3.6 and 7.1. In exploring the possibility of fitting the data using the negative binomial distribution, we would be interested in the negative binomial distribution with this mean and variance. Then plugging these into $(14)$ produces the negative binomial distribution with $r=3.7$ and $p=0.507$.
________________________________________________________________________

The Poisson-Gamma Mixture
One important application of the negative binomial distribution is that it is a mixture of a family of Poisson distributions with Gamma mixing weights. Thus the negative binomial distribution can be viewed as a generalization of the Poisson distribution. The negative binomial distribution can be viewed as a Poisson distribution where the Poisson parameter is itself a random variable, distributed according to a Gamma distribution. Thus the negative binomial distribution is known as a Poisson-Gamma mixture.

In an insurance application, the negative binomial distribution can be used as a model for claim frequency when the risks are not homogeneous. Let $N$ has a Poisson distribution with parameter $\theta$, which can be interpreted as the number of claims in a fixed period of time from an insured in a large pool of insureds. There is uncertainty in the parameter $\theta$, reflecting the risk characteristic of the insured. Some insureds are poor risks (with large $\theta$) and some are good risks (with small $\theta$). Thus the parameter $\theta$ should be regarded as a random variable $\Theta$. The following is the conditional distribution of the random variable $N$ (conditional on $\Theta=\theta$):

$\displaystyle (15) \ \ \ \ \ P(N=n \lvert \Theta=\theta)=\frac{e^{-\theta} \ \theta^n}{n!} \ \ \ \ \ \ \ \ \ \ n=0,1,2,\cdots$

Suppose that $\Theta$ has a Gamma distribution with scale parameter $\alpha$ and shape parameter $\beta$. The following is the probability density function of $\Theta$.

$\displaystyle (16) \ \ \ \ \ g(\theta)=\frac{\alpha^\beta}{\Gamma(\beta)} \theta^{\beta-1} e^{-\alpha \theta} \ \ \ \ \ \ \ \ \ \ \theta>0$

Then the joint density of $N$ and $\Theta$ is:

$\displaystyle (17) \ \ \ \ \ P(N=n \lvert \Theta=\theta) \ g(\theta)=\frac{e^{-\theta} \ \theta^n}{n!} \ \frac{\alpha^\beta}{\Gamma(\beta)} \theta^{\beta-1} e^{-\alpha \theta}$

The unconditional distribution of $N$ is obtained by summing out $\theta$ in $(17)$.

\displaystyle \begin{aligned}(18) \ \ \ \ \ P(N=n)&=\int_0^\infty P(N=n \lvert \Theta=\theta) \ g(\theta) \ d \theta \\&=\int_0^\infty \frac{e^{-\theta} \ \theta^n}{n!} \ \frac{\alpha^\beta}{\Gamma(\beta)} \ \theta^{\beta-1} \ e^{-\alpha \theta} \ d \theta \\&=\int_0^\infty \frac{\alpha^\beta}{n! \ \Gamma(\beta)} \ \theta^{n+\beta-1} \ e^{-(\alpha+1) \theta} d \theta \\&=\frac{\alpha^\beta}{n! \ \Gamma(\beta)} \ \frac{\Gamma(n+\beta)}{(\alpha+1)^{n+\beta}} \int_0^\infty \frac{(\alpha+1)^{n+\beta}}{\Gamma(n+\beta)} \theta^{n+\beta-1} \ e^{-(\alpha+1) \theta} d \theta \\&=\frac{\alpha^\beta}{n! \ \Gamma(\beta)} \ \frac{\Gamma(n+\beta)}{(\alpha+1)^{n+\beta}} \\&=\frac{\Gamma(n+\beta)}{\Gamma(n+1) \ \Gamma(\beta)} \ \biggl( \frac{\alpha}{\alpha+1}\biggr)^\beta \ \biggl(\frac{1}{\alpha+1}\biggr)^n \\&=\binom{n+\beta-1}{n} \ \biggl( \frac{\alpha}{\alpha+1}\biggr)^\beta \ \biggl(\frac{1}{\alpha+1}\biggr)^n \ \ \ \ \ \ \ \ \ n=0,1,2,\cdots \end{aligned}

Note that the integral in the fourth step in $(18)$ is 1.0 since the integrand is the pdf of a Gamma distribution. The above probability function is that of a negative binomial distribution. It is of the same form as $(6a)$. Equivalently, it is also of the form $(6)$ with parameter $r=\beta$ and $p=\frac{\alpha}{\alpha+1}$.

The variance of the negative binomial distribution is greater than the mean. In a Poisson distribution, the mean equals the variance. Thus the unconditional distribution of $N$ is more dispersed than its conditional distributions. This is a characteristic of mixture distributions. The uncertainty in the parameter variable $\Theta$ has the effect of increasing the unconditional variance of the mixture distribution of $N$. The variance of a mixture distribution has two components, the weighted average of the conditional variances and the variance of the conditional means. The second component represents the additional variance introduced by the uncertainty in the parameter $\Theta$ (see The variance of a mixture).

________________________________________________________________________

The Poisson Distribution as Limit of Negative Binomial
There is another connection to the Poisson distribution, that is, the Poisson distribution is a limiting case of the negative binomial distribution. We show that the generating function of the Poisson distribution can be obtained by taking the limit of the negative binomial generating function as $r \rightarrow \infty$. Interestingly, the Poisson distribution is also the limit of the binomial distribution.

In this section, we use the negative binomial parametrization of $(6a)$. By replacing $\frac{\alpha}{\alpha+1}$ for $p$, the following are the mean, variance, and the generating function for the probability function in $(6a)$:

\displaystyle \begin{aligned}(19) \ \ \ \ \ \ &E(X)=\frac{r}{\alpha} \\&\text{ }\\&Var(X)=\frac{\alpha+1}{\alpha} \ \frac{r}{\alpha}=\frac{r(\alpha+1)}{\alpha^2} \\&\text{ } \\&g(z)=\frac{1}{[1-\frac{1}{\alpha}(z-1)]^r} \ \ \ \ \ \ \ z<\alpha+1 \end{aligned}

Let r goes to infinity and $\displaystyle \frac{1}{\alpha}$ goes to zero and at the same time keeping their product constant. Thus $\displaystyle \mu=\frac{r}{\alpha}$ is constant (this is the mean of the negative binomial distribution). We show the following:

$\displaystyle (20) \ \ \ \ \ \lim \limits_{r \rightarrow \infty} [1-\frac{\mu}{r}(z-1)]^{-r}=e^{\mu (z-1)}$

The right-hand side of $(20)$ is the generating function of the Poisson distribution with mean $\mu$. The generating function in the left-hand side is that of a negative binomial distribution with mean $\displaystyle \mu=\frac{r}{\alpha}$. The following is the derivation of $(20)$.

\displaystyle \begin{aligned}(21) \ \ \ \ \ \lim \limits_{r \rightarrow \infty} [1-\frac{\mu}{r}(z-1)]^{-r}&=\lim \limits_{r \rightarrow \infty} e^{\displaystyle \biggl(ln[1-\frac{\mu}{r}(z-1)]^{-r}\biggr)} \\&=\lim \limits_{r \rightarrow \infty} e^{\displaystyle \biggl(-r \ ln[1-\frac{\mu}{r}(z-1)]\biggr)} \\&=e^{\displaystyle \biggl(\lim \limits_{r \rightarrow \infty} -r \ ln[1-\frac{\mu}{r}(z-1)]\biggr)} \end{aligned}

We now focus on the limit in the exponent.

\displaystyle \begin{aligned}(22) \ \ \ \ \ \lim \limits_{r \rightarrow \infty} -r \ ln[1-\frac{\mu}{r}(z-1)]&=\lim \limits_{r \rightarrow \infty} \frac{ln(1-\frac{\mu}{r} (z-1))^{-1}}{r^{-1}} \\&=\lim \limits_{r \rightarrow \infty} \frac{(1-\frac{\mu}{r} (z-1)) \ \mu (z-1) r^{-2}}{r^{-2}} \\&=\mu (z-1) \end{aligned}

The middle step in $(22)$ uses the L’Hopital’s Rule. The result in $(20)$ is obtained by combining $(21)$ and $(22)$.

________________________________________________________________________

Reference

1. Klugman S.A., Panjer H. H., Wilmot G. E. Loss Models, From Data to Decisions, Second Edition., Wiley-Interscience, a John Wiley & Sons, Inc., New York, 2004

________________________________________________________________________
$\copyright \ \text{2011-2015 by Dan Ma}$

# Splitting a Poisson Distribution

We consider a remarkable property of the Poisson distribution that has a connection to the multinomial distribution. We start with the following examples.

Example 1
Suppose that the arrivals of customers in a gift shop at an airport follow a Poisson distribution with a mean of $\alpha=5$ per 10 minutes. Furthermore, suppose that each arrival can be classified into one of three distinct types – type 1 (no purchase), type 2 (purchase under $20), and type 3 (purchase over$20). Records show that about 25% of the customers are of type 1. The percentages of type 2 and type 3 are 60% and 15%, respectively. What is the probability distribution of the number of customers per hour of each type?

Example 2
Roll a fair die $N$ times where $N$ is random and follows a Poisson distribution with parameter $\alpha$. For each $i=1,2,3,4,5,6$, let $N_i$ be the number of times the upside of the die is $i$. What is the probability distribution of each $N_i$? What is the joint distribution of $N_1,N_2,N_3,N_4,N_5,N_6$?

In Example 1, the stream of customers arrive according to a Poisson distribution. It can be shown that the stream of each type of customers also has a Poisson distribution. One way to view this example is that we can split the Poisson distribution into three Poisson distributions.

Example 2 also describes a splitting process, i.e. splitting a Poisson variable into 6 different Poisson variables. We can also view Example 2 as a multinomial distribution where the number of trials is not fixed but is random and follows a Poisson distribution. If the number of rolls of the die is fixed in Example 2 (say 10), then each $N_i$ would be a binomial distribution. Yet, with the number of trials being Poisson, each $N_i$ has a Poisson distribution with mean $\displaystyle \frac{\alpha}{6}$. In this post, we describe this Poisson splitting process in terms of a “random” multinomial distribution (the view point of Example 2).

________________________________________________________________________

Suppose we have a multinomial experiment with parameters $N$, $r$, $p_1, \cdots, p_r$, where

• $N$ is the number of multinomial trials,
• $r$ is the number of distinct possible outcomes in each trial (type 1 through type $r$),
• the $p_i$ are the probabilities of the $r$ possible outcomes in each trial.

Suppose that $N$ follows a Poisson distribution with parameter $\alpha$. For each $i=1, \cdots, r$, let $N_i$ be the number of occurrences of the $i^{th}$ type of outcomes in the $N$ trials. Then $N_1,N_2,\cdots,N_r$ are mutually independent Poisson random variables with parameters $\alpha p_1,\alpha p_2,\cdots,\alpha p_r$, respectively.

The variables $N_1,N_2,\cdots,N_r$ have a multinomial distribution and their joint probability function is:

$\displaystyle (1) \ \ \ \ P(N_1=n_1,N_2=n_2,\cdots,N_r=n_r)=\frac{N!}{n_1! n_2! \cdots n_r!} \ p_1^{n_1} p_2^{n_2} \cdots p_r^{n_r}$

where $n_i$ are nonnegative integers such that $N=n_1+n_2+\cdots+n_r$.

Since the total number of multinomial trials $N$ is not fixed and is random, $(1)$ is not the end of the story. The probability in $(1)$ is only a conditional probability. The following is the joint probability function of $N_1,N_2,\cdots,N_r$:

$\displaystyle (2) \ \ \ \ P(N_1=n_1,N_2=n_2,\cdots,N_r=n_r)$

\displaystyle \begin{aligned}&=P(N_1=n_1,N_2=n_2,\cdots,N_r=n_r \lvert N=\sum \limits_{k=0}^r n_k) \\&\ \ \ \ \ \times P(N=\sum \limits_{k=0}^r n_k) \\&\text{ } \\&=\frac{(\sum \limits_{k=0}^r n_k)!}{n_1! \ n_2! \ \cdots \ n_r!} \ p_1^{n_1} \ p_2^{n_2} \ \cdots \ p_r^{n_r} \ \times \frac{e^{-\alpha} \alpha^{\sum \limits_{k=0}^r n_k}}{(\sum \limits_{k=0}^r n_k!)} \\&\text{ } \\&=\frac{e^{-\alpha p_1} \ (\alpha p_1)^{n_1}}{n_1!} \ \frac{e^{-\alpha p_2} \ (\alpha p_2)^{n_2}}{n_2!} \ \cdots \ \frac{e^{-\alpha p_r} \ (\alpha p_r)^{n_r}}{n_r!} \end{aligned}

To obtain the marginal probability function of $N_j$, $j=1,2,\cdots,r$, we sum out the other variables $N_k=n_k$ ($k \ne j$) in $(2)$ and obtain the following:

$\displaystyle (3) \ \ \ \ P(N_j=n_j)=\frac{e^{-\alpha p_j} \ (\alpha p_j)^{n_j}}{n_j!}$

Thus we can conclude that $N_j$, $j=1,2,\cdots,r$, has a Poisson distribution with parameter $\alpha p_j$. Furrthermore, the joint probability function of $N_1,N_2,\cdots,N_r$ is the product of the marginal probability functions. Thus we can conclude that $N_1,N_2,\cdots,N_r$ are mutually independent.

________________________________________________________________________
Example 1
Let $N_1,N_2,N_3$ be the number of customers per hour of type 1, type 2, and type 3, respectively. Here, we attempt to split a Poisson distribution with mean 30 per hour (based on 5 per 10 minutes). Thus $N_1,N_2,N_3$ are mutually independent Poisson variables with means $30 \times 0.25=7.5$, $30 \times 0.60=18$, $30 \times 0.15=4.5$, respectively.

Example 2
As indicated earlier, each $N_i$, $i=1,2,3,4,5,6$, has a Poisson distribution with mean $\frac{\alpha}{6}$. According to $(2)$, the joint probability function of $N_1,N_2,N_3,N_4,N_5,N_6$ is simply the product of the six marginal Poisson probability functions.

# The Poisson Distribution

Let $\alpha$ be a positive constant. Consider the following probability distribution:

$\displaystyle (1) \ \ \ \ \ P(X=j)=\frac{e^{-\alpha} \alpha^j}{j!} \ \ \ \ \ j=0,1,2,\cdots$

The above distribution is said to be a Poisson distribution with parameter $\alpha$. The Poisson distribution is usually used to model the random number of events occurring in a fixed time interval. As will be shown below, $E(X)=\alpha$. Thus the parameter $\alpha$ is the rate of occurrence of the random events; it indicates on average how many events occur per unit of time. Examples of random events that may be modeled by the Poisson distribution include the number of alpha particles emitted by a radioactive substance counted in a prescribed area during a fixed period of time, the number of auto accidents in a fixed period of time or the number of losses arising from a group of insureds during a policy period.

Each of the above examples can be thought of as a process that generates a number of arrivals or changes in a fixed period of time. If such a counting process leads to a Poisson distribution, then the process is said to be a Poisson process.

We now discuss some basic properties of the Poisson distribution. Using the Taylor series expansion of $e^{\alpha}$, the following shows that $(1)$ is indeed a probability distribution.

$\displaystyle . \ \ \ \ \ \ \ \sum \limits_{j=0}^\infty \frac{e^{-\alpha} \alpha^j}{j!}=e^{-\alpha} \sum \limits_{j=0}^\infty \frac{\alpha^j}{j!}=e^{-\alpha} e^{\alpha}=1$

The generating function of the Poisson distribution is $g(z)=e^{\alpha (z-1)}$ (see The generating function). The mean and variance can be calculated using the generating function.

\displaystyle \begin{aligned}(2) \ \ \ \ \ &E(X)=g'(1)=\alpha \\&\text{ } \\&E[X(X-1)]=g^{(2)}(1)=\alpha^2 \\&\text{ } \\&Var(X)=E[X(X-1)]+E(X)-E(X)^2=\alpha^2+\alpha-\alpha^2=\alpha \end{aligned}

The Poisson distribution can also be interpreted as an approximation to the binomial distribution. It is well known that the Poisson distribution is the limiting case of binomial distributions (see [1] or this post).

$\displaystyle (3) \ \ \ \ \ \lim \limits_{n \rightarrow \infty} \binom{n}{j} \biggl(\frac{\alpha}{n}\biggr)^j \biggl(1-\frac{\alpha}{n}\biggr)^{n-j}=\frac{e^{-\alpha} \alpha^j}{j!}$

One application of $(3)$ is that we can use Poisson probabilities to approximate Binomial probabilities. The approximation is reasonably good when the number of trials $n$ in a binomial distribution is large and the probability of success $p$ is small. The binomial mean is $n p$ and the variance is $n p (1-p)$. When $p$ is small, $1-p$ is close to 1 and the binomial variance is approximately $np \approx n p (1-p)$. Whenever the mean of a discrete distribution is approximately equaled to the mean, the Poisson approximation is quite good. As a rule of thumb, we can use Poisson to approximate binomial if $n \le 100$ and $p \le 0.01$.

As an example, we use the Poisson distribution to estimate the probability that at most 1 person out of 1000 will have a birthday on the New Year Day. Let $n=1000$ and $p=365^{-1}$. So we use the Poisson distribution with $\alpha=1000 365^{-1}$. The following is an estimate using the Poisson distribution.

$\displaystyle . \ \ \ \ \ \ \ P(X \le 1)=e^{-\alpha}+\alpha e^{-\alpha}=(1+\alpha) e^{-\alpha}=0.2415$

Another useful property is that the independent sum of Poisson distributions also has a Poisson distribution. Specifically, if each $X_i$ has a Poisson distribution with parameter $\alpha_i$, then the independent sum $X=X_1+\cdots+X_n$ has a Poisson distribution with parameter $\alpha=\alpha_1+\cdots+\alpha_n$. One way to see this is that the product of Poisson generating functions has the same general form as $g(z)=e^{\alpha (z-1)}$ (see The generating function). One interpretation of this property is that when merging several arrival processes, each of which follow a Poisson distribution, the result is still a Poisson distribution.

For example, suppose that in an airline ticket counter, the arrival of first class customers follows a Poisson process with a mean arrival rate of 8 per 15 minutes and the arrival of customers flying coach follows a Poisson distribution with a mean rate of 12 per 15 minutes. Then the arrival of customers of either types has a Poisson distribution with a mean rate of 20 per 15 minutes or 80 per hour.

A Poisson distribution with a large mean can be thought of as an independent sum of Poisson distributions. For example, a Poisson distribution with a mean of 50 is the independent sum of 50 Poisson distributions each with mean 1. Because of the central limit theorem, when the mean is large, we can approximate the Poisson using the normal distribution.

In addition to merging several Poisson distributions into one combined Poisson distribution, we can also split a Poisson into several Poisson distributions. For example, suppose that a stream of customers arrives according to a Poisson distribution with parameter $\alpha$ and each customer can be classified into one of two types (e.g. no purchase vs. purchase) with probabilities $p_1$ and $p_2$, respectively. Then the number of “no purchase” customers and the number of “purchase” customers are independent Poisson random variables with parameters $\alpha p_1$ and $\alpha p_2$, respectively. For more details on the splitting of Poisson, see Splitting a Poisson Distribution.

Reference

1. Feller W. An Introduction to Probability Theory and Its Applications, Third Edition, John Wiley & Sons, New York, 1968

# The generating function

Consider the function $g(z)=\displaystyle e^{\alpha (z-1)}$ where $\alpha$ is a positive constant. The following shows the derivatives of this function.

\displaystyle \begin{aligned}. \ \ \ \ \ \ &g(z)=e^{\alpha (z-1)} \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ g(0)=e^{-\alpha} \\&\text{ } \\&g'(z)=e^{\alpha (z-1)} \ \alpha \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ g'(0)=e^{-\alpha} \ \alpha \\&\text{ } \\&g^{(2)}(z)=e^{\alpha (z-1)} \ \alpha^2 \ \ \ \ \ \ \ \ \ \ \ \ \ \ g^{(2)}(0)=2! \ \frac{e^{-\alpha} \ \alpha^2}{2!} \\&\text{ } \\&g^{(3)}(z)=e^{\alpha (z-1)} \ \alpha^3 \ \ \ \ \ \ \ \ \ \ \ \ \ \ g^{(3)}(0)=3! \ \frac{e^{-\alpha} \ \alpha^3}{3!} \\&\text{ } \\&\ \ \ \ \ \ \ \ \cdots \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \cdots \\&\text{ } \\&g^{(n)}(z)=e^{\alpha (z-1)} \ \alpha^n \ \ \ \ \ \ \ \ \ \ \ \ \ \ g^{(n)}(0)=n! \ \frac{e^{-\alpha} \ \alpha^n}{n!} \end{aligned}

Note that the derivative of $g(z)$ at each order is a multiple of a Poisson probability. Thus the Poisson distribution is coded by the function $g(z)=\displaystyle e^{\alpha (z-1)}$. Because of this reason, such a function is called a generating function (or probability generating function). This post discusses some basic facts about the generating function (gf) and its cousin, the moment generating function (mgf). One important characteristic is that these functions generate probabilities and moments. Another important characteristic is that there is a one-to-one correspondence between a probability distribution and its generating function and moment generating function, i.e. two random variables with different cumulative distribution functions cannot have the same gf or mgf. In some situations, this fact is useful in working with independent sum of random variables.

________________________________________________
The Generating Function
Suppose that $X$ is a random variable that takes only nonegative integer values with the probability function given by

$(1) \ \ \ \ \ \ P(X=j)=a_j, \ \ \ \ j=0,1,2,\cdots$

The idea of the generating function is that we use a power series to capture the entire probability distribution. The following defines the generating function that is associated with the above sequence $a_j$, .

$(2) \ \ \ \ \ \ g(z)=a_0+a_1 \ z+a_2 \ z^2+ \cdots=\sum \limits_{j=0}^\infty a_j \ z^j$

Since the elements of the sequence $a_j$ are probabilities, we can also call $g(z)$ the generating function of the probability distribution defined by the sequence in $(1)$. The generating function $g(z)$ is defined wherever the power series converges. It is clear that at the minimum, the power series in $(2)$ converges for $\lvert z \lvert \le 1$.

We discuss the following three properties of generating functions:

1. The generating function completely determines the distribution.
2. The moments of the distribution can be derived from the derivatives of the generating function.
3. The generating function of a sum of independent random variables is the product of the individual generating functions.

The Poisson generating function at the beginning of the post is an example demonstrating property 1 (see Example 0 below for the derivation of the generating function). In some cases, the probability distribution of an independent sum can be deduced from the product of the individual generating functions. Some examples are given below.

________________________________________________
Generating Probabilities
We now discuss the property 1 indicated above. To see that $g(z)$ generates the probabilities, let’s look at the derivatives of $g(z)$:

\displaystyle \begin{aligned}(3) \ \ \ \ \ \ &g'(z)=a_1+2 a_2 \ z+3 a_3 \ z^2+\cdots=\sum \limits_{j=1}^\infty j a_j \ z^{j-1} \\&\text{ } \\&g^{(2)}(z)=2 a_2+6 a_3 \ z+ 12 a_4 \ z^2=\sum \limits_{j=2}^\infty j (j-1) a_j \ z^{j-2} \\&\text{ } \\&g^{(3)}(z)=6 a_3+ 24 a_4 \ z+60 a_5 \ z^2=\sum \limits_{j=3}^\infty j (j-1)(j-2) a_j \ z^{j-3} \\&\text{ } \\&\ \ \ \ \ \ \ \ \cdots \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \cdots \\&\text{ } \\&g^{(n)}(z)=\sum \limits_{j=n}^\infty j(j-1) \cdots (j-n+1) a_j \ z^{j-n}=\sum \limits_{j=n}^\infty \binom{j}{n} n! \ a_j \ z^{j-n} \end{aligned}

By letting $z=0$ above, all the terms vanishes except for the constant term. We have:

$(4) \ \ \ \ \ \ g^{(n)}(0)=n! \ a_n=n! \ P(X=n)$

Thus the generating function is a compact way of encoding the probability distribution. The probability distribution determines the generating function as seen in $(2)$. On the other hand, $(3)$ and $(4)$ demonstrate that the generating function also determines the probability distribution.

________________________________________________
Generating Moments
The generating function also determines the moments (property 2 indicated above). For example, we have:

\displaystyle \begin{aligned}(5) \ \ \ \ \ \ &g'(1)=0 \ a_0+a_1+2 a_2+3 a_3+\cdots=\sum \limits_{j=0}^\infty j a_j=E(X) \\&\text{ } \\&g^{(2)}(1)=0 a_0 + 0 a_1+2 a_2+6 a_3+ 12 a_4+\cdots=\sum \limits_{j=0}^\infty j (j-1) a_j=E[X(X-1)] \\&\text{ } \\&E(X)=g'(1) \\&\text{ } \\&E(X^2)=g'(1)+g^{(2)}(1) \end{aligned}

Note that $g^{(n)}(1)=E[X(X-1) \cdots (X-(n-1))]$. Thus the higher moment $E(X^n)$ can be expressed in terms of $g^{(n)}(1)$ and $g^{(k)}(1)$ where $k.
________________________________________________
More General Definitions
Note that the definition in $(2)$ can also be interpreted as the mathematical expectation of $z^X$, i.e., $g(z)=E(z^X)$. This provides a way to define the generating function for random variables that may take on values outside of the nonnegative integers. The following is a more general definition of the generating function of the random variable $X$, which is defined for all $z$ where the expectation exists.

$(6) \ \ \ \ \ \ g(z)=E(z^X)$

________________________________________________
The Generating Function of Independent Sum
Let $X_1,X_2,\cdots,X_n$ be independent random variables with generating functions $g_1,g_2,\cdots,g_n$, respectively. Then the generating function of $X_1+X_2+\cdots+X_n$ is given by the product $g_1 \cdot g_2 \cdots g_n$.

Let $g(z)$ be the generating function of the independent sum $X_1+X_2+\cdots+X_n$. The following derives $g(z)$. Note that the general form of generating function $(6)$ is used.

\displaystyle \begin{aligned}(7) \ \ \ \ \ \ g(z)&=E(z^{X_1+\cdots+X_n}) \\&\text{ } \\&=E(z^{X_1} \cdots z^{X_n}) \\&\text{ } \\&=E(z^{X_1}) \cdots E(z^{X_n}) \\&\text{ } \\&=g_1(z) \cdots g_n(z) \end{aligned}

The probability distribution of a random variable is uniquely determined by its generating function. In particular, the generating function $g(z)$ of the independent sum $X_1+X_2+\cdots+X_n$ that is derived in $(7)$ is unique. So if the generating function is of a particular distribution, we can deduce that the distribution of the sum must be of the same distribution. See the examples below.

________________________________________________
Example 0
In this example, we derive the generating function of the Poisson distribution. Based on the definition, we have:

\displaystyle \begin{aligned}. \ \ \ \ \ \ g(z)&=\sum \limits_{j=0}^\infty \frac{e^{-\alpha} \alpha^j}{j!} \ z^j \\&\text{ } \\&=\sum \limits_{j=0}^\infty \frac{e^{-\alpha} (\alpha z)^j}{j!} \\&\text{ } \\&=\frac{e^{-\alpha}}{e^{- \alpha z}} \sum \limits_{j=0}^\infty \frac{e^{-\alpha z} (\alpha z)^j}{j!} \\&\text{ } \\&=e^{\alpha (z-1)} \end{aligned}

Example 1
Suppose that $X_1,X_2,\cdots,X_n$ are independent random variables where each $X_i$ has a Bernoulli distribution with probability of success $p$. Let $q=1-p$. The following is the generating function for each $X_i$.

$. \ \ \ \ \ \ g(z)=q+p z$

Then the generating function of the sum $X=X_1+\cdots+X_n$ is $g(z)^n=(q+p z)^n$. The following is the binomial expansion:

\displaystyle \begin{aligned}(8) \ \ \ \ \ \ g(z)^n&=(q+p z)^n \\&\text{ } \\&=\sum \limits_{j=0}^n \binom{n}{j} q^{n-j} \ p^j \ z^j \end{aligned}

By definition $(2)$, the generating function of $X=X_1+\cdots+X_n$ is:

$(9) \ \ \ \ \ \ g(z)^n=\sum \limits_{j=0}^\infty P(X=j) \ z^j$

Comparing $(8)$ and $(9)$, we have

$\displaystyle (10) \ \ \ \ \ \ P(X=j)=\left\{\begin{matrix}\displaystyle \binom{n}{j} p^j \ q^{n-j}&\ 0 \le j \le n\\{0}&\ j>n \end{matrix}\right.$

The probability distribution indicated by $(8)$ and $(10)$ is that of a binomial distribution. Since the probability distribution of a random variable is uniquely determined by its generating function, the independent sum of Bernoulli distributions must ave a Binomial distribution.

Example 2
Suppose that $X_1,X_2,\cdots,X_n$ are independent and have Poisson distributions with parameters $\alpha_1,\alpha_2,\cdots,\alpha_n$, respectively. Then the independent sum $X=X_1+\cdots+X_n$ has a Poisson distribution with parameter $\alpha=\alpha_1+\cdots+\alpha_n$.

Let $g(z)$ be the generating function of $X=X_1+\cdots+X_n$. For each $i$, the generating function of $X_i$ is $g_i(z)=e^{\alpha_i (z-1)}$. The key to the proof is that the product of the $g_i$ has the same general form as the individual $g_i$.

\displaystyle \begin{aligned}(11) \ \ \ \ \ \ g(z)&=g_1(z) \cdots g_n(z) \\&\text{ } \\&=e^{\alpha_1 (z-1)} \cdots e^{\alpha_n (z-1)} \\&\text{ } \\&=e^{(\alpha_1+\cdots+\alpha_n)(z-1)} \end{aligned}

The generating function in $(11)$ is that of a Poisson distribution with mean $\alpha=\alpha_1+\cdots+\alpha_n$. Since the generating function uniquely determines the distribution, we can deduce that the sum $X=X_1+\cdots+X_n$ has a Poisson distribution with parameter $\alpha=\alpha_1+\cdots+\alpha_n$.

Example 3
In rolling a fair die, let $X$ be the number shown on the up face. The associated generating function is:

$\displaystyle. \ \ \ \ \ \ g(z)=\frac{1}{6}(z+z^2+z^3+z^4+z^5+z^6)=\frac{z(1-z^6)}{6(1-z)}$

The generating function can be further reduced as:

\displaystyle \begin{aligned}. \ \ \ \ \ \ g(z)&=\frac{z(1-z^6)}{6(1-z)} \\&\text{ } \\&=\frac{z(1-z^3)(1+z^3)}{6(1-z)} \\&\text{ } \\&=\frac{z(1-z)(1+z+z^2)(1+z^3)}{6(1-z)} \\&\text{ } \\&=\frac{z(1+z+z^2)(1+z^3)}{6} \end{aligned}

Suppose that we roll the fair dice 4 times. Let $W$ be the sum of the 4 rolls. Then the generating function of $Z$ is

$\displaystyle. \ \ \ \ \ \ g(z)^4=\frac{z^4 (1+z^3)^4 (1+z+z^2)^4}{6^4}$

The random variable $W$ ranges from 4 to 24. Thus the probability function ranges from $P(W=4)$ to $P(W=24)$. To find these probabilities, we simply need to decode the generating function $g(z)^4$. For example, to find $P(W=12)$, we need to find the coefficient of the term $z^{12}$ in the polynomial $g(z)^4$. To help this decoding, we can expand two of the polynomials in $g(z)^4$.

\displaystyle \begin{aligned}. \ \ \ \ \ \ g(z)^4&=\frac{z^4 (1+z^3)^4 (1+z+z^2)^4}{6^4} \\&\text{ } \\&=\frac{z^4 \times A \times B}{6^4} \\&\text{ } \\&A=(1+z^3)^4=1+4z^3+6z^6+4z^9+z^{12} \\&\text{ } \\&B=(1+z+z^2)^4=1+4z+10z^2+16z^3+19z^4+16z^5+10z^6+4z^7+z^8 \end{aligned}

Based on the above polynomials, there are three ways of forming $z^{12}$. They are: $(z^4 \times 1 \times z^8)$, $(z^4 \times 4z^3 \times 16z^5)$, $(z^4 \times 6z^6 \times 10z^2)$. Thus we have:

$\displaystyle. \ \ \ \ \ \ P(W=12)=\frac{1}{6^4}(1+4 \times 16+6 \times 10)=\frac{125}{6^4}$

To find the other probabilities, we can follow the same decoding process.

________________________________________________
Remark
The probability distribution of a random variable is uniquely determined by its generating function. This fundamental property is useful in determining the distribution of an independent sum. The generating function of the independent sum is simply the product of the individual generating functions. If the product is of a certain distributional form (as in Example 1 and Example 2), then we can deduce that the sum must be of the same distribution.

We can also decode the product of generating functions to obtain the probability function of the independent sum (as in Example 3). The method in Example 3 is quite tedious. But one advantage is that it is a “machine process”, a pretty fool proof process that can be performed mechanically.

The machine process is this: Code the individual probability distribution in a generating function $g(z)$. Then raise it to $n$. After performing some manipulation to $g(z)^n$, decode the probabilities from $g(z)^n$.

As long as we can perform the algebraic manipulation carefully and correctly, this process will be sure to provide the probability distribution of an independent sum.

________________________________________________
The Moment Generating Function
The moment generating function of a random variable $X$ is $M_X(t)=E(e^{tX})$ on all real numbers $t$ for which the expected value exists. The moments can be computed more directly using an mgf. From the theory of mathematical analysis, it can be shown that if $M_X(t)$ exists on some interval $-a, then the derivatives of $M_X(t)$ of all orders exist at $t=0$. Furthermore, it can be show that $E(X^n)=M_X^{(n)}(0)$.

Suppose that $g(z)$ is the generating function of a random variable. The following relates the generating function and the moment generating function.

\displaystyle \begin{aligned}. \ \ \ \ \ \ &M_X(t)=g(e^t) \\&\text{ } \\&g(z)=M_X(ln z) \end{aligned}

________________________________________________

Reference

1. Feller W. An Introduction to Probability Theory and Its Applications, Third Edition, John Wiley & Sons, New York, 1968