Simulation of a Poisson Process using R and ggplot2

Introduction

A Poisson process is a stochastic process that represents the number of events occurring in a fixed interval of time or space, where these events occur independently and at a constant average rate. The Poisson distribution is commonly used to model the number of arrivals (events) in a given time period. In this article, we will explore how to simulate a Poisson process using R and ggplot2.

Background

The Poisson distribution is defined as follows:

P(k; λ) = (e^(-λ) * (λ^k)) / k!

where P(k; λ) is the probability of observing exactly k events, λ is the average rate of events, and e is the base of the natural logarithm.

To simulate a Poisson process, we need to generate random numbers that follow this distribution. One way to do this is by using the inverse transform method, which involves generating a uniform random variable U and then transforming it into a Poisson distributed random variable using the cumulative distribution function (CDF) of the Poisson distribution.

However, in this article, we will focus on using R’s built-in rpois function to generate Poisson distributed random variables directly.

Creating a Simulation Function

To simulate a Poisson process, we need to create a simulation function that generates arrival times and counts. We can use the make_sample_df function provided in the answer to achieve this.

Here is an example implementation of a Poisson process simulation function:

## make_sample_df function
make_sample_df <- function(run, tmax, lambda) {
  # Initialize variables
  x <- 0
  
  # Grow vector by adding random samples until sum reaches target duration
  while (sum(x) < tmax) {
    x <- c(x, rexp(1, lambda))
  }
  
  # Return data frame with cumulative times, counts, and run grouping variable
  data.frame(
    t = cumsum(x), 
    N = seq_along(x), 
    run = rep(run, length(x))
  )
}

Plotting the Poisson Process

To plot the Poisson process, we can use the plot_poisson function provided in the answer:

## plot_poisson function
plot_poisson <- function(runs, tmax, lambda) {
  # Create data frames for each run and combine them into a single data frame
  df <- do.call("rbind", lapply(seq(runs), make_sample_df, tmax, lambda))
  
  # Plot Poisson process with ggplot2
  ggplot2::ggplot(df, aes(t, N, group = run)) + 
    geom_step(alpha = 0.25) +
    labs(title = paste(runs, "runs of Poisson process with lambda", lambda)) +
    theme(legend.position = "none") +
    coord_cartesian(xlim = c(0, tmax))
}

Example Usage

We can now use these functions to simulate a Poisson process and plot it:

# Simulate 10 runs of Poisson process with lambda = 0.7
plot_poisson(runs = 10, tmax = 100, lambda = 0.7)

# Simulate 100 runs of Poisson process with lambda = 0.7
plot_poisson(runs = 100, tmax = 100, lambda = 0.7)

This will generate two plots showing the simulated Poisson processes.

Conclusion

In this article, we have explored how to simulate a Poisson process using R and ggplot2. We created a simulation function make_sample_df that generates arrival times and counts, and used it in a plotting function plot_poisson to visualize the Poisson process. The example usage demonstrates how to use these functions to simulate multiple runs of the Poisson process and plot them.

Last modified on 2024-06-22