Tutorials › Exam ASTAM - Advanced Short-Term Actuarial Mathematics › Poisson Sampling By Recursive Inversion

Home  /  Tutorials  /  Exam ASTAM / Random Variate Generation

Exam ASTAM · Random Variate Generation · Tutorial 594 of 979

Poisson Sampling By Recursive Inversion

Use recursive Poisson probabilities and a precise cumulative stopping rule to turn a uniform draw into a simulated count.

Beginner 12 min read

What You'll Learn

  • Initialize the Poisson probability recursion from the probability of zero.
  • Update each probability mass without repeatedly evaluating factorials.
  • Accumulate masses and apply the first-crossing stopping rule.
  • Handle cumulative-probability endpoints consistently with the comparison rule.
  • Check a simulated count and identify common recurrence and rounding errors.

Why Poisson Sampling Uses a Recursion

A Poisson count is used to represent the number of events in a specified period or exposure when the count model has a fixed mean and independent event occurrence in disjoint portions of exposure. In actuarial work, such a count might represent hypothetical claims during a policy period. To simulate it, we need a method that assigns the correct probability to every possible count.

The direct probability formula contains a factorial for each possible count. Recursive inversion avoids repeatedly calculating factorials: calculate the probability at zero, obtain each next probability from the one before it, and keep a running cumulative total. Stop as soon as that total reaches the uniform draw.

Definition: Let \(N\sim\operatorname{Poisson}(\lambda)\), where \(\lambda>0\) is the expected count over the stated exposure. Its probability mass function is \(P(N=k)=e^{-\lambda}\lambda^k/k!\), for \(k=0,1,2,\ldots\). Define \(p_k=P(N=k)\) and \(C_k=\sum_{j=0}^{k}p_j\).

Generate \(U\sim\operatorname{Uniform}[0,1)\). Starting at \(k=0\), calculate the Poisson masses and cumulative probabilities in order. Return the first \(k\) for which \(C_k\ge U\). The value of \(\lambda\) must correspond to the exposure being simulated; for example, a mean count per year cannot be used unchanged for a different exposure period unless the model specifies that mean for the new period.

The Recursive Formula and Stopping Rule

The initial mass is \(p_0=e^{-\lambda}\). For \(k\ge1\), divide the mass at \(k\) by the mass at \(k-1\):

$$ \frac{p_k}{p_{k-1}} = \frac{e^{-\lambda}\lambda^k/k!}{e^{-\lambda}\lambda^{k-1}/(k-1)!} = \frac{\lambda}{k}, \qquad p_k=\frac{\lambda}{k}p_{k-1}. $$

Thus, after initializing \(p_0\), each new mass requires one multiplication and one division. The running cumulative probability is updated by \(C_k=C_{k-1}+p_k\), with \(C_0=p_0\). In an implementation, keep the current mass and cumulative total at full available precision; do not round intermediate values before making the stopping comparison.

Stopping rule: With the convention “return the first \(k\) such that \(C_k\ge U\),” return 0 when \(0\le U\le C_0\), and return \(k\ge1\) when \(C_{k-1}<U\le C_k\). The endpoint convention agrees with the stated comparison rule. Since the Poisson probabilities sum to 1, \(C_k\) tends to 1 as \(k\) increases, so every \(U\in[0,1)\) is covered after finitely many steps.
1
Set the parameter and draw.
Specify the Poisson mean \(\lambda>0\) for the exposure, and generate one \(U\in[0,1)\).
2
Initialize at zero.
Set \(k=0\), \(p_0=e^{-\lambda}\), and \(C_0=p_0\). If \(C_0\ge U\), return 0.
3
Advance the mass and total.
Increase \(k\) by 1, calculate \(p_k=(\lambda/k)p_{k-1}\), and set \(C_k=C_{k-1}+p_k\).
4
Check and stop.
If \(C_k\ge U\), return \(k\). Otherwise repeat Step 3.

The reason this procedure produces the desired distribution is the length of each interval assigned to a count. Count 0 is returned on an interval of length \(C_0=p_0\). For \(k\ge1\), count \(k\) is returned when \(C_{k-1}<U\le C_k\), an interval of length \(C_k-C_{k-1}=p_k\). Because \(U\) is uniform, the probability of falling in that interval is \(p_k\), exactly the required Poisson probability.

Worked Applications

Worked Example: Finding a Count for a Uniform Draw

A hypothetical insurer models the number of small equipment incidents in a defined service period as \(N\sim\operatorname{Poisson}(2)\). A uniform draw is \(U=0.70\). Use recursive inversion to generate the count. Report intermediate probabilities to seven decimal places, but make the stopping decision using the unrounded values.

Method. Initialize with \(p_0=e^{-2}=0.1353352832\), so \(C_0=0.1353352832\). Since \(C_0<0.70\), continue.

For \(k=1\), \(p_1=(2/1)p_0=0.2706705665\), and \(C_1=p_0+p_1=0.4060058497\). Since \(C_1<0.70\), continue.

For \(k=2\), \(p_2=(2/2)p_1=0.2706705665\), and \(C_2=C_1+p_2=0.6766764162\). Since \(C_2<0.70\), continue.

For \(k=3\), \(p_3=(2/3)p_2=0.1804470443\), and \(C_3=C_2+p_3=0.8571234605\). Now \(C_3\ge0.70\), so the generated count is 3.

Check. The draw lies between \(C_2=0.6766764162\) and \(C_3=0.8571234605\), so the first-crossing rule returns 3. The calculation avoids computing separate factorials for each mass.

Worked Example: Checking the First Crossing Near a Cumulative Boundary

A hypothetical monitoring program models its event count as \(N\sim\operatorname{Poisson}(4)\). For \(U=0.30\), find the generated count. Show the masses and cumulative probabilities through the first crossing, reporting cumulative probabilities to ten decimal places.

The starting mass is \(p_0=e^{-4}=0.0183156389\), so \(C_0=0.0183156389\). This is below 0.30.

Next, \(p_1=(4/1)p_0=0.0732625556\), giving \(C_1=0.0915781944\). Then \(p_2=(4/2)p_1=0.1465251111\), giving \(C_2=0.2381033056\). This cumulative probability is still below 0.30.

For \(k=3\), \(p_3=(4/3)p_2=0.1953668148\), and \(C_3=C_2+p_3=0.4334701204\). This is the first cumulative probability at least 0.30, so return 3.

Check. The unrounded cumulative values satisfy \(C_2<0.30\le C_3\), so the count is 3. The displayed values are rounded to ten decimal places; the crossing decision is consistent with the unrounded values.

Worked Example: Continuing Until a Later Crossing

For a hypothetical portfolio segment, suppose the claim count over the stated period is \(N\sim\operatorname{Poisson}(3)\). A uniform draw is \(U=0.90\). Use recursive inversion and report cumulative probabilities to ten decimal places.

Start with \(p_0=e^{-3}=0.0497870684\), so \(C_0=0.0497870684\). The next masses and cumulative values follow from \(p_k=(3/k)p_{k-1}\):

$$ \begin{aligned} p_1&=0.1493612051, & C_1&=0.1991482735,\\ p_2&=0.2240418077, & C_2&=0.4231900811,\\ p_3&=0.2240418077, & C_3&=0.6472318888,\\ p_4&=0.1680313557, & C_4&=0.8152632445,\\ p_5&=0.1008188134, & C_5&=0.9160820580. \end{aligned} $$

The cumulative values through \(k=4\) are all below 0.90, while \(C_5=0.9160820580\) exceeds it. Therefore, the generated count is 5.

Check. The first-crossing condition is \(C_4<0.90\le C_5\), which holds. This confirms that the algorithm does not stop merely because the cumulative probability is large; it stops at the first crossing of the actual draw.

Common Mistakes and Numerical Checks

A frequent error is to use \(p_k=\lambda p_{k-1}\) instead of \(p_k=(\lambda/k)p_{k-1}\). The missing division by \(k\) changes the probabilities and does not generate a Poisson distribution. For example, with \(\lambda=2\), the correct values are \(p_1=2p_0\) and \(p_2=(2/2)p_1=p_1\). The incorrect recurrence would give \(p_2=2p_1\), which overstates the mass at 2. The correction is to include the current count \(k\) in the denominator at every step.

Another trap is mixing the stopping comparison with its interval endpoints. If the algorithm returns the first \(k\) with \(C_k\ge U\), then \(U=C_k\) returns \(k\), while \(U=C_{k-1}\) has already caused a return at \(k-1\). Therefore, the interval for \(k\ge1\) is \(C_{k-1}<U\le C_k\), not \(C_{k-1}\le U<C_k\). Endpoint events have probability zero under a continuous uniform draw, but consistent endpoints matter when explaining or implementing the stated algorithm.

Rounding can also create a false crossing when a draw is close to a cumulative boundary. Keep unrounded values for the comparison whenever possible. If an implementation has limited precision and a very large \(\lambda\), direct initialization with \(e^{-\lambda}\) can underflow numerically; a more stable computational method may then be required. The simple recursion is appropriate when its intermediate values can be represented accurately.

Key takeaway: For \(\lambda>0\), start with \(p_0=e^{-\lambda}\), update with \(p_k=(\lambda/k)p_{k-1}\), and accumulate \(C_k\) until the first \(C_k\ge U\). If \(\lambda=0\), the Poisson count is identically zero, so return 0 without using the recursion.

Check Your Understanding

For each problem, use the first-crossing rule and the stated uniform draw. Treat displayed parameter and draw values as exact.

  1. A hypothetical event count has \(N\sim\operatorname{Poisson}(1)\), and \(U=0.20\). Which count is returned? (Answer: 0)
  2. A hypothetical service count has \(N\sim\operatorname{Poisson}(2)\), and \(U=0.50\). Which count is returned? (Answer: 2)
  3. A hypothetical count has \(N\sim\operatorname{Poisson}(4)\), and \(U=0.24\). Which count is returned? (Answer: 3)
  4. A hypothetical annual count has \(N\sim\operatorname{Poisson}(3)\), and \(U=0.80\). Which count is returned? (Answer: 4)
  5. A hypothetical count over a specified exposure has \(N\sim\operatorname{Poisson}(5)\), and \(U=0.95\). Which count is returned? (Answer: 9)

Previous

Next

Back to Tutorials
Clinical Biostats