Poisson distribution probabilities using R

Poisson distribution probabilities using R

In this tutorial, you will learn about how to use dpois(), ppois(), qpois() and rpois() functions in R programming language to compute the individual probabilities, cumulative probabilities, quantiles and to generate random sample for Poisson distribution.

Before we discuss R functions for Poisson distribution, let us see what is Poisson distribution.

Poisson Distribution

Poisson distribution helps to describe the probability of occurrence of a number of events in some given time interval or in a specified region. The time interval may be of any length, such as a minutes, a day, a week etc.

Let $X\sim P(\lambda)$. Then the probability distribution of $X$ is

$$ \begin{aligned} P(X=x)= \left\{ \begin{array}{ll} \frac{e^{-\lambda}\lambda^x}{x!} , & \hbox{$x=0,1,2,\cdots; \lambda>0$;} \\ 0, & \hbox{Otherwise.} \end{array} \right. \end{aligned} $$

where $\lambda$ is the parameter of Poisson distribution. $\lambda$ is the average of random variable $X$.

Read more about the theory and results of Poisson distribution here.

Poisson probabilities using dpois() function in R

For discrete probability distribution, density is the probability of getting exactly the value $x$ (i.e., $P(X=x)$).

The syntax to compute the probability at $x$ for Poisson distribution using R is

dpois(x,lambda)

where

  • x : the value(s) of the variable and,
  • lambda : average.

The dpois() function gives the probability for given value(s) x and lambda.

Numerical Problem for Poisson Distribution

To understand the four functions dpois(), ppois(), qpois() and rpois(), let us take the following numerical problem.

Poisson Distribution Example

An old bus breaks down an average of 3 times per month. Use the Poisson probability distribution.

(a) Find the probability that exactly 2 breakdowns during next month.
(b) Plot the graph of Poisson probability distribution.
(c) What is the probability that at most one breakdown during next month?
(d) What is the probability at least 3 breakdowns during next month?
(e) What is the probability that 2 to 4 (inclusive) breakdowns during next month?
(f) Plot the graph of cumulative Poisson probabilities.
(g) What is the value of $c$, if $P(X\leq c) \geq 0.60$?
(h) Simulate 100 Poisson distributed random variables with $\lambda = 3$.

Example 1: How to use dpois() function in R?

To find the probability that exactly 2 breakdowns during next month, we need to use dpois() function.

Let $X$ denote the number of breakdowns during a month. is $3$. Then $X\sim P(3)$.

First let us define the given terms as

# assign the average value
lambda <- 3
# x is the possible values of random variable x
x <- 0:16

The probability mass function of $X$ is

$$ \begin{aligned} P(X=x) &= \frac{e^{-3}(3)^x}{x!},\\ & \quad x=0,1,2,\cdots \end{aligned} $$

For part (a), we need to find the probability $P(X = 2)$.

First I will show you how to calculate this probability using manual calculation, then I will show you how to compute the same probability using dpois() function in R.

(a) The probability that exactly 2 breakdowns during next month is

$$ \begin{aligned} P(X = 2) & =\frac{e^{-3}3^2}{2!} \\ & = 0.2240418\\ \end{aligned} $$

The above probability can be calculated using dpois(2,3) function in R.

# Compute Poisson probability
result1 <- dpois(2,lambda)
result1
[1] 0.2240418

Example 2 Visualize Poisson probability distribution

Using dpois() function we can compute Poisson distribution probabilities and make a table of it.

## Compute the Poisson probabilities 
px<-dpois(x,lambda)
# make a table 
b_table <- cbind(x,px)
# specify the column names
colnames(b_table) <- c("x", "P(X=x)")
b_table
       x       P(X=x)
 [1,]  0 4.978707e-02
 [2,]  1 1.493612e-01
 [3,]  2 2.240418e-01
 [4,]  3 2.240418e-01
 [5,]  4 1.680314e-01
 [6,]  5 1.008188e-01
 [7,]  6 5.040941e-02
 [8,]  7 2.160403e-02
 [9,]  8 8.101512e-03
[10,]  9 2.700504e-03
[11,] 10 8.101512e-04
[12,] 11 2.209503e-04
[13,] 12 5.523758e-05
[14,] 13 1.274713e-05
[15,] 14 2.731529e-06
[16,] 15 5.463057e-07
[17,] 16 1.024323e-07

Using kable() function from knitr package, we can create table in LaTeX, HTML, Markdown and reStructured Text.

# to make table
library(knitr)
kable(b_table)
x P(X=x)
0 0.0497871
1 0.1493612
2 0.2240418
3 0.2240418
4 0.1680314
5 0.1008188
6 0.0504094
7 0.0216040
8 0.0081015
9 0.0027005
10 0.0008102
11 0.0002210
12 0.0000552
13 0.0000127
14 0.0000027
15 0.0000005
16 0.0000001

(b) Visualizing Poisson Distribution with dpois() function and plot() function in R:

The probability mass function of Poisson distribution with given lambda can be visualized using dpois() function in plot() function as follows:

## Plot the Poisson probability dist
plot(x,px,type="h",xlim=c(0,11),ylim=c(0,max(px)),
     lwd=10, col="blue",ylab="P(X=x)")
title("PMF of Poisson (lambda= 3)")
PMF of Poisson Distribution
PMF of Poisson Distribution

Poisson cumulative probability using ppois() function in R

The syntax to compute the cumulative probability distribution function (CDF) for Poisson distribution using R is

ppois(q,lambda)

where

  • q : the value(s) of the variable,
  • lambda : the mean of $X$.

This function is very useful for calculating the cumulative Poisson probabilities for given value(s) of q (value of the variable x), lambda (mean of $X$).

Example 3: How to use ppois() function in R?

In the above example, for part (c), we need to find the probability $P(X\leq 1)$.

First I will show you how to calculate this probability using manual calculation, then I will show you how to compute the same probability using ppois() and dpois() function in R.

(c) The probability that at most 1breakdown during next month

$$ \begin{aligned} P(X\leq1) &= P(X=0)+ P(X=1)\\ &= \frac{e^{-3}3^{0}}{0!}+\frac{e^{-3}3^{1}}{1!}\\ &= 0.0497871+0.1493612\\ &= 0.1991483 \end{aligned} $$

## Compute cumulative Poisson probability
result2 <- ppois(1,lambda)
result2
[1] 0.1991483

Above probability can also be calculated using dpois() function and the sum() function as follows:

sum(dpois(0:1,lambda))
[1] 0.1991483

Example 4: How to use ppois() function in R?

In the above example, for part (d), we need to find the probability $P(X \geq 3)$.

Numerically the probability that at least 3 breakdowns during next month can be calculated as

$$ \begin{aligned} P(X \geq 3) & =1-P(X\leq 2)\\ &=1-\sum_{x=0}^{2} P(X=x)\\ & = 1- (P(X=0)+P(X=1)+P(X=2))\\ &= 1- \big(0.0497871+0.1493612\big)\\ &\quad +0.2240418\\ & = 0.5768099\\ \end{aligned} $$

To calculate the probability that a random variable $X$ is greater than a given number you can use the option lower.tail=FALSE in ppois() function.

Above probability can be calculated easily using ppois() function with argument lower.tail=FALSE as

$P(X \geq 3) =$ ppois(2,lambda,lower.tail=FALSE)

or by using complementary event as

$P(X \geq 3) = 1- P(X\leq 2)$= 1- ppois(2,lambda)

# compute cumulative Poisson probabilities
# with lower.tail False
ppois(2,lambda,lower.tail=FALSE)
[1] 0.5768099
1-ppois(2,lambda)
[1] 0.5768099

Example 5: How to use ppois() function in R?

One can also use ppois() function to calculate the probability that the random variable $X$ is between two values.

(e) The probability that 2 to 4 (inclusive) breakdowns during next month is

$$ \begin{aligned} P(2 \leq X \leq 4) &= P(X=2)+P(X=3)+P(X=4)\\ &=\frac{e^{-3}3^{2}}{2!}+\frac{e^{-3}3^{3}}{3!}\\ &\quad +\frac{e^{-3}3^{3}}{3!}\\ &= 0.2240418+0.2240418+0.1680314\\ &= 0.616115 \end{aligned} $$

Above event can also be written as

$$ \begin{aligned} P(2 \leq X \leq 4) &= P(X\leq 4) -P(X\leq 1)\\ &= 0.8152632 - 0.1991483\\ &= 0.616115 \end{aligned} $$

The above probability can be calculated using ppois() function as follows:

result3 <- ppois(4,lambda)-ppois(1,lambda)
result3
[1] 0.616115

The above probability can also be calculated using dpois() function along with sum() function.

result4 <- sum(dpois(2:4,lambda))
result4
[1] 0.616115

The first command compute the Poisson probability for $x=2$, $x=3$ and $x=4$. Then add all the probabilities using sum() function and store the result in result4.

Example 6: Visualize the cumulative Poisson probability distribution

# the value of x
x <- 0:16
## Compute the Poisson probabilities 
px <- dpois(x,lambda)
## Compute the cumulative Poisson probabilities 
Fx <- ppois(x,lambda)
## make a table 
b_table2 <- cbind(x,px,Fx)
## assign column names
colnames(b_table2) <- c("x", "P(X=x)","P(X<=x)")
# display result
b_table2
       x       P(X=x)    P(X<=x)
 [1,]  0 4.978707e-02 0.04978707
 [2,]  1 1.493612e-01 0.19914827
 [3,]  2 2.240418e-01 0.42319008
 [4,]  3 2.240418e-01 0.64723189
 [5,]  4 1.680314e-01 0.81526324
 [6,]  5 1.008188e-01 0.91608206
 [7,]  6 5.040941e-02 0.96649146
 [8,]  7 2.160403e-02 0.98809550
 [9,]  8 8.101512e-03 0.99619701
[10,]  9 2.700504e-03 0.99889751
[11,] 10 8.101512e-04 0.99970766
[12,] 11 2.209503e-04 0.99992861
[13,] 12 5.523758e-05 0.99998385
[14,] 13 1.274713e-05 0.99999660
[15,] 14 2.731529e-06 0.99999933
[16,] 15 5.463057e-07 0.99999988
[17,] 16 1.024323e-07 0.99999998
kable(b_table2)
x P(X=x) P(X<=x)
0 0.0497871 0.0497871
1 0.1493612 0.1991483
2 0.2240418 0.4231901
3 0.2240418 0.6472319
4 0.1680314 0.8152632
5 0.1008188 0.9160821
6 0.0504094 0.9664915
7 0.0216040 0.9880955
8 0.0081015 0.9961970
9 0.0027005 0.9988975
10 0.0008102 0.9997077
11 0.0002210 0.9999286
12 0.0000552 0.9999839
13 0.0000127 0.9999966
14 0.0000027 0.9999993
15 0.0000005 0.9999999
16 0.0000001 1.0000000

The cumulative probability distribution of Poisson distribution with given lambda can be visualized using plot() function with argument type="s" (step function) as follows:

# Plot the cumulative Poisson dist
plot(x,Fx,type="s",lwd=2,col="blue",
     ylab=expression(P(X<=x)),
main="Distribution Function of P(lambda = 3)")
CDF of Poisson Distribution
CDF of Poisson Distribution

Poisson Distribution Quantiles using qpois() in R

The syntax to compute the quantiles of Poisson distribution using R is

qpois(p,lambda)

where

  • p : the value(s) of the probabilities,
  • lambda : the mean of $X$

The function qpois(p,lambda) gives $100*p^{th}$ quantile of Poisson distribution for given value of p and lambda.

The $p^{th}$ quantile is the smallest value of Poisson random variable $X$ such that $P(X\leq x) \geq p$.

It is the inverse of ppois() function. That is, inverse cumulative probability distribution function for Poisson distribution.

Example 7: How to use qpois() function in R?

In part (g), we need to find the value of $c$ such a that $P(X\leq c) \geq 0.60$. That is we need to find the $60^{th}$ quantile of given Poisson distribution.

# compute the quantile for Poisson dist
qpois(0.60,lambda)
[1] 3

From the above table of Poisson probabilities and cumulative probabilities, it is clear that $60^{th}$ percentile is 3.

Visualize the quantiles of Poisson Distribution

The quantiles of Poisson distribution with given p and lambda can be visualized using plot() function as follows:

p <- seq(0,1,by=0.02)
qx <- qpois(p,lambda)
# Plot the quantiles of Poisson dist
plot(p,qx,type="s",lwd=2,col="darkred",
     ylab="quantiles",
main="Quantiles of P(lambda = 4)")
Quantiles of Poisson Distribution
Quantiles of Poisson Distribution

Simulating Poisson random variable using rpois() function in R

The general R function to generate random numbers from Poisson distribution is

rpois(n,lambda)

where,

  • n is the sample size,
  • lambda is the mean of $X$.

The function rpois(n,lambda) generates n random numbers from Poisson distribution with the average lambda.

Example 8: How to use rpois() function in R?

In part (h), we need to generate 100 random numbers from Poisson distribution with average $\lambda$.

We can use rpois() function to generate random numbers from Poisson distribution.

## initialize sample size to generate
n <- 100
# Simulate 100 values From Poisson dist
x_sim <- rpois(n,lambda)
# print values at console
x_sim  
  [1] 2 4 2 5 6 0 3 5 3 3 6 3 4 3 1 5 2 0 2 6 5 4 3 8 4 4 3 3 2 1 6 5 4 4 0 3 4
 [38] 2 2 2 1 2 2 2 1 1 2 3 2 5 0 3 4 1 3 2 1 4 5 2 4 1 2 2 4 3 4 4 4 3 4 3 4 0
 [75] 3 2 2 3 2 1 2 4 2 4 1 3 7 5 5 1 1 4 2 4 2 1 4 1 3 3

The frequency table for Poisson simulated data x_sim can be obtained using table() command.

## Print the frequency table
table(x_sim)
x_sim
 0  1  2  3  4  5  6  7  8 
 5 14 24 20 22  9  4  1  1 
## Plot the simulated data
plot(table(x_sim),xlab="x",ylab="frequency",
     lwd=10,col="skyblue",
     main="Simulated data from P(3) dist")
Random Sample from Poisson Distribution
Random Sample from Poisson Distribution

If you use same function again, R will generate another set of random numbers from $P(3)$.

# Simulate 100 values From Poisson dist
x_sim_2 <- rpois(n,lambda)
# print values at console
x_sim_2
  [1] 4 4 1 4 1 2 1 2 2 2 2 4 3 3 5 6 3 1 3 4 2 2 2 3 4 2 3 6 1 2 2 4 0 5 2 0 4
 [38] 3 1 4 2 2 1 5 1 0 4 3 1 5 2 3 0 4 6 1 5 1 4 8 4 0 5 3 3 2 4 7 1 1 3 5 3 4
 [75] 6 3 4 2 1 2 2 2 3 2 2 5 2 3 2 1 3 3 3 4 2 5 6 2 2 1

The frequency table of simulated data from Poisson distribution is as follow:

## Print the frequency table
table(x_sim_2)
x_sim_2
 0  1  2  3  4  5  6  7  8 
 5 16 27 19 17  9  5  1  1 
## Plot the simulated data
plot(table(x_sim_2),xlab="x",ylab="frequency",
     lwd=10,col="skyblue",
     main="Simulated data from P(3) dist")
Random Sample from Poisson Distribution 2
Random Sample from Poisson Distribution 2

For the simulation purpose to reproduce same set of random numbers, one can use set.seed() function.

# set seed for reproducibility
set.seed(1457)
# Simulate 100 values From Poisson dist
x_sim_3 <- rpois(n,lambda)
# print values at console
x_sim_3
  [1] 4 3 3 2 3 1 3 2 4 2 4 2 1 0 2 3 2 3 1 4 1 2 3 4 3 2 4 6 2 2 5 1 1 2 1 4 4
 [38] 1 3 2 4 3 1 1 1 2 3 6 2 3 1 3 6 2 6 1 2 0 5 4 2 2 2 2 3 5 6 4 4 2 4 1 1 3
 [75] 2 4 2 7 3 0 2 0 2 1 4 3 4 4 1 4 3 3 4 1 1 3 2 5 3 3

The frequency table of x_sim_3 is as follows:

## Print the frequency table
table(x_sim_3)
x_sim_3
 0  1  2  3  4  5  6  7 
 4 19 26 22 19  4  5  1 
plot(table(x_sim_3),xlab="x",ylab="frequency",
     lwd=10,col="purple",
     main="Simulated data from P(3) dist")
Random Sample from Poisson Distribution 3
Random Sample from Poisson Distribution 3
set.seed(1457)
# Simulate 100 values From Poisson dist
x_sim_4 <- rpois(n,lambda)
# print values at console
x_sim_4
  [1] 4 3 3 2 3 1 3 2 4 2 4 2 1 0 2 3 2 3 1 4 1 2 3 4 3 2 4 6 2 2 5 1 1 2 1 4 4
 [38] 1 3 2 4 3 1 1 1 2 3 6 2 3 1 3 6 2 6 1 2 0 5 4 2 2 2 2 3 5 6 4 4 2 4 1 1 3
 [75] 2 4 2 7 3 0 2 0 2 1 4 3 4 4 1 4 3 3 4 1 1 3 2 5 3 3

The frequency table of x_sim_4 is as follows:

## Print the frequency table
table(x_sim_4)
x_sim_4
 0  1  2  3  4  5  6  7 
 4 19 26 22 19  4  5  1 
plot(table(x_sim_4),xlab="x",ylab="frequency",
     lwd=10,col="purple",
     main="Simulated data from P(3) dist")
Random Sample from Poisson Distribution 3
Random Sample from Poisson Distribution 4

Since we have used set.seed(1457) function for both the simulation, the x_sim_3 and x_sim_4 are same.

To learn more about other discrete and continuous probability distributions using R, go through the following tutorials:

Binomial distribution in R
Geometric distribution in R
Negative Binomial distribution in R
Hypergeometric distribution in R

Continuous Distributions Using R

Uniform distribution in R
Exponential distribution in R
Normal distribution in R
Log-Normal distribution in R
Beta distribution in R
Gamma distribution in R
Cauchy distribution in R
Laplace distribution in R
Logistic distribution in R
Weibull distribution in R

Endnote

In this tutorial, you learned about how to compute the probabilities, cumulative probabilities and quantiles of Poisson distribution in R programming. You also learned about how to simulate a Poisson distribution using R programming.

To learn more about R code for discrete and continuous probability distributions, please refer to the following tutorials:

Probability Distributions using R

Let me know in the comments below, if you have any questions on Poisson Distribution using R and your thought on this article.

Leave a Comment