{
"cells": [
{
"cell_type": "markdown",
"metadata": {
"button": false,
"deletable": true,
"id": "xxc43bkDjrvO",
"new_sheet": false,
"run_control": {
"read_only": false
}
},
"source": [
"Peter Norvig, 12 Feb 2016
Revised 17 Feb 2018
Python 3.8 updates 2024
\n",
"\n",
"# A Concrete Introduction to Probability (using Python)\n",
"\n",
"In 1814, Pierre-Simon Laplace provided a classic [definition](https://en.wikipedia.org/wiki/Classical_definition_of_probability) of probability:\n",
"\n",
">*Probability theory is nothing but common sense reduced to calculation. ... [Probability] is thus simply a fraction whose numerator is the number of favorable cases and whose denominator is the number of all the cases possible ... when nothing leads us to expect that any one of these cases should occur more than any other.*\n",
"\n",
"| |\n",
"|--|\n",
"||\n",
"|Pierre-Simon Laplace (1749–1827)|\n",
"\n",
"\n",
"Laplace nailed it. To untangle a probability problem, all you have to do is define exactly what the cases are, and carefully count the favorable and total cases. Let's be clear on our vocabulary words:\n",
"\n",
"\n",
"- **[Trial](https://en.wikipedia.org/wiki/Experiment_%28probability_theory%29):**\n",
" A single occurrence with an outcome that is uncertain until it happens.\n",
"
*For example, rolling a single die.*\n",
"- **[Outcome](https://en.wikipedia.org/wiki/Outcome_%28probability%29):**\n",
" A possible result of a trial; one particular state of the world. What Laplace calls a **case.**\n",
"
*For example: the die comes up as* `4`.\n",
"- **[Sample Space](https://en.wikipedia.org/wiki/Sample_space):**\n",
" The set of all possible outcomes for the trial.\n",
"
*For example,* `{1, 2, 3, 4, 5, 6}`.\n",
"- **[Event](https://en.wikipedia.org/wiki/Event_%28probability_theory%29):**\n",
" A subset of the sample space, a set of outcomes that together have some property we are interested in.\n",
"
*For example, the event \"even die roll\" is the set of outcomes* `{2, 4, 6}`.\n",
"- **[Probability](https://en.wikipedia.org/wiki/Probability_theory):**\n",
" As Laplace said, the probability of an event with respect to a sample space is the \"number of favorable cases\" (outcomes from the sample space that are in the event) divided by the \"number of all the cases\" in the sample space, assuming \"nothing leads us to expect that any one of these cases should occur more than any other.\" Since this is a proper fraction, probability will always be a number between 0 (representing an impossible event) and 1 (representing a certain event).\n",
"
*For example, the probability of an even die roll is 3/6 = 1/2.*\n",
"\n",
"This notebook will explore these concepts in a concrete way using Python code. The code is meant to be succint and explicit, and fast enough to handle sample spaces with millions of outcomes. If you need to handle trillions, you'll want a more efficient implementation. I also have [another notebook](http://nbviewer.jupyter.org/url/norvig.com/ipython/ProbabilityParadox.ipynb) that covers paradoxes in Probability Theory.\n",
"\n",
"First some imports we will need later, and some type definitions:"
]
},
{
"cell_type": "code",
"execution_count": 1,
"metadata": {},
"outputs": [],
"source": [
"from fractions import Fraction\n",
"from itertools import combinations, product\n",
"from typing import *\n",
"import math\n",
"import random\n",
"\n",
"Space = set # A sample space is a set of all possible outcomes\n",
"Event = set # An event is a subset of the sample space"
]
},
{
"cell_type": "markdown",
"metadata": {
"button": false,
"deletable": true,
"id": "Lcnp5JR3jrvP",
"new_sheet": false,
"run_control": {
"read_only": false
}
},
"source": [
"# P is for Probability\n",
"\n",
"The code below implements Laplace's quote directly: *Probability is thus simply a fraction whose numerator is the number of favorable cases and whose denominator is the number of all the cases possible.*"
]
},
{
"cell_type": "code",
"execution_count": 2,
"metadata": {},
"outputs": [],
"source": [
"def P(event: Event, space: Space) -> Fraction:\n",
" \"\"\"The probability of an event, given a sample space:\n",
" the number of favorable cases divided by the number of all the cases possible.\"\"\"\n",
" return Fraction(number_cases(favorable(event, space)),\n",
" number_cases(space))\n",
"\n",
"favorable = set.intersection # Favorable cases are in the event and also in the sample space\n",
"number_cases = len # The number of cases is the length, or size, of a set"
]
},
{
"cell_type": "markdown",
"metadata": {
"button": false,
"deletable": true,
"id": "JTy695a1jrvQ",
"new_sheet": false,
"run_control": {
"read_only": false
}
},
"source": [
"\n",
"# Warm-up Problem: Die Roll"
]
},
{
"cell_type": "markdown",
"metadata": {
"button": false,
"deletable": true,
"id": "3FK0zMz6jrvQ",
"new_sheet": false,
"run_control": {
"read_only": false
}
},
"source": [
"What's the probability of an even number with a single roll of a six-sided fair die? \n",
"\n",
"Mathematicians traditionally use a single capital letter to denote a sample space; I'll use `D` for the die:"
]
},
{
"cell_type": "code",
"execution_count": 3,
"metadata": {
"button": false,
"deletable": true,
"id": "p42G7bpSjrvQ",
"new_sheet": false,
"outputId": "0d464c7f-fe16-45e0-bd0d-9cb8898eb344",
"run_control": {
"read_only": false
}
},
"outputs": [],
"source": [
"D = {1, 2, 3, 4, 5, 6} # the sample space for the die"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"I can then define the event of rolling an even number, and ask for the probability of that event:"
]
},
{
"cell_type": "code",
"execution_count": 4,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"Fraction(1, 2)"
]
},
"execution_count": 4,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"even = {2, 4, 6} # the event of an even roll\n",
"\n",
"P(even, D)"
]
},
{
"cell_type": "markdown",
"metadata": {
"id": "2YxbcGMCjrvR"
},
"source": [
"Good to confirm what we already knew. We can explore some other events:"
]
},
{
"cell_type": "code",
"execution_count": 5,
"metadata": {
"id": "7Uuwz8b6jrvR"
},
"outputs": [
{
"data": {
"text/plain": [
"Fraction(1, 2)"
]
},
"execution_count": 5,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"odd = {1, 3, 5, 7, 9, 11, 13}\n",
"\n",
"P(odd, D)"
]
},
{
"cell_type": "code",
"execution_count": 6,
"metadata": {
"id": "L5GOF5nKjrvR",
"outputId": "2058c35e-86a2-4b95-c248-7fe5a197b7ec"
},
"outputs": [
{
"data": {
"text/plain": [
"Fraction(5, 6)"
]
},
"execution_count": 6,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"prime = {2, 3, 5, 7, 11, 13}\n",
"\n",
"P((even | prime), D) # The probability of an even or prime die roll"
]
},
{
"cell_type": "code",
"execution_count": 7,
"metadata": {
"id": "jatr0HDFjrvR",
"outputId": "894b269a-b2a3-4d89-8500-712e587a78c0"
},
"outputs": [
{
"data": {
"text/plain": [
"Fraction(1, 3)"
]
},
"execution_count": 7,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"P((odd & prime), D) # The probability of an odd prime die roll"
]
},
{
"cell_type": "markdown",
"metadata": {
"button": false,
"id": "CCmQ-m1vjrvR",
"new_sheet": false,
"run_control": {
"read_only": false
}
},
"source": [
"# Card Problems\n",
"\n",
"Consider a deck of playiong cards. An individual card has a rank and suit, and will be represented as a string, like `'A♥'` for the Ace of Hearts. There are 4 suits and 13 ranks, so there are 52 cards in a deck:"
]
},
{
"cell_type": "code",
"execution_count": 8,
"metadata": {
"button": false,
"id": "GakwLmAKjrvR",
"new_sheet": false,
"outputId": "866a4a08-4b37-4d3d-d23a-c504f9ff27bc",
"run_control": {
"read_only": false
}
},
"outputs": [
{
"data": {
"text/plain": [
"52"
]
},
"execution_count": 8,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"suits = '♥♠♦♣'\n",
"ranks = 'AKQJT98765432'\n",
"deck = [r + s for r in ranks for s in suits]\n",
"len(deck)"
]
},
{
"cell_type": "markdown",
"metadata": {
"id": "fmlleKe5jrvR"
},
"source": [
"Now I want to define `Hands` as the sample space of all possible 5-card hands that could be dealt from a deck. The function `itertools.combinations` does most of the work; we than concatenate the combinations into space-separateds string using `joins`:\n"
]
},
{
"cell_type": "code",
"execution_count": 9,
"metadata": {
"button": false,
"id": "MfQeXBn5jrvS",
"new_sheet": false,
"outputId": "744a98ef-1ea6-48e7-bab1-ad73ba8e6994",
"run_control": {
"read_only": false
}
},
"outputs": [
{
"data": {
"text/plain": [
"2598960"
]
},
"execution_count": 9,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"def joins(strings) -> Set[str]: return {' '.join(s) for s in strings}\n",
"\n",
"Hands = joins(combinations(deck, 5))\n",
"\n",
"len(Hands)"
]
},
{
"cell_type": "markdown",
"metadata": {
"id": "VpEXxDZgjrvS"
},
"source": [
"There are over 2.5 million hands; too many to look at them all, but we can sample:"
]
},
{
"cell_type": "code",
"execution_count": 10,
"metadata": {
"id": "KMLal8AijrvS",
"outputId": "09ca6e36-1cc0-481b-df56-b3dea7bd34b8"
},
"outputs": [
{
"data": {
"text/plain": [
"['T♥ 7♣ 6♦ 4♠ 2♥',\n",
" 'K♣ Q♦ 8♣ 6♦ 3♥',\n",
" 'A♦ 9♣ 7♥ 7♦ 4♣',\n",
" 'A♥ T♥ 8♠ 6♣ 5♥',\n",
" 'A♠ K♦ Q♥ 7♥ 3♦',\n",
" 'A♦ K♦ J♣ 3♦ 2♦',\n",
" 'K♦ 8♥ 6♠ 6♣ 2♣']"
]
},
"execution_count": 10,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"random.sample(list(Hands), 7)"
]
},
{
"cell_type": "markdown",
"metadata": {
"button": false,
"id": "c34QOuuMjrvS",
"new_sheet": false,
"run_control": {
"read_only": false
}
},
"source": [
"Now we can answer questions like the probability of being dealt a flush (5 cards of the same suit):"
]
},
{
"cell_type": "code",
"execution_count": 11,
"metadata": {
"button": false,
"id": "QMeYmswHjrvS",
"new_sheet": false,
"outputId": "f25a0634-dcce-43cf-8db5-c1cc88aa76ac",
"run_control": {
"read_only": false
}
},
"outputs": [
{
"data": {
"text/plain": [
"Fraction(33, 16660)"
]
},
"execution_count": 11,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"flush = {hand for hand in Hands if any(hand.count(suit) == 5 for suit in suits)}\n",
"\n",
"P(flush, Hands)"
]
},
{
"cell_type": "markdown",
"metadata": {
"button": false,
"id": "LsSlVHUcjrvS",
"new_sheet": false,
"run_control": {
"read_only": false
}
},
"source": [
"Or the probability of four of a kind:"
]
},
{
"cell_type": "code",
"execution_count": 12,
"metadata": {
"button": false,
"id": "9EhqiYrsjrvS",
"new_sheet": false,
"outputId": "38314ccf-2bb5-4439-9d43-3b80afc78975",
"run_control": {
"read_only": false
}
},
"outputs": [
{
"data": {
"text/plain": [
"Fraction(1, 4165)"
]
},
"execution_count": 12,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"four_kind = {hand for hand in Hands if any(hand.count(rank) == 4 for rank in ranks)}\n",
"\n",
"P(four_kind, Hands)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"To make these calculations we need to go through all 2.5 million hands, so this is not the fastest way to compute probabilities, but it is straightforward."
]
},
{
"cell_type": "markdown",
"metadata": {
"button": false,
"deletable": true,
"id": "ww9hGNBVjrvS",
"new_sheet": false,
"run_control": {
"read_only": false
}
},
"source": [
"# Urn Problems\n",
"\n",
"Around 1700, Jacob Bernoulli wrote about removing colored balls from an urn in his landmark treatise *[Ars Conjectandi](https://en.wikipedia.org/wiki/Ars_Conjectandi)*, and ever since then, explanations of probability have relied on [urn problems](https://www.google.com/search?q=probability+ball+urn). (You'd think the urns would be empty by now.)\n",
"\n",
"| |\n",
"|---|\n",
"||\n",
"|Jacob Bernoulli (1655–1705)|\n",
"\n",
"For example, here is a three-part problem [adapted](http://mathforum.org/library/drmath/view/69151.html) from mathforum.org:\n",
"\n",
"> *An urn contains 6 blue, 9 red, and 8 white balls. We select six balls at random. What is the probability of each of these outcomes:*\n",
">\n",
"> - *All balls are red*.\n",
"> - *3 are blue, and 1 is red, and 2 are white*.\n",
"> - *Exactly 4 balls are white*.\n",
"\n",
"We'll start by defining the contents of the urn. We will need a way to name the balls so that we call tell blue from red, but also tell one red ball from another red ball. I'll use the function `names('R', 9)` to create names for 9 red balls:"
]
},
{
"cell_type": "code",
"execution_count": 13,
"metadata": {
"button": false,
"deletable": true,
"id": "fqxshx8yjrvS",
"new_sheet": false,
"outputId": "28605592-fbef-4b8f-968e-d04d07426f52",
"run_control": {
"read_only": false
}
},
"outputs": [
{
"data": {
"text/plain": [
"['B1',\n",
" 'B2',\n",
" 'B3',\n",
" 'B4',\n",
" 'B5',\n",
" 'B6',\n",
" 'R1',\n",
" 'R2',\n",
" 'R3',\n",
" 'R4',\n",
" 'R5',\n",
" 'R6',\n",
" 'R7',\n",
" 'R8',\n",
" 'R9',\n",
" 'W1',\n",
" 'W2',\n",
" 'W3',\n",
" 'W4',\n",
" 'W5',\n",
" 'W6',\n",
" 'W7',\n",
" 'W8']"
]
},
"execution_count": 13,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"def names(name: str, n) -> List[str]:\n",
" \"\"\"A list of `n` distinct names.\"\"\"\n",
" return [name + str(i) for i in range(1, n + 1)]\n",
"\n",
"urn = names('B', 6) + names('R', 9) + names('W', 8)\n",
"urn"
]
},
{
"cell_type": "markdown",
"metadata": {
"button": false,
"deletable": true,
"id": "xkJWJNeGjrvS",
"new_sheet": false,
"run_control": {
"read_only": false
}
},
"source": [
"Now we can define the sample space, `U6`, as the set of all combinations of 6 balls: "
]
},
{
"cell_type": "code",
"execution_count": 14,
"metadata": {
"button": false,
"deletable": true,
"id": "7JeJXu7ujrvS",
"new_sheet": false,
"outputId": "8c17ce77-5a51-44c0-96fa-7831941b0aff",
"run_control": {
"read_only": false
}
},
"outputs": [
{
"data": {
"text/plain": [
"['B2 B5 R9 W3 W6 W8',\n",
" 'B4 R4 R7 W1 W2 W4',\n",
" 'B1 B6 R1 W1 W4 W7',\n",
" 'B1 B3 B4 R1 R5 W3',\n",
" 'B1 B3 B5 R5 W3 W6']"
]
},
"execution_count": 14,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"U6 = joins(combinations(urn, 6))\n",
"\n",
"random.sample(list(U6), 5) # A sample from the U6 set"
]
},
{
"cell_type": "markdown",
"metadata": {
"id": "JDUY1Te6jrvS"
},
"source": [
"Now I will define `select('R', 6)` to be the event of picking exactly 6 red balls from the urn:"
]
},
{
"cell_type": "code",
"execution_count": 15,
"metadata": {
"id": "YBIbRWikjrvS"
},
"outputs": [],
"source": [
"def select(color, n, space=U6) -> set:\n",
" \"The subset of the sample space with exactly `n` balls of given `color`.\"\n",
" return {s for s in space if s.count(color) == n}"
]
},
{
"cell_type": "markdown",
"metadata": {
"id": "GaKtWY7BjrvT"
},
"source": [
"Now I can answer the three questions:"
]
},
{
"cell_type": "code",
"execution_count": 16,
"metadata": {
"id": "33QapmtBjrvT",
"outputId": "c47ae0e4-ac43-4c69-849f-cb4314f6d1b5"
},
"outputs": [
{
"data": {
"text/plain": [
"Fraction(4, 4807)"
]
},
"execution_count": 16,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"P(select('R', 6), U6) # Probability that all 6 balls are red."
]
},
{
"cell_type": "code",
"execution_count": 17,
"metadata": {
"id": "VdN2Mzy7jrvT",
"outputId": "c1c0c76e-0cda-4ad9-b63c-2639b347ea57"
},
"outputs": [
{
"data": {
"text/plain": [
"Fraction(240, 4807)"
]
},
"execution_count": 17,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"P(select('B', 3) & select('R', 1) & select('W', 2), U6) # Probability that 3 are blue, and 1 is red, and 2 are white"
]
},
{
"cell_type": "code",
"execution_count": 18,
"metadata": {
"id": "RhexMme5jrvT",
"outputId": "0c80f443-e8a4-45cf-93ad-9225aca9f8a5"
},
"outputs": [
{
"data": {
"text/plain": [
"Fraction(350, 4807)"
]
},
"execution_count": 18,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"P(select('W', 4), U6) # Probability that xxactly 4 balls are white"
]
},
{
"cell_type": "markdown",
"metadata": {
"button": false,
"deletable": true,
"id": "5vAFCkojjrvU",
"new_sheet": false,
"run_control": {
"read_only": false
}
},
"source": [
"## Urn problems via arithmetic\n",
"\n",
"Let's verify these calculations using basic arithmetic, rather than exhaustive counting. First, how many ways can I choose 6 out of 9 red balls? It could be any of the 9 for the first ball, any of 8 remaining for the second, and so on down to any of the remaining 4 for the sixth and final ball. But we don't care about the *order* of the six balls, so divide that product by the number of permutations of 6 things, which is 6!, giving us\n",
"9 × 8 × 7 × 6 × 5 × 4 / 6! = 84. In general, the number of ways of choosing *k* out of *n* items is (*n* choose *k*) = *n*! / ((*n* - *c*)! × c!).\n",
"In Python 3.8+ that is provided as the `math.comb` function."
]
},
{
"cell_type": "markdown",
"metadata": {
"id": "E9OsEwL8jrvU"
},
"source": [
"Now we can verify the answers to the three problems. (Since `P` computes a ratio and `choose` computes a count,\n",
"I multiply the left-hand-side by `N`, the length of the sample space, to make both sides be counts.)"
]
},
{
"cell_type": "code",
"execution_count": 19,
"metadata": {
"id": "_mruTBRfjrvU",
"outputId": "394ac4fe-9805-4e18-e770-1b4f58505cbf"
},
"outputs": [],
"source": [
"N = len(U6)\n",
"\n",
"assert math.comb(9, 6) == 84\n",
"\n",
"assert (N * P(select('R', 6), U6) ==\n",
" math.comb(9, 6))\n",
"\n",
"assert (N * P(select('B', 3) & select('W', 2) & select('R', 1), U6) ==\n",
" math.comb(6, 3) * math.comb(8, 2) * math.comb(9, 1))\n",
"\n",
"assert (N * P(select('W', 4), U6) == \n",
" math.comb(8, 4) * math.comb(6 + 9, 2)) # (6 + 9 non-white balls)"
]
},
{
"cell_type": "code",
"execution_count": 20,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"(Fraction(4, 4807), 0.0008321198252548367)"
]
},
"execution_count": 20,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"P(select('R', 6), U6) , math.comb(9, 6) / N"
]
},
{
"cell_type": "markdown",
"metadata": {
"button": false,
"id": "D7lXV37vjrvY",
"new_sheet": false,
"run_control": {
"read_only": false
}
},
"source": [
"# Non-Equiprobable Outcomes\n",
"\n",
"So far, we have accepted Laplace's assumption that *nothing leads us to expect that any one of these cases should occur more than any other*.\n",
"In real life, we often get outcomes that are not equiprobable--for example, a loaded die favors one side over the others. We will introduce three more vocabulary items:\n",
"\n",
"* [Frequency](https://en.wikipedia.org/wiki/Frequency_%28statistics%29): a non-negative number describing how often an outcome occurs. It Can be a count like 5, or a ratio like 1/6.\n",
"* [Distribution](http://mathworld.wolfram.com/StatisticalDistribution.html): A mapping from outcome to frequency of that outcome. We will allow a sample spaces to be either a set (of equi-probable outcomes) or a distribution.\n",
"* [Probability Distribution](https://en.wikipedia.org/wiki/Probability_distribution): A probability distribution\n",
"is a distribution whose frequencies sum to 1.\n",
"\n",
"\n",
"I'll implemet `Dist` as a subclass of `Counter`, and re-define the type `Space` to be a set or a `Dist`."
]
},
{
"cell_type": "code",
"execution_count": 21,
"metadata": {
"id": "15XeSxtwjrvY"
},
"outputs": [],
"source": [
"class Dist(Counter):\n",
" \"A Distribution of {outcome: frequency} pairs.\"\n",
"\n",
"Space = Union[set, Dist]"
]
},
{
"cell_type": "markdown",
"metadata": {
"id": "v1UgD04cjrvY"
},
"source": [
"Because a `Dist` is a `Counter`, we can initialize it three different ways:\n",
"- With a collection of outcomes (equiprobable or not).\n",
"- With a mapping of `{outcome: frequency}` pairs.\n",
"- With keyword arguments, each being assigned a frequency number.\n",
"\n",
"You can get the same result with any of the three ways:"
]
},
{
"cell_type": "code",
"execution_count": 22,
"metadata": {
"id": "lF98Mu_XjrvY",
"outputId": "fd3e1a79-c821-4500-9f56-c1369076e4c2"
},
"outputs": [],
"source": [
"assert Dist('THHHTTHHT') == Dist({'H': 5, 'T': 4}) == Dist(H=5, T=4)"
]
},
{
"cell_type": "markdown",
"metadata": {
"id": "GzAMe7ZIjrvY"
},
"source": [
"Now I will modify the code to handle distributions.\n",
"Here's my plan:\n",
"\n",
"- The function `P` is unchanged. Laplace's advice still stands!\n",
"- A sample space can be either a set or a distribution, so I will redefine three helper functions:\n",
" - `number_cases` now sums the frequencies in a distribution (it previously counted the length).\n",
" - `favorable` now returns a `Dist` of favorable outcomes and their frequencies (not a `set`).\n",
" - `Fraction` now uses `\"/\"`, not `fractions.Fraction`, because frequencies might be floats.\n",
"- The new function `cast(object, type)` converts an object to the given type (if it is not already of that type)."
]
},
{
"cell_type": "code",
"execution_count": 23,
"metadata": {
"id": "Tt0iAITXjrvY"
},
"outputs": [],
"source": [
"def number_cases(outcomes) -> float:\n",
" \"\"\"The total frequency of all the outcomes.\"\"\"\n",
" return sum(cast(outcomes, Dist).values())\n",
"\n",
"def favorable(event: Event, space: Space) -> Dist:\n",
" \"\"\"A distribution of outcomes from the sample space that are in the event.\"\"\"\n",
" space = cast(space, Dist)\n",
" return Dist({x: space[x] for x in space if x in event})\n",
"\n",
"def Fraction(n, d): return n / d\n",
"\n",
"def cast(object, typ: Type) -> object:\n",
" \"\"\"Convert `object` to `typ`, unless it is already of type `typ`.\"\"\"\n",
" return object if isinstance(object, typ) else typ(object)"
]
},
{
"cell_type": "markdown",
"metadata": {
"id": "HUAFGxQpjrvY"
},
"source": [
"For example, here's the probability of rolling an even number with a crooked die that is loaded to prefer 6:"
]
},
{
"cell_type": "code",
"execution_count": 24,
"metadata": {
"id": "Rwdr1hOOjrvY",
"outputId": "302150cb-66cf-47ac-ec4a-6bc86c515833"
},
"outputs": [
{
"data": {
"text/plain": [
"0.7"
]
},
"execution_count": 24,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"Crooked = Dist({1: 0.1, 2: 0.1, 3: 0.1, 4: 0.1, 5: 0.1, 6: 0.5})\n",
"\n",
"P(even, Crooked)"
]
},
{
"cell_type": "markdown",
"metadata": {
"id": "W5Oqm-f1jrvY"
},
"source": [
"## Example: Birth Counts\n",
"\n",
"As another example, an [article](http://people.kzoo.edu/barth/math105/moreboys.pdf) gives the following counts for two-child families in Denmark, where `GB` means a family where the first-born child is a girl and the second a boy (I'm aware that not all births can be classified as the binary \"boy\" or \"girl,\" but this particular data set was reported that way):"
]
},
{
"cell_type": "code",
"execution_count": 25,
"metadata": {
"button": false,
"id": "bizO1Ic4jrvY",
"new_sheet": false,
"run_control": {
"read_only": false
}
},
"outputs": [],
"source": [
"DK = Dist(GG=121801, GB=126840,\n",
" BG=127123, BB=135138)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"I can define some events:"
]
},
{
"cell_type": "code",
"execution_count": 26,
"metadata": {
"id": "wSgMO5VYjrvZ",
"outputId": "0e10ee5c-6a81-45b4-e916-57af163802de"
},
"outputs": [],
"source": [
"first_girl = {'GG', 'GB'}\n",
"second_girl = {'GG', 'BG'}\n",
"first_boy = {'BB', 'BG'}\n",
"second_boy = {'BB', 'GB'}\n",
"same = {'GG', 'BB'}"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"And ask for the probability that, say, the first or second child is a girl, or that the two children have the same sex:"
]
},
{
"cell_type": "code",
"execution_count": 27,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"0.48667063350701306"
]
},
"execution_count": 27,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"P(first_girl, DK)"
]
},
{
"cell_type": "code",
"execution_count": 28,
"metadata": {
"id": "alGEqnOxjrvZ",
"outputId": "fb9b2cbc-4e25-4e2c-a828-fd3a563f7901"
},
"outputs": [
{
"data": {
"text/plain": [
"0.4872245557856497"
]
},
"execution_count": 28,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"P(second_girl, DK)"
]
},
{
"cell_type": "code",
"execution_count": 29,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"0.5029124959385557"
]
},
"execution_count": 29,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"P(same, DK)"
]
},
{
"cell_type": "markdown",
"metadata": {
"button": false,
"id": "-sxc-yYxjrvZ",
"new_sheet": false,
"run_control": {
"read_only": false
}
},
"source": [
"The numbers say that you are slighltly more likely to have a second child of the same sex, but only by about 0.3%."
]
},
{
"cell_type": "markdown",
"metadata": {
"button": false,
"deletable": true,
"id": "hDmhQHsvjrvZ",
"new_sheet": false,
"run_control": {
"read_only": false
}
},
"source": [
"# Predicates as events\n",
"\n",
"To calculate the probability of an even die roll, I originally said\n",
"\n",
" even = {2, 4, 6}\n",
" \n",
"But that's inelegant—I had to explicitly enumerate all the even numbers from one to six. If I ever wanted to deal with a twelve or twenty-sided die, I would have to go back and redefine `even`. I would prefer to define `even` once and for all like this:"
]
},
{
"cell_type": "code",
"execution_count": 30,
"metadata": {
"button": false,
"deletable": true,
"id": "nb3OdU1BjrvZ",
"new_sheet": false,
"run_control": {
"read_only": false
}
},
"outputs": [],
"source": [
"def even(n: int) -> bool: return n % 2 == 0"
]
},
{
"cell_type": "markdown",
"metadata": {
"button": false,
"deletable": true,
"id": "LXdsmLnfjrvZ",
"new_sheet": false,
"run_control": {
"read_only": false
}
},
"source": [
"Now in order to make `P(even, D)` work, I'll allow an `Event` to be either a collection of outcomes or a `callable` predicate (that is, a function that returns true for outcomes that are part of the event). I don't need to modify `P`, but `favorable` will have to convert a callable `event` to a `set`:"
]
},
{
"cell_type": "code",
"execution_count": 31,
"metadata": {
"button": false,
"deletable": true,
"id": "5ZWWAjJnjrvZ",
"new_sheet": false,
"run_control": {
"read_only": false
}
},
"outputs": [],
"source": [
"Event = Union[set, Callable]\n",
"\n",
"def favorable(event, space):\n",
" \"\"\"A distribution of outcomes from the sample space that are in the event.\"\"\"\n",
" if callable(event):\n",
" event = {x for x in space if event(x)}\n",
" space = cast(space, Dist)\n",
" return Dist({x: space[x] for x in space if x in event})"
]
},
{
"cell_type": "code",
"execution_count": 32,
"metadata": {
"id": "L5cVog-ojrvZ",
"outputId": "65622e94-38d0-44ea-ce07-1010cd3a4cd9"
},
"outputs": [
{
"data": {
"text/plain": [
"Dist({2: 1, 4: 1, 6: 1})"
]
},
"execution_count": 32,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"favorable(even, D)"
]
},
{
"cell_type": "code",
"execution_count": 33,
"metadata": {
"id": "ng1tYqPHjrvZ",
"outputId": "76c65bf4-c022-4eec-c7bd-c87bd87f017d"
},
"outputs": [
{
"data": {
"text/plain": [
"0.5"
]
},
"execution_count": 33,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"P(even, D)"
]
},
{
"cell_type": "markdown",
"metadata": {
"id": "A4NNZhJXjrvZ"
},
"source": [
"I'll define `die(s)` to make a sample space for an *s*-sided die, and `roll(r, s)` to make a sample space for the sum of rolling *r* *s*-sided dice:"
]
},
{
"cell_type": "code",
"execution_count": 34,
"metadata": {
"id": "r_OXRoIujrvZ"
},
"outputs": [],
"source": [
"def die(s: int) -> Space: \n",
" \"\"\"The sample space for an s-sided die.\"\"\"\n",
" return set(range(1, s + 1))\n",
"\n",
"def roll(r: int, s: int) -> Space:\n",
" \"\"\"The sample space for rolling `r` s-sided dice and summing them.\"\"\"\n",
" return Dist(map(sum, product(die(s), repeat=r)))"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"For example, here's the distribution for the sum of two six-sided dice:"
]
},
{
"cell_type": "code",
"execution_count": 35,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"Dist({2: 1, 3: 2, 4: 3, 5: 4, 6: 5, 7: 6, 8: 5, 9: 4, 10: 3, 11: 2, 12: 1})"
]
},
"execution_count": 35,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"roll(2, 6)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"We can check if various dice-related sample spaces are odd or even, or prime:"
]
},
{
"cell_type": "code",
"execution_count": 36,
"metadata": {},
"outputs": [],
"source": [
"def is_prime(n) -> bool: return (n > 1 and not any(n % i == 0 for i in range(2, n)))"
]
},
{
"cell_type": "code",
"execution_count": 37,
"metadata": {
"id": "FxLyaZ5BjrvZ",
"outputId": "3c77cbe8-99a7-4854-9e7b-38a6e9ed4199"
},
"outputs": [
{
"data": {
"text/plain": [
"0.5"
]
},
"execution_count": 37,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"P(even, die(12))"
]
},
{
"cell_type": "code",
"execution_count": 38,
"metadata": {
"id": "EqoaurrSjrva",
"outputId": "aa1a002d-fe3b-427b-9b65-7b7a70163cca"
},
"outputs": [
{
"data": {
"text/plain": [
"0.45454545454545453"
]
},
"execution_count": 38,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"P(even, die(11))"
]
},
{
"cell_type": "code",
"execution_count": 39,
"metadata": {
"id": "8rlgl6sGjrva",
"outputId": "245acf34-0f5a-47f8-b1db-44cd349b5613"
},
"outputs": [
{
"data": {
"text/plain": [
"0.5"
]
},
"execution_count": 39,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"P(prime, die(6))"
]
},
{
"cell_type": "code",
"execution_count": 40,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"0.4166666666666667"
]
},
"execution_count": 40,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"P(prime, roll(2, 6))"
]
},
{
"cell_type": "code",
"execution_count": 41,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"0.4166666666666667"
]
},
"execution_count": 41,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"P(prime, die(12))"
]
},
{
"cell_type": "code",
"execution_count": 42,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"0.3"
]
},
"execution_count": 42,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"P(prime, die(20))"
]
},
{
"cell_type": "code",
"execution_count": 43,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"0.0875"
]
},
"execution_count": 43,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"P(prime, roll(2, 20))"
]
},
{
"cell_type": "code",
"execution_count": 44,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"0.002275"
]
},
"execution_count": 44,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"P(prime, roll(4, 20))"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# Conditional Probability\n",
"\n",
"Conditional Probability is used to answer questions of the form \"What is the probability of event *X*, given that event *Y* has occurred?\" The \"given *Y* has occurred\" part makes a new sample space, one that is *favorable* to the event *Y*. So I'll define `given` to be the same function as `favorable`:"
]
},
{
"cell_type": "code",
"execution_count": 45,
"metadata": {},
"outputs": [],
"source": [
"given = favorable"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"For example, what's the probability that the second child is a girl, given that the first is a girl (with the Denmark data)? And how does that compare to the case where the first is a boy, or is unspecified?"
]
},
{
"cell_type": "code",
"execution_count": 46,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"0.4898669165584115"
]
},
"execution_count": 46,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"P(second_girl, given(first_girl, DK))"
]
},
{
"cell_type": "code",
"execution_count": 47,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"0.48471942072973107"
]
},
"execution_count": 47,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"P(second_girl, given(first_boy, DK))"
]
},
{
"cell_type": "code",
"execution_count": 48,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"0.4872245557856497"
]
},
"execution_count": 48,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"P(second_girl, DK)"
]
},
{
"cell_type": "markdown",
"metadata": {
"id": "Ywi2NzMFjrva"
},
"source": [
"# Fermat and Pascal: The Unfinished Game\n",
"\n",
"\n",
"\n",
"Consider a two-player gambling game consisting of tossing a coin repeatedly. Player H wins as soon as a total of 10 heads come up, but player T wins if a total of 10 tails come up first. If the game is interrupted when H has 8 heads and T has 7 tails, how should the pot of money (which happens to be 100 Francs) be split? Here are some proposals, and arguments against them:\n",
"- It is uncertain, so just split the pot 50-50.\n",
"
*No, because surely H is more likely to win and thus deserves more.*\n",
"- In proportion to each player's current score, so H gets a 8/(8+7) share.\n",
"
*No, because if the score was 0 heads to 1 tail, H should surely deserve more than 0/1.*\n",
"- In proportion to how many tosses the opponent needs to win, so H gets 3/(3+2).\n",
"
*No, because if H is 9 away and T is only 1 away from winning, then it seems that giving H a 1/10 share is too generous.*\n",
"\n",
"In 1654, Blaise Pascal and Pierre de Fermat corresponded on this problem, with Fermat [writing](http://mathforum.org/isaac/problems/prob1.html):\n",
"\n",
">Dearest Blaise,\n",
">\n",
">As to the problem of how to divide the 100 Francs, I think I have found a solution that you will find to be fair. Seeing as I needed only two points to win the game, and you needed 3, I think we can establish that after four more tosses of the coin, the game would have been over. For, in those four tosses, if you did not get the necessary 3 points for your victory, this would imply that I had in fact gained the necessary 2 points for my victory. In a similar manner, if I had not achieved the necessary 2 points for my victory, this would imply that you had in fact achieved at least 3 points and had therefore won the game. Thus, I believe the following list of possible endings to the game is exhaustive. I have denoted 'heads' by an 'h', and tails by a 't.' I have starred the outcomes that indicate a win for myself.\n",
">\n",
"> h h h h * h h h t * h h t h * h h t t *\n",
"> h t h h * h t h t * h t t h * h t t t\n",
"> t h h h * t h h t * t h t h * t h t t\n",
"> t t h h * t t h t t t t h t t t t\n",
"> \n",
">\n",
">I think you will agree that all of these outcomes are equally likely. Thus I believe that we should divide the stakes by the ratio 11:5 in my favor, that is, I should receive (11/16)×100 = 68.75 Francs, while you should receive 31.25 Francs.\n",
">\n",
">\n",
">I hope all is well in Paris,\n",
">\n",
">Your friend and colleague,\n",
">\n",
">Pierre\n",
"\n",
"Pascal agreed with this solution, and [replied](http://mathforum.org/isaac/problems/prob2.html) with a generalization that made use of his previous invention, Pascal's Triangle. There's even [a book](https://smile.amazon.com/Unfinished-Game-Pascal-Fermat-Seventeenth-Century/dp/0465018963?sa-no-redirect=1) about it.\n",
"\n",
"We can solve the problem with the tools we have:"
]
},
{
"cell_type": "code",
"execution_count": 49,
"metadata": {
"id": "FzdjxDZejrva"
},
"outputs": [],
"source": [
"def at_least(n, item) -> Event:\n",
" \"The event of getting at least n instances of item in an outcome.\"\n",
" return lambda outcome: outcome.count(item) >= n\n",
"\n",
"def all_finishes(tosses: int) -> Set[tuple]:\n",
" \"All finishes of a game with `tosses` more tosses.\"\n",
" return joins(product(*['ht'] * tosses))"
]
},
{
"cell_type": "markdown",
"metadata": {
"id": "Mib-C0yBjrva"
},
"source": [
"We can generate the 16 equiprobable 4-toss finishes that Pierre wrote about:"
]
},
{
"cell_type": "code",
"execution_count": 50,
"metadata": {
"id": "x3kP8N0Vjrva",
"outputId": "e14967c2-9a09-4d39-cd9e-af911aed12a6"
},
"outputs": [
{
"data": {
"text/plain": [
"{'h h h h',\n",
" 'h h h t',\n",
" 'h h t h',\n",
" 'h h t t',\n",
" 'h t h h',\n",
" 'h t h t',\n",
" 'h t t h',\n",
" 'h t t t',\n",
" 't h h h',\n",
" 't h h t',\n",
" 't h t h',\n",
" 't h t t',\n",
" 't t h h',\n",
" 't t h t',\n",
" 't t t h',\n",
" 't t t t'}"
]
},
"execution_count": 50,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"all_finishes(4)"
]
},
{
"cell_type": "markdown",
"metadata": {
"id": "sK2gLIBGjrva"
},
"source": [
"And we can find the 11 of them that are favorable to player `H`:"
]
},
{
"cell_type": "code",
"execution_count": 51,
"metadata": {
"id": "Cirm14DBjrva",
"outputId": "d87bc76e-a2a5-4110-991c-29843256059b"
},
"outputs": [
{
"data": {
"text/plain": [
"Dist({'h h t h': 1,\n",
" 'h t t h': 1,\n",
" 't h t h': 1,\n",
" 'h h h h': 1,\n",
" 'h h t t': 1,\n",
" 't t h h': 1,\n",
" 't h h t': 1,\n",
" 't h h h': 1,\n",
" 'h t h t': 1,\n",
" 'h t h h': 1,\n",
" 'h h h t': 1})"
]
},
"execution_count": 51,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"favorable(at_least(2, 'h'), all_finishes(4))"
]
},
{
"cell_type": "markdown",
"metadata": {
"collapsed": false,
"id": "s7RI7cp4jrva",
"jupyter": {
"outputs_hidden": false
}
},
"source": [
"Finally, we can answer the question:"
]
},
{
"cell_type": "code",
"execution_count": 52,
"metadata": {
"id": "gH5UuyDYjrva",
"outputId": "a7fd44c2-e7e1-47d0-e6ef-a9988e383f07"
},
"outputs": [
{
"data": {
"text/plain": [
"68.75"
]
},
"execution_count": 52,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"100 * P(at_least(2, 'h'), all_finishes(4))"
]
},
{
"cell_type": "markdown",
"metadata": {
"id": "iS97I1shjrva"
},
"source": [
"Blaise deserves 68.75 francs. We agree with Pascal and Fermat; we're in good company!"
]
},
{
"cell_type": "markdown",
"metadata": {
"id": "Xf-FTePajrva"
},
"source": [
"# Newton's Answer to a Problem by Pepys\n",
"\n",
"\n",
"\n",
"Let's jump ahead from 1654 all the way to 1693, [when](http://fermatslibrary.com/s/isaac-newton-as-a-probabilist) Samuel Pepys wrote to Isaac Newton posing the problem:\n",
"\n",
"> Which of the following three propositions has the greatest chance of success?\n",
" 1. Six fair dice are tossed independently and at least one “6” appears.\n",
" 2. Twelve fair dice are tossed independently and at least two “6”s appear.\n",
" 3. Eighteen fair dice are tossed independently and at least three “6”s appear.\n",
" \n",
"Newton was able to answer the question correctly (although his reasoning was not quite right); let's see how we can do. Since we're only interested in whether a die comes up as \"6\" or not, we can define a single die like this:"
]
},
{
"cell_type": "code",
"execution_count": 53,
"metadata": {
"id": "4FdHEVLVjrva"
},
"outputs": [],
"source": [
"die6 = Dist({'6': 1/6, '.': 5/6})"
]
},
{
"cell_type": "markdown",
"metadata": {
"id": "qLn1d8tUjrva"
},
"source": [
"Next we can define the joint distribution formed by combining two independent distribution like this:"
]
},
{
"cell_type": "code",
"execution_count": 54,
"metadata": {
"id": "UfhfuDQ6jrvb",
"outputId": "bfdb7d58-5ab8-4c1e-aacd-5503fcfa8ed7"
},
"outputs": [],
"source": [
"def joint2(A: Dist, B: Dist, format_fn='{}{}'.format) -> Dist:\n",
" \"\"\"The joint distribution of two independent distributions.\n",
" Result is all entries of the form {'ab': frequency(a) * frequency(b)}\"\"\"\n",
" return Dist({format_fn(a, b): A[a] * B[b]\n",
" for a in A for b in B})"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"(Note we pass in a format function to say how to combine the two outcomes. The default is to just concatenate them together, but it could be any function.)\n",
"\n",
"Here's the joint distribution of outcomes (not the sum) from rolling two dice:"
]
},
{
"cell_type": "code",
"execution_count": 55,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"Dist({'66': 0.027777777777777776,\n",
" '6.': 0.1388888888888889,\n",
" '.6': 0.1388888888888889,\n",
" '..': 0.6944444444444445})"
]
},
"execution_count": 55,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"joint2(die6, die6)"
]
},
{
"cell_type": "markdown",
"metadata": {
"collapsed": false,
"id": "asYQEoI-jrvb",
"jupyter": {
"outputs_hidden": false
}
},
"source": [
"And here is the joint distribution for `n` copies of the same distribution:"
]
},
{
"cell_type": "code",
"execution_count": 56,
"metadata": {
"id": "gDJ26cfGjrvb",
"outputId": "599ff7fe-e6e4-4dcf-eaa3-c8476001d50e"
},
"outputs": [],
"source": [
"def joint(n, dist: Dist, format_fn='{}{}'.format) -> Dist:\n",
" \"Joint probability distribution from rolling `n` dice.\"\n",
" if n == 1:\n",
" return dist\n",
" else:\n",
" return joint2(dist, joint(n - 1, dist, format_fn))"
]
},
{
"cell_type": "code",
"execution_count": 57,
"metadata": {
"id": "gDJ26cfGjrvb",
"outputId": "599ff7fe-e6e4-4dcf-eaa3-c8476001d50e"
},
"outputs": [
{
"data": {
"text/plain": [
"Dist({'6666': 0.0007716049382716049,\n",
" '666.': 0.0038580246913580245,\n",
" '66.6': 0.0038580246913580245,\n",
" '66..': 0.019290123456790126,\n",
" '6.66': 0.0038580246913580245,\n",
" '6.6.': 0.019290123456790126,\n",
" '6..6': 0.019290123456790126,\n",
" '6...': 0.09645061728395063,\n",
" '.666': 0.0038580246913580245,\n",
" '.66.': 0.019290123456790122,\n",
" '.6.6': 0.019290123456790122,\n",
" '.6..': 0.09645061728395063,\n",
" '..66': 0.019290123456790122,\n",
" '..6.': 0.09645061728395063,\n",
" '...6': 0.09645061728395063,\n",
" '....': 0.48225308641975323})"
]
},
"execution_count": 57,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"joint(4, die6)"
]
},
{
"cell_type": "markdown",
"metadata": {
"id": "xMgdWFuhjrvb"
},
"source": [
"Now we are ready to determine which proposition is more likely to have the required number of sixes:"
]
},
{
"cell_type": "code",
"execution_count": 58,
"metadata": {
"id": "70egSPGTjrvb",
"outputId": "ce5f10d8-2ae7-4a9f-aeeb-63bae635c5be"
},
"outputs": [
{
"data": {
"text/plain": [
"0.665102023319616"
]
},
"execution_count": 58,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"P(at_least(1, '6'), joint(6, die6))"
]
},
{
"cell_type": "code",
"execution_count": 59,
"metadata": {
"id": "wawDwjKWjrvb",
"outputId": "baf788bc-531e-4bf8-c1cd-37fee6820ca7"
},
"outputs": [
{
"data": {
"text/plain": [
"0.61866737373231"
]
},
"execution_count": 59,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"P(at_least(2, '6'), joint(12, die6))"
]
},
{
"cell_type": "code",
"execution_count": 60,
"metadata": {
"id": "H_RtixGLjrvb",
"outputId": "06eab287-0557-46fb-8814-b94a24ad6e45"
},
"outputs": [
{
"data": {
"text/plain": [
"0.5973456859477678"
]
},
"execution_count": 60,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"P(at_least(3, '6'), joint(18, die6))"
]
},
{
"cell_type": "markdown",
"metadata": {
"id": "ZC6eE6XYjrvb"
},
"source": [
"We reach the same conclusion Newton did, that the best chance is rolling six dice."
]
},
{
"cell_type": "markdown",
"metadata": {
"button": false,
"id": "khve5o3fjrvb",
"new_sheet": false,
"run_control": {
"read_only": false
}
},
"source": [
"# More Urn Problems: M&Ms and Bayes\n",
"\n",
"Here's another urn problem (actually a \"bag\" problem) [from](http://allendowney.blogspot.com/2011/10/my-favorite-bayess-theorem-problems.html) prolific Python/Probability pundit [Allen Downey ](http://allendowney.blogspot.com/):\n",
"\n",
"> The blue M&M was introduced in 1995. Before then, the color mix in a bag of plain M&Ms was (30% Brown, 20% Yellow, 20% Red, 10% Green, 10% Orange, 10% Tan). Afterward it was (24% Blue , 20% Green, 16% Orange, 14% Yellow, 13% Red, 13% Brown).\n",
"A friend of mine has two bags of M&Ms, and he tells me that one is from 1994 and one from 1996. He won't tell me which is which, but he randomly selects one M&M from each bag. One is yellow and one is green. What is the probability that the yellow M&M came from the 1994 bag?\n",
"\n",
"To solve this problem, we'll first create distributions for each bag: `bag94` and `bag96`:"
]
},
{
"cell_type": "code",
"execution_count": 61,
"metadata": {
"button": false,
"id": "4AapS83ujrvb",
"new_sheet": false,
"run_control": {
"read_only": false
}
},
"outputs": [],
"source": [
"bag94 = Dist(brown=30, yellow=20, red=20, green=10, orange=10, tan=10)\n",
"bag96 = Dist(blue=24, green=20, orange=16, yellow=14, red=13, brown=13)"
]
},
{
"cell_type": "markdown",
"metadata": {
"button": false,
"id": "MM6M8MN3jrvb",
"new_sheet": false,
"run_control": {
"read_only": false
}
},
"source": [
"Next, define `MM` as the joint distribution–the sample space for picking one M&M from each bag. The outcome `'94:yellow 96:green'` means that a yellow M&M was selected from the 1994 bag and a green one from the 1996 bag. In this problem we don't get to see the actual outcome; we just see some evidence about the outcome, that it contains a yellow and a green."
]
},
{
"cell_type": "code",
"execution_count": 62,
"metadata": {
"button": false,
"id": "XWiKZWKrjrvb",
"new_sheet": false,
"outputId": "db411dcf-9651-44f9-af86-f0331939a46e",
"run_control": {
"read_only": false
}
},
"outputs": [
{
"data": {
"text/plain": [
"Dist({'94:brown 96:blue': 720,\n",
" '94:brown 96:green': 600,\n",
" '94:brown 96:orange': 480,\n",
" '94:brown 96:yellow': 420,\n",
" '94:brown 96:red': 390,\n",
" '94:brown 96:brown': 390,\n",
" '94:yellow 96:blue': 480,\n",
" '94:yellow 96:green': 400,\n",
" '94:yellow 96:orange': 320,\n",
" '94:yellow 96:yellow': 280,\n",
" '94:yellow 96:red': 260,\n",
" '94:yellow 96:brown': 260,\n",
" '94:red 96:blue': 480,\n",
" '94:red 96:green': 400,\n",
" '94:red 96:orange': 320,\n",
" '94:red 96:yellow': 280,\n",
" '94:red 96:red': 260,\n",
" '94:red 96:brown': 260,\n",
" '94:green 96:blue': 240,\n",
" '94:green 96:green': 200,\n",
" '94:green 96:orange': 160,\n",
" '94:green 96:yellow': 140,\n",
" '94:green 96:red': 130,\n",
" '94:green 96:brown': 130,\n",
" '94:orange 96:blue': 240,\n",
" '94:orange 96:green': 200,\n",
" '94:orange 96:orange': 160,\n",
" '94:orange 96:yellow': 140,\n",
" '94:orange 96:red': 130,\n",
" '94:orange 96:brown': 130,\n",
" '94:tan 96:blue': 240,\n",
" '94:tan 96:green': 200,\n",
" '94:tan 96:orange': 160,\n",
" '94:tan 96:yellow': 140,\n",
" '94:tan 96:red': 130,\n",
" '94:tan 96:brown': 130})"
]
},
"execution_count": 62,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"MM = joint2(bag94, bag96, '94:{} 96:{}'.format)\n",
"MM"
]
},
{
"cell_type": "markdown",
"metadata": {
"button": false,
"id": "e80N-6s0jrvb",
"new_sheet": false,
"run_control": {
"read_only": false
}
},
"source": [
"We are given the observation that \"One is yellow and one is green\":"
]
},
{
"cell_type": "code",
"execution_count": 63,
"metadata": {
"button": false,
"id": "2GMvAOHQjrvb",
"new_sheet": false,
"outputId": "c6f461b8-0269-41bc-b922-d6485a87627a",
"run_control": {
"read_only": false
}
},
"outputs": [
{
"data": {
"text/plain": [
"Dist({'94:yellow 96:green': 400, '94:green 96:yellow': 140})"
]
},
"execution_count": 63,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"def yellow_and_green(outcome: str) -> bool: return 'yellow' in outcome and 'green' in outcome\n",
"\n",
"given(yellow_and_green, MM)"
]
},
{
"cell_type": "markdown",
"metadata": {
"id": "zb1IQU_ejrvb"
},
"source": [
"We want to know \"What is the probability that the yellow M&M came from the 1994 bag, given the observation?\""
]
},
{
"cell_type": "code",
"execution_count": 64,
"metadata": {
"button": false,
"id": "FCtLAPQfjrvb",
"new_sheet": false,
"outputId": "9c67a8e5-d521-4f39-de5c-54e893de388c",
"run_control": {
"read_only": false
}
},
"outputs": [
{
"data": {
"text/plain": [
"0.7407407407407407"
]
},
"execution_count": 64,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"def yellow94(outcome: str) -> bool: return '94:yellow' in outcome\n",
"\n",
"P(yellow94, given(yellow_and_green, MM))"
]
},
{
"cell_type": "markdown",
"metadata": {
"button": false,
"id": "SSxTlftyjrvb",
"new_sheet": false,
"run_control": {
"read_only": false
}
},
"source": [
"So there is a 74% chance that the yellow M&M comes from the 1994 bag.\n",
"\n",
"Answering this question was straightforward: just like all the other probability problems, we simply create a sample space, and use `P` to pick out the probability of the event in question, given what we know about the outcome.\n",
"\n",
"But in a sense it is curious that we were able to solve this problem with the same methodology as the others: this problem comes from a section of Downey's book titled **My favorite Bayes's Theorem Problems**, so one would expect that we'd need to invoke Bayes Theorem to solve it. The computation above shows that that is not necessary.\n",
"\n",
"| |\n",
"|---|\n",
"||\n",
"|Rev. Thomas Bayes (1701-1761)|\n",
"\n",
"Of course, we *could* solve it using Bayes Theorem. Why is Bayes Theorem recommended? Because we are asked about the probability of an outcome given the evidence—the probability the yellow came from the 94 bag, given that there is a yellow and a green. But the problem statement doesn't directly tell us the probability of that outcome given the evidence; it just tells us the probability of the evidence given the outcome.\n",
"\n",
"Before we see the colors of the M&Ms, there are two hypotheses, `A` and `B`, both with equal probability:\n",
"\n",
" A: first M&M from 94 bag, second from 96 bag\n",
" B: first M&M from 96 bag, second from 94 bag\n",
" P(A) = P(B) = 0.5\n",
" \n",
"Then we get some evidence:\n",
" \n",
" E: first M&M yellow, second green\n",
" \n",
"We want to know the probability of hypothesis `A`, given the evidence:\n",
" \n",
" P(A | E)\n",
" \n",
"That's not easy to calculate (except by enumerating the sample space, which our `P` function does). But Bayes Theorem says:\n",
" \n",
" P(A | E) = P(E | A) * P(A) / P(E)\n",
" \n",
"The quantities on the right-hand-side are easier to calculate:\n",
" \n",
" P(E | A) = 0.20 * 0.20 = 0.04\n",
" P(E | B) = 0.10 * 0.14 = 0.014\n",
" P(A) = 0.5\n",
" P(B) = 0.5\n",
" P(E) = P(E | A) * P(A) + P(E | B) * P(B)\n",
" = 0.04 * 0.5 + 0.014 * 0.5 = 0.027\n",
" \n",
"And we can get a final answer:\n",
" \n",
" P(A | E) = P(E | A) * P(A) / P(E)\n",
" = 0.04 * 0.5 / 0.027\n",
" = 0.7407407407\n",
" \n",
"You have a choice: Bayes Theorem allows you to do less calculation at the cost of more algebra; that is a great trade-off if you are working with pencil and paper. \n",
"\n",
"Enumerating the sample space allows you to do less algebra at the cost of more calculation; usually a good trade-off if you have a computer. But regardless of the approach you use, it is important to understand Bayes theorem and how it works.\n",
"\n",
"There is one important question that Allen Downey does not address: *would you eat thirty-year-old M&Ms*?\n",
"😨"
]
},
{
"cell_type": "markdown",
"metadata": {
"button": false,
"deletable": true,
"id": "sIT0znXcjrvd",
"new_sheet": false,
"run_control": {
"read_only": false
}
},
"source": [
"# Conclusion\n",
"\n",
"We can solve all these problems just by counting! \n",
"\n",
"All you ever needed to know about probability problems you learned from Sesame Street:\n",
"\n",
"| |\n",
"|---|\n",
"||\n",
"|The Count (1972—)|\n",
"\n",
"We've had an interesting tour of probability and met Laplace, Bernoulli, Fermat, Pascal, Bayes, Newton, and The Count.\n",
"\n",
"The conclusion is: be methodical in defining the sample space and the event(s) of interest, and be careful in counting the number of outcomes in the numerator and denominator. and you can't go wrong. Easy as 1-2-3."
]
}
],
"metadata": {
"colab": {
"provenance": []
},
"kernelspec": {
"display_name": "Python 3 (ipykernel)",
"language": "python",
"name": "python3"
},
"language_info": {
"codemirror_mode": {
"name": "ipython",
"version": 3
},
"file_extension": ".py",
"mimetype": "text/x-python",
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython3",
"version": "3.9.12"
}
},
"nbformat": 4,
"nbformat_minor": 4
}
![](\"https://upload.wikimedia.org/wikipedia/commons/f/f3/Pierre_de_Fermat.jpg\")
![](\"https://upload.wikimedia.org/wikipedia/commons/thumb/9/98/Blaise_Pascal_Versailles.JPG/440px-Blaise_Pascal_Versailles.JPG\")
![](\"http://scienceworld.wolfram.com/biography/pics/Newton.jpg\")
![](\"https://upload.wikimedia.org/wikipedia/commons/thumb/f/f8/Samuel_Pepys_portrait.jpg/148px-Samuel_Pepys_portrait.jpg\")