Let us apply the Poisson distribution in form of an example. We will look at a one-hundred-year flood, a concept that is widely used in river engineering to design flood control measures.

Recall the mathematical notation of a Poisson random variable:

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

where $$\lambda$$ is a positive real number, which represents the average number of events occurring during a fixed time interval, and $$e \approx 2.7182818$$.

A 100-year flood is a shorthand term for a flood with an annual exceedance probability of 1 % and an average recurrence interval of 100 years. The term can be confusing as is might suggest that a given flood occurs once every 100-years. However, this is not true. A flood with an annual exceedance probability of 1 % indicates that with a probability of 0.01 a flood of a magnitude corresponding to a 100-year flood occurs in any given year.

Putting this in the context of a Poisson distribution the expected value, $$E(x) = \lambda$$, of such a flood event during the fixed interval of 100 years is set to $$\lambda= 100\times 0.01 = 1$$. Consequently, the Poisson random variable $$X$$ is the number of occurrences, which of course may take several values depending on the question. We may be interested in the probability that such a flood event will not occur during the 100-year interval $$P(x=0)$$. Or, we want to know the probability of such a flood event occurring exactly once during the 100-year interval, thus $$P(x=1)$$. We might also want to know the probability that two or more of such flood events occur during the 100-year interval, thus $$P(x \ge 2)$$. Plugging these values in the equation from above yields:

$$\lambda = 1, x = 0,1,2,..,n$$

$P(X = 0) = e^{-1}\frac{1 \times 0}{0!}, \qquad \text{for } x = 0$ $P(X = 1) = e^{-1}\frac{1 \times 1}{1!}, \qquad \text{for } x = 1$ $P(X \ge 2) = \sum_{i=2}^n e^{-1}\frac{1 \times x_i}{x_i!}, \qquad \text{for } x_i = 2,3,...,n$

We turn to R to do the calculations. We will make use of the dpois() and the ppois() functions.

lambda <- 1 # set expected value
x <- c(0, 1) # set Poisson random variables of interest

dpois(x = x[1], lambda = lambda)
## [1] 0.3678794
dpois(x = x[2], lambda = lambda)
## [1] 0.3678794
# recall: lower.tail = TRUE, probabilities are P[X <= x], otherwise P[X > x]
ppois(q = x[2], lambda = lambda, lower.tail = FALSE)
## [1] 0.2642411

The results indicate, that the probability of no flood $$P(X = 0)$$ with a magnitude corresponding to a 100-year flood occurring during a period of 100 years is 0.37. Interestingly, this is as likely as the occurrence of exactly one flood $$P(X = 1)$$. The probability that two or more $$(P(X \ge 2))$$ of such flood events occur during the 100-year interval is 0.26 and thus lower. However, please note that the probability that two or more $$(P(X \ge 2))$$ of such flood events occur during the 100-year interval is approximately 26 %!

For a sanity check we sum up the probabilities $$P(X = 0)$$, $$P(X = 1)$$ and $$P(X \ge 2)$$, which should yield 1:

dpois(x = x[1], lambda = lambda) +
dpois(x = x[2], lambda = lambda) +
ppois(q = x[2], lambda = lambda, lower.tail = FALSE)
## [1] 1

In order to gain a better intuition we plot the probabilities of the Poisson random variable $$x = 0,1,2,3,4,\ge 5$$.

lambda <- 1 # set expected value
out <- NULL
out[1] <- dpois(x = 0, lambda = lambda)
out[2] <- dpois(x = 1, lambda = lambda)
out[3] <- dpois(x = 2, lambda = lambda)
out[4] <- dpois(x = 3, lambda = lambda)
out[5] <- dpois(x = 4, lambda = lambda)
out[6] <- ppois(q = 5, lambda = lambda, lower.tail = FALSE)

# plot
barplot(out,
ylab = "Probability P(X = x)",
xlab = "Number of floods (x)",
names = paste(c("0", "1", "2", "3", "4", ">= 5")),
main = 'Probability of the occurence of "100-year floods" \nduring a period of 100 years'
)

Exercise: At the bus stop Emmichstrasse near FU Campus Lankwitz the transport company schedules 6 buses per hour. Thus, we may regard 6 buses per hour as the expected value lambda. Let us look at this declared intention as a Poisson process and calculate the probability for X = 0, X = 1, …, X = 12 buses/h. Provide a graph for the results!

### your code here
Show code
lambda <- 6
prob_buses <- NA
for (x in 0:12) {
prob_buses[x + 1] <- dpois(x = x, lambda = lambda)
}

barplot(prob_buses,
ylab = "Probability P(X = x)", xlab = "Number of buses per hour (x)",
col = "#f7f90e", names = 0:12,
main = "Probability for the expected number of buses\n at the bus stop Emmichstrasse"
)

Citation

The E-Learning project SOGA-R was developed at the Department of Earth Sciences by Kai Hartmann, Joachim Krois and Annette Rudolph. You can reach us via mail by soga[at]zedat.fu-berlin.de.

You may use this project freely under the Creative Commons Attribution-ShareAlike 4.0 International License.

Please cite as follow: Hartmann, K., Krois, J., Rudolph, A. (2023): Statistics and Geodata Analysis using R (SOGA-R). Department of Earth Sciences, Freie Universitaet Berlin.