Four ways to calculate the probability that if it rains, I have an umbrella, using R

This post shows four ways to calculate the conditional probability that I have an umbrella, given that it’s raining, \(P(umbrella|raining)\). The first way uses the ratio formula, \(P(raining \land umbrella)/P(raining)\). The second uses a three-valued logical connective, the conditional event, which can be traced back to the work of Bruno de Finetti in the 1930s (see, e.g., Baratgin, 2021). The third filters the dataset to rows where it’s raining before counting the proportion of those where I had an umbrella. The fourth uses logistic regression.

Here’s a simulated dataset about rain and umbrellas, where each row is a day.

library(conflicted)
library(tidyverse)
conflicts_prefer(dplyr::filter)

days_n <- 20
set.seed(5)

sim_dat <- data.frame(day = 1:days_n, raining = rbinom(days_n, 1, .5)) |>
  mutate(umbrella = rbinom(days_n, 1, ifelse(raining, .8, .2))) |>
  mutate(across(raining:umbrella, as.logical))
sim_dat

The variable raining denotes whether it was raining and umbrella denotes whether I had an umbrella with me.

Here are the joint probabilities:

joints <- sim_dat |>
  group_by(raining, umbrella) |>
  tally() |>
  ungroup() |>
  mutate(p = n / sum(n)) |>
  select(-n) |>
  arrange(desc(raining), desc(umbrella))
joints

Then

\[P(umbrella|raining) = \frac{P(raining \land umbrella)}{P(raining)}\]

\(P(raining \land umbrella)\), i.e., the probability that it’s raining and I have an umbrella, is:

raining_and_umbrella_p <- joints |>
  filter(raining & umbrella) |>
  pull(p)
raining_and_umbrella_p 

[1] 0.45

\(P(raining)\) is:

raining_p <- joints |>
  filter(raining) |>
  pull(p) |>
  sum()
raining_p

[1] 0.5

So, \(P(umbrella|raining)\) is:

raining_and_umbrella_p / raining_p

[1] 0.9

Another way to calculate this is using the conditional event, \(B \dashv A\) (“B given A”), here using Hailperin’s (1996) formulation (he calls it the suppositional):

dont_care <- function(x)
  ifelse(x, NA, x)

`%-|%` <- function(y, x)
  dont_care(!x) | (x & y)

The don’t_care function has one parameter and returns the missing value, NA, if that parameter is TRUE. \(B \dashv A\) is then defined using this so it gets the value NA if \(A\) is FALSE, otherwise it’s equal to \(A \land B\). Here’s a table:

vals <- c(TRUE, FALSE)
expand.grid(A = vals, B = vals) |>
  mutate(`B -| A` = B %-|% A) |>
  arrange(desc(A), desc(B))

We can apply this to the original dataset:

sim_dat <- sim_dat |>
  mutate(`umbrella -| raining` = umbrella %-|% raining)
sim_dat

Now for \(P(umbrella|raining)\) we drop the NAs and count the proportion of the time that \(umbrella \dashv raining\) is TRUE:

sim_dat |>
  pull(`umbrella -| raining`) |>
  mean(na.rm = TRUE)

[1] 0.9

An easier, though less fun, way would be just:

sim_dat |>
  filter(raining) |>
  pull(umbrella) |>
  mean()

[1] 0.9

Finally, we could use logistic regresion:

logistic <- glm(umbrella ~ raining, data = sim_dat, family = binomial)
predict(logistic, newdat = data.frame(raining = TRUE), type = "response")

  1 
0.9 

References

Baratgin, J. (2021). Discovering early de Finetti’s writings on trivalent theory of conditionals. Argumenta, 6, 267–291.

Hailperin, T. (1996). Sentential Probability Logic. Origins, Development, Current Status, and Technical Applications. Lehigh University Press.