Home About
Courses Concepts Tools
GitHub

Simulating a Poisson Distribution

Simulating a Poisson distribution (starting with a uniform RV)

This is a simple example of simulating a Poisson distribution.

import numpy as np
from matplotlib import pyplot as plt

Recall that a Poisson distribution is defined as having exponential interarrival times. This is useful in simulating a Poisson distribution, as we can easily get an exponential RV from uniforms with inversion.

def gen_exp(size=1):
    uniforms = np.random.uniform(0,1,size=size)
    return -np.log(uniforms)
plt.hist(gen_exp(1000), bins=30)

Now that we have our exponential RV we can use it to generate interarrival times for our Poisson RV.

def poisson(time):
    arrivals = 0
    total = 0
    while total < time:
        total += gen_exp()
        arrivals += 1
    return arrivals - 1

The below histograms illustrate that our Poisson RV does indeed behave as expected, given various expected numbers of arrivals \lambda

# expected Poisson values, for use with overlaying onto histograms
def pmf(lam, k):
    return (lam**k * np.exp(-lam))/np.math.factorial(k)
sample = [poisson(50) for _ in range(10_000)]

exp_space = np.linspace(25, 80, 100)
expected = pmf(exp_space, 50)
plt.plot(exp_space, expected)
plt.title(r'$\lambda = {}$'.format(50))
plt.hist(sample, density=True, bins=30)

sample = [poisson(10) for _ in range(10_000)]

exp_space = np.linspace(0, 25, 100)
expected = pmf(exp_space, 10)
plt.plot(exp_space, expected)
plt.title(r'$\lambda = {}$'.format(10))
plt.hist(sample, density=True, bins=20)

sample = [poisson(5) for _ in range(10_000)]

exp_space = np.linspace(0, 16, 100)
expected = pmf(exp_space, 5)
plt.plot(exp_space, expected)
plt.title(r'$\lambda = {}$'.format(5))
plt.hist(sample, density=True, bins=15)

sample = [poisson(1) for _ in range(10_000)]

exp_space = np.linspace(0, 7, 100)
expected = pmf(exp_space, 1)
plt.plot(exp_space, expected)
plt.title(r'$\lambda = {}$'.format(1))
plt.hist(sample, density=True, bins=8)