diff --git a/ipynb/war.ipynb b/ipynb/war.ipynb
index f688571..0d786c5 100644
--- a/ipynb/war.ipynb
+++ b/ipynb/war.ipynb
@@ -10,7 +10,7 @@
"\n",
"The [538 Riddler Classic for 28 August 2020](https://fivethirtyeight.com/features/can-you-cover-the-globe/) asks about the chance of winning the children's [card game **War**](https://en.wikipedia.org/wiki/War_%28card_game%29) in a **sweep**: a game where player **A** wins 26 turns in a row against player **B**. (In **War**, players are dealt 26 cards each and on each turn they flip the top card in their respective hands; higher rank card wins. Other rules are not relevant to this problem.)\n",
"\n",
- "We'll analyze this problem and come up with a program to solve it; first let's get the imports out of the way:"
+ "We'll analyze this problem and come up with a program to solve it. First let's get the imports out of the way:"
]
},
{
@@ -20,8 +20,8 @@
"outputs": [],
"source": [
"import random\n",
- "import itertools\n",
"import collections\n",
+ "from itertools import permutations, combinations\n",
"from statistics import mean"
]
},
@@ -29,46 +29,62 @@
"cell_type": "markdown",
"metadata": {},
"source": [
- "In my analysis I considered four different approaches to the problem, which I will describe in the order I considered them (although I actually implemented the final one first, and then came back to fill in the details of the earlier approaches).\n",
+ "Here are the four approaches I considered:\n",
"\n",
"# Approach 1: Simple Arithmetic?\n",
"\n",
- "A naive approach reasons as follows: there are 13 ranks, so perhaps the outcome probabilities for the first turn are:\n",
+ "There is no strategy in **War**, so on every turn, both players have an equal chance of winning the turn. If player **A** wins each turn with probability 1/2, than the probability of sweeping would be $(1/2)^{26}$ or about 15 in a billion. But that's **wrong** because there is also the possibility of a **tie**. After **A** picks a card, there are 51 cards remaining, and 3 have the same rank as **A**'s card, so the possible outcomes for the first turn are: \n",
"\n",
- " A wins=6/13; B wins=6/13; tie=1/13\n",
- " \n",
- "And thus the probability that **A** wins 26 turns in a row would be:"
- ]
- },
- {
- "cell_type": "code",
- "execution_count": 2,
- "metadata": {},
- "outputs": [
- {
- "data": {
- "text/plain": [
- "1.8595392516568175e-09"
- ]
- },
- "execution_count": 2,
- "metadata": {},
- "output_type": "execute_result"
- }
- ],
- "source": [
- "(6/13) ** 26"
+ " tie: 3/51\n",
+ " A wins: 24/51\n",
+ " B wins: 24/51\n",
+ " \n",
+ "If every subsequent turn were the same, the probability that **A** sweeps would be $(24/51)^{26}$ or about 3 in a billion. But that's not quite right because the outcome of each turn depends on the cards picked in the previous turns. So simple arithmetic isn't enough."
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
- "or about 2 in a billion. Unfortunately, that reasoning is **wrong**, even for one turn. The probability of a tie (that is, both players flipping cards of the same rank) on the first turn is actually 3/51 (not 1/13 or 4/52) because after player **A** flips a card, there are 51 equiprobable cards remaining for player **B**, and 3 have the same rank. So we actually have for the first turn:\n",
+ "# Approach 2: Brute Force Enumeration?\n",
"\n",
- " A wins=24/51; B wins=24/51; tie=3/51\n",
- " \n",
- "And if every subsequent turn were the same, the probability that **A** wins 26 turns in a row would be:"
+ "Brute force enumeration means:\n",
+ "- Consider every possible permutation of the deck of cards.\n",
+ "- For each permutation, deal out the cards and see whether or not **A** sweeps.\n",
+ "- The probability that **A** sweeps is the number of sweeps divided by the number of permutations.\n",
+ "\n",
+ "Easy-peasy; here's the code:"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 2,
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "Probability = float # Type: Probability is a number between 0.0 and 1.0\n",
+ "Deck = list # Type: Deck is a list of card ranks\n",
+ "\n",
+ "def make_deck(ranks=13, suits=4) -> Deck: \n",
+ " \"\"\"A Deck is a list of ranks, each rank repeated `suits` times.\"\"\"\n",
+ " return [r + 1 for r in range(ranks) for _ in range(suits)] \n",
+ "\n",
+ "def sweep(deck) -> bool: \n",
+ " \"\"\"Upon dealing this deck, does player A win every turn?\"\"\"\n",
+ " return all(deck[i] > deck[i + 1] for i in range(0, len(deck), 2))\n",
+ "\n",
+ "def p_sweeps(decks) -> Probability:\n",
+ " \"\"\"The probability that A sweeps, averaged over the decks.\"\"\"\n",
+ " return mean(sweep(deck) for deck in decks)"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "(*Implementation details:* in Python `False` is equivalent to `0` and `True` is equivalent to `1`, so the `mean` of an iterable of `bool` values is the same as the proportion (or probability) of `True` values. Also, it just felt wrong to have a card with rank `0`, thus I add 1 so that ranks start at `1`.)\n",
+ "\n",
+ "Here are two 8-card decks:"
]
},
{
@@ -79,7 +95,7 @@
{
"data": {
"text/plain": [
- "3.0808297965386556e-09"
+ "[1, 2, 3, 4, 5, 6, 7, 8]"
]
},
"execution_count": 3,
@@ -87,157 +103,22 @@
"output_type": "execute_result"
}
],
- "source": [
- "(24/51) ** 26"
- ]
- },
- {
- "cell_type": "markdown",
- "metadata": {},
- "source": [
- "But actually the probabilities on subsequent turns would be slightly different, depending on the cards picked in the previous turns. So simple arithmetic doesn't give us the answer.\n",
- "\n",
- "# Approach 2: Brute Force Enumeration?\n",
- "\n",
- "Brute force enumeration means:\n",
- "- Consider every possible permutation of the deck of cards.\n",
- "- For each permutation, deal out the cards and see whether or not **A** sweeps.\n",
- "- The probability that **A** sweeps is the number of sweeps divided by the number of permutations.\n",
- "\n",
- "Easy-peasy; here's the code. First the types:"
- ]
- },
- {
- "cell_type": "code",
- "execution_count": 4,
- "metadata": {},
- "outputs": [],
- "source": [
- "Probability = float # Type: Probability is a number between 0.0 and 1.0\n",
- "Deck = list # Type: Deck is a list of card ranks"
- ]
- },
- {
- "cell_type": "markdown",
- "metadata": {},
- "source": [
- "Now the functions:"
- ]
- },
- {
- "cell_type": "code",
- "execution_count": 5,
- "metadata": {},
- "outputs": [],
- "source": [
- "def make_deck(ranks=13, suits=4) -> Deck: \n",
- " \"\"\"Make a list of ranks, each rank repeated `suits` times.\"\"\"\n",
- " return list(range(ranks)) * suits\n",
- "\n",
- "def A_sweeps(deck) -> bool: \n",
- " \"\"\"Upon dealing this deck, does player A win every turn?\"\"\"\n",
- " return all(deck[i] > deck[i + 1] for i in range(0, len(deck), 2))\n",
- "\n",
- "def brute_force_p_sweep(deck) -> Probability:\n",
- " \"\"\"The probability that A sweeps, considering every permutation of deck.\"\"\"\n",
- " return mean(A_sweeps(d) for d in itertools.permutations(deck))"
- ]
- },
- {
- "cell_type": "markdown",
- "metadata": {},
- "source": [
- "(Note that in Python the `bool` value `False` is equivalent to `0` and `True` is equivalent to `1`, so the `mean` of an iterable of `bool` values is the same as the proportion or probability of `True` values.)\n",
- "\n",
- "Here we run the code on a tiny deck of just 8 cards, all with different ranks:"
- ]
- },
- {
- "cell_type": "code",
- "execution_count": 6,
- "metadata": {},
- "outputs": [
- {
- "data": {
- "text/plain": [
- "[0, 1, 2, 3, 4, 5, 6, 7]"
- ]
- },
- "execution_count": 6,
- "metadata": {},
- "output_type": "execute_result"
- }
- ],
"source": [
"make_deck(8, 1)"
]
},
{
"cell_type": "code",
- "execution_count": 7,
+ "execution_count": 4,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
- "0.0625"
+ "[1, 1, 2, 2, 3, 3, 4, 4]"
]
},
- "execution_count": 7,
- "metadata": {},
- "output_type": "execute_result"
- }
- ],
- "source": [
- "brute_force_p_sweep(make_deck(8, 1))"
- ]
- },
- {
- "cell_type": "markdown",
- "metadata": {},
- "source": [
- "There can be no ties, so this is four turns where on each turn **A** and **B** have equal chances of winning and the probability of **A** sweeping is: "
- ]
- },
- {
- "cell_type": "code",
- "execution_count": 8,
- "metadata": {},
- "outputs": [
- {
- "data": {
- "text/plain": [
- "0.0625"
- ]
- },
- "execution_count": 8,
- "metadata": {},
- "output_type": "execute_result"
- }
- ],
- "source": [
- "(1/2) ** 4"
- ]
- },
- {
- "cell_type": "markdown",
- "metadata": {},
- "source": [
- "If the 8-card deck had 4 ranks and 2 suits, then the probability of **A** sweeping would be less, because there could be ties:"
- ]
- },
- {
- "cell_type": "code",
- "execution_count": 9,
- "metadata": {},
- "outputs": [
- {
- "data": {
- "text/plain": [
- "[0, 1, 2, 3, 0, 1, 2, 3]"
- ]
- },
- "execution_count": 9,
+ "execution_count": 4,
"metadata": {},
"output_type": "execute_result"
}
@@ -246,9 +127,36 @@
"make_deck(4, 2)"
]
},
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "And here are the probabilities of sweeping a game played with those decks:"
+ ]
+ },
{
"cell_type": "code",
- "execution_count": 10,
+ "execution_count": 5,
+ "metadata": {},
+ "outputs": [
+ {
+ "data": {
+ "text/plain": [
+ "0.0625"
+ ]
+ },
+ "execution_count": 5,
+ "metadata": {},
+ "output_type": "execute_result"
+ }
+ ],
+ "source": [
+ "p_sweeps(permutations(make_deck(8, 1)))"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 6,
"metadata": {},
"outputs": [
{
@@ -257,13 +165,13 @@
"0.03571428571428571"
]
},
- "execution_count": 10,
+ "execution_count": 6,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
- "brute_force_p_sweep(make_deck(4, 2))"
+ "p_sweeps(permutations(make_deck(4, 2)))"
]
},
{
@@ -281,29 +189,25 @@
"source": [
"# Approach 3: Simulation?\n",
"\n",
- "It would take too long to look at **all** the permutations, but we can **randomly sample** the space of possible permutations; we call that a **simulation**:"
+ "It takes too long to look at **all** the permutations, but we can **randomly sample** the permutations. We call that a **simulation**:"
]
},
{
"cell_type": "code",
- "execution_count": 11,
+ "execution_count": 7,
"metadata": {},
"outputs": [],
"source": [
- "def deals(deck, N) -> Deck: \n",
+ "def sample(deck, N) -> Deck: \n",
" \"\"\"Yield N randomly shuffled deals of deck.\"\"\"\n",
" for _ in range(N):\n",
" random.shuffle(deck)\n",
- " yield deck\n",
- "\n",
- "def simulate(deck, N=10000) -> Probability:\n",
- " \"\"\"Simulate N games of War, and return the estimated probability of A sweeping.\"\"\"\n",
- " return mean(A_sweeps(d) for d in deals(deck, N))"
+ " yield deck"
]
},
{
"cell_type": "code",
- "execution_count": 12,
+ "execution_count": 8,
"metadata": {},
"outputs": [
{
@@ -312,20 +216,20 @@
"0"
]
},
- "execution_count": 12,
+ "execution_count": 8,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
- "simulate(make_deck())"
+ "p_sweeps(sample(make_deck(), 10_000))"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
- "Sampling 10,000 deals with the full 52-card deck we got zero sweeps. Since the probability of a sweep is probably somewhere around one in a billion, we would need to look at trillions of deals to get a reliable estimate. That would require roughly a day of run time: much better than millions of years, but still not good enough."
+ "Sampling 10,000 deals with the full 52-card deck we got zero sweeps. Approach 1 suggests the probability of a sweep is somewhere near 3 in a billion, which means we would need to sample trillions of deals to get a reliable estimate of the probability. That would require roughly a day of run time: much better than millions of years, but still not good enough."
]
},
{
@@ -334,93 +238,76 @@
"source": [
"# Approach 4: Abstract Incremental Enumeration!\n",
"\n",
- "The first three approaches didn't pan out. What next? \n",
- "\n",
- "It is feasible to consider every possible deal if we are careful. As discussed in my [How To Count Things](How%20To%20Count%20Things.ipynb) notebook, the idea is to represent a deck not as a concrete permutation of 52 cards but rather with a representation that is:\n",
- "- **Abstract**: What really matters is not the exact rank of every card in the deck, but rather whether **A**'s next card is higher, lower, or the same as **B**'s next card. For example, if there are only two cards remaining and we know they have different ranks, then the probability of **A** winning is 1/2; it doesn't matter whether the cards are [10, 8] or [2, 5] or any of the 52 × 51 possibilities for two cards.\n",
- "- **Incremental**: First we'll consider the possibilities for the two cards in the first turn, and only if **A** wins will we then move on to consider possible cards for the next turn. For cases where **A** loses or ties the first turn, there is no need to consider the 50! permutations of the remaining cards.\n",
"\n",
"\n",
- "The key to the abstract representation is knowing the cards of the same rank. I will define an abstract deck (or `AbDeck`) as a tuple where `deck[i]` gives the number of different ranks that have exactly `i` cards remaining in the deck. That's all the information you need to determine the probability of winning. (This representation is a bit wasteful, because `deck[0]` is always a redundant `0`, but the alternative would be to sprinkle the code with `[i - 1]` indexes, risking an off-by-one error).\n",
+ "As discussed in my [**How To Count Things**](How%20To%20Count%20Things.ipynb) notebook, in problems where brute force enumeration is not feasible it is often possible to consider every outcome if we use a representation that is:\n",
+ "- **Abstract**: What matters is not the exact rank of every card in the deck, but rather whether **A**'s next card is higher, lower, or the same as **B**'s next card. The representation should tell us enough to compute the odds of those three outcomes, and no more.\n",
+ "- **Incremental**: First we'll consider the possibilities for the two cards in the first turn, and only for the outcomes in which **A** wins will we then move on to consider possible cards for the next turn. For outcomes in which **A** loses or ties, there is no need to consider the permutations of the remaining cards.\n",
"\n",
- "Consider these sample abstract decks:"
+ "For example, if there are two cards remaining then there are 52 × 51 possible concrete decks, but only two abstract cases: \n",
+ "* The two cards are **the same** rank, in which case the probability of **A** winning is 0, because it will be a tie.
It doesn't matter whether the deck is `[5, 5]` or `[9, 9]`, or any other pair of cards. \n",
+ "* The two cards are **different** ranks, in which case the probability of **A** winning is 1/2.
It doesn't matter whether the deck is `[10, 8]` or `[2, 5]` or any other two distinct ranks.\n",
+ " \n",
+ "# Abstract Deck\n",
+ "\n",
+ "The representation I will use, which I call `AbDeck` for **abstract deck**, is as follows:\n",
+ "- An `AbDeck` is a tuple of counts of the number of cards of each rank, without specifying the rank. \n",
+ "- Ranks with counts of zero are dropped from the representation.\n",
+ "- The counts are put in sorted order.\n",
+ "- The tuple represents the deck without regard to the order of the cards (that is, not a specific permutation).\n",
+ "- The tuple tells you just what you need (no more, no less) to compute the probability of a sweep.\n",
+ "- The initial 52-card deck is 13 fours: `(4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4)`\n",
+ "- The possible two-card decks are `(1, 1)` and `(2,)`.\n",
+ "- The possible four-card decks are `(1, 1, 1, 1)`, `(1, 1, 2)`, `(1, 3)`, `(2, 2)`, and `(4,)`.\n",
+ "\n"
]
},
{
"cell_type": "code",
- "execution_count": 13,
- "metadata": {},
- "outputs": [],
- "source": [
- "AbDeck = tuple # The Abstract Deck Type\n",
- "\n",
- "deck = AbDeck((0, 0, 0, 0, 13))\n",
- "tie1 = AbDeck((0, 0, 1, 0, 12))\n",
- "dif1 = AbDeck((0, 0, 0, 2, 11))\n",
- "end4 = AbDeck((0, 4))"
- ]
- },
- {
- "cell_type": "markdown",
- "metadata": {},
- "source": [
- "- `deck` is the normal 52 card full deck with 13 ranks, each with 4 cards of that rank, so `deck[4]` is 13.\n",
- "- `tie1` results when both players flip a card of the same rank on the first turn.
There are 12 ranks with 4 cards each and 1 rank with 2 cards.\n",
- "- `dif1` results when the players flip cards of different ranks on the first turn.
There are 11 ranks each with 4 cards and 2 ranks with 3 cards. \n",
- "- `end4` is a deck you would see near the end of a game, with just 4 remaining cards, each of a different rank. \n",
- "\n",
- "Note that `tie1` and `dif1` are the only two possible abstract deck outcomes for turn 1. That's a lot better than the 52! possibilities that we would need to consider on turn 1 with a concrete deck. \n",
- "\n",
- "We can create abstract decks and count the number of cards as follows:"
- ]
- },
- {
- "cell_type": "code",
- "execution_count": 14,
- "metadata": {},
- "outputs": [],
- "source": [
- "def make_abdeck(ranks=13, suits=4) -> AbDeck: \n",
- " \"\"\"Make an abstract deck.\"\"\"\n",
- " return (0,) * suits + (ranks,)\n",
- " \n",
- "def ncards(deck: AbDeck) -> int:\n",
- " \"\"\"Number of cards remaining in an abstract deck.\"\"\"\n",
- " return sum(i * deck[i] for i in range(len(deck)))"
- ]
- },
- {
- "cell_type": "code",
- "execution_count": 15,
+ "execution_count": 9,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
- "[52, 50, 50, 4]"
+ "(4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4)"
]
},
- "execution_count": 15,
+ "execution_count": 9,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
- "[ncards(d) for d in (deck, tie1, dif1, end4)]"
+ "AbDeck = tuple # Type: An abstract deck\n",
+ "\n",
+ "def make_abdeck(ranks=13, suits=4) -> AbDeck: \n",
+ " \"\"\"Make an abstract deck.\"\"\"\n",
+ " return AbDeck([suits] * ranks)\n",
+ "\n",
+ "deck = make_abdeck(13, 4) # The initial 52-card abstract deck\n",
+ "deck"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
- "# Probability Distributions\n",
+ ">***Note:*** I originally thought of *two* possible abstract representations for decks: `AbDeck` as described above, and an alternative in which: \n",
+ "- A deck was a tuple where `deck[i]` gives the number of different ranks that have exactly `i` cards remaining in the deck.\n",
+ "- The initial deck is represented as `(0, 0, 0, 0, 13)`: 13 ranks with 4-of-a-kind.\n",
+ "- The two possible two-card decks are `(0, 2)` (two ranks with 1 card each) and `(0, 0, 1)` (one rank with two cards).\n",
+ ">\n",
+ ">That representation was more compact, and ran faster; I used it in the first version of this notebook. But the very nice [implementation](https://laurentlessard.com/bookproofs/flawless-war/) shared by [Laurent Lessard](https://laurentlessard.com/) convinced me to switch to the `AbDeck` representation becausse it is simpler to code and understand.\n",
"\n",
- "We want to keep track of the probability of **A** winning, over a variety of possible events (flips of cards). We'll do that with the help of a class called `PDist` (or **probability distribution**), a mapping that lists all the possible abstract decks that might occur on a turn, each mapped to their probability of occurrance."
+ "# Probability Distribution\n",
+ "\n",
+ "We want to keep track of the possible outcomes of each turn. We'll do that with the help of a class called `PDist` (for **probability distribution**), a mapping of `{deck: p, ...}` where `deck` is an abstract deck that is the outcome after some number of turns, each of which **A** has won, and `p` is the probability of that outcome."
]
},
{
"cell_type": "code",
- "execution_count": 16,
+ "execution_count": 10,
"metadata": {},
"outputs": [],
"source": [
@@ -434,82 +321,46 @@
"source": [
"# Abstract Incremental Enumeration Strategy\n",
"\n",
- "Now we're ready to outline a strategy to exactly and efficiently compute `p_sweep(deck)`, the probability that **A** sweeps a game of **War**:\n",
+ "The probability that **A** sweeps a game of **War** given an abstract deck, `p_sweep(deck)`, is computed as follows:\n",
"\n",
- "- Start with `P` being a **probability distribution** of deck outcomes after 0 turns: `{deck: 1.0}`.\n",
- "- **for** each of the 26 turns (or in general ncards(deck) / 2 turns):\n",
- " - Update `P` to be the result of playing a turn, which is given by:\n",
- " - **for** each of the entries in `P`:\n",
- " - See what possible outcomes can arise from selecting 2 cards from the deck\n",
- " - Update the new `PDist` `P1` to include each outcome where **A** might sweep, with appropriate probability\n",
- "- Sum the probabilities in `P`"
+ "- Start with `P` being a **probability distribution** of outcomes after 0 turns: `{deck: 1.0}`.\n",
+ "- **for** each turn of the game (26 turns for a 52-card deck):\n",
+ " - Update `P` to be the probability distribution that results from playing a turn, `P = play_turn(P)`:\n",
+ " - **for** each `deck` in `P`:\n",
+ " - **for** each way of removing two cards of rank `i` and rank `j` from `deck` to yield a resulting deck:\n",
+ " - Increment the probability of the resulting deck by (probability of `deck` × probability of `i` × probability of `j`).\n",
+ "- **return** the sum of the probabilities of the winning decks in `P` from the final turn\n",
+ "\n",
+ "Note that (for example) `combinations(range(3), 2)` yields `(1, 2)`, but not `(1, 1)` nor `(2, 1)`. That's just what we want: if the two selected cards are `i, j = (1, 1)` then the result will be a tie, which we don't want to put in the probability distribution of decks where **A** sweeps. And out of the two possibilities `i, j = (1, 2)` and `i, j = (2, 1)`, exactly one will result in a win for **A**, so we only want to include one of them."
]
},
{
"cell_type": "code",
- "execution_count": 17,
+ "execution_count": 11,
"metadata": {},
"outputs": [],
"source": [
"def p_sweep(deck: AbDeck) -> Probability:\n",
" \"\"\"The probability that player A sweeps a game of War.\"\"\"\n",
" P = PDist({deck: 1.0}) # The initial probability distribution\n",
- " for turn in range(ncards(deck) // 2):\n",
+ " for turn in range(sum(deck) // 2):\n",
" P = play_turn(P)\n",
" return sum(P.values())\n",
"\n",
"def play_turn(P) -> PDist:\n",
" \"\"\"Play one turn with all possible card choices for players A and B; return\n",
" the probability distribution of outcomes where A still might sweep.\"\"\"\n",
- " P1 = PDist() # The probability distribution after this turn\n",
+ " P1 = PDist() # The probability distribution of {deck: prob} after this turn\n",
" for deck in P:\n",
- " for deck2, p_Awin in select2(deck):\n",
- " P1[deck2] += P[deck] * p_Awin\n",
- " return P1"
- ]
- },
- {
- "cell_type": "markdown",
- "metadata": {},
- "source": [
- "Now the key is figuring out:\n",
- "- All the possible ways to select two cards from the deck.\n",
- "- The probability of each selection.\n",
- "- Who won (or was it a tie) with that selection."
- ]
- },
- {
- "cell_type": "code",
- "execution_count": 18,
- "metadata": {},
- "outputs": [],
- "source": [
- "def select1(deck) -> list:\n",
- " \"\"\"All ways to pick one card from deck, returning a list of tuples:\n",
- " (remaining deck, probability of pick, index of pick in deck)\"\"\"\n",
- " return [(remove(i, deck), i * deck[i] / ncards(deck), i)\n",
- " for i in range(len(deck)) if deck[i]]\n",
+ " n = sum(deck) # number of cards in deck\n",
+ " for i, j in combinations(range(len(deck)), 2):\n",
+ " P1[remove(deck, i, j)] += P[deck] * deck[i] / n * deck[j] / (n - 1)\n",
+ " return P1\n",
"\n",
- "def select2(deck) -> list:\n",
- " \"\"\"All ways to select two cards from deck, returning a list of tuples:\n",
- " (remaining deck, probability that deck occurs and player A wins).\"\"\"\n",
- " return [(deck2, p_i * p_j * p_Awin(deck1, i, j))\n",
- " for deck1, p_i, i in select1(deck)\n",
- " for deck2, p_j, j in select1(deck1)\n",
- " if p_Awin(deck1, i, j) > 0]\n",
- " \n",
- "def p_Awin(deck1, i, j) -> Probability:\n",
- " \"\"\"Probability that A wins turn when A selects i giving deck1 and B selects j.\"\"\"\n",
- " p_tie = 1 / deck1[j] if j == i - 1 else 0\n",
- " return (1 - p_tie) / 2\n",
- "\n",
- "def remove(i, deck) -> AbDeck:\n",
- " \"\"\"Remove one card from deck[i].\"\"\"\n",
- " deck1 = list(deck)\n",
- " deck1[i] -= 1\n",
- " if i - 1 != 0:\n",
- " deck1[i - 1] += 1\n",
- " return AbDeck(deck1)"
+ "def remove(deck, *indexes) -> AbDeck:\n",
+ " \"\"\"Remove cards at indexes from deck.\"\"\"\n",
+ " counts = [c - indexes.count(i) for i, c in enumerate(deck)]\n",
+ " return AbDeck(sorted(c for c in counts if c > 0))"
]
},
{
@@ -523,16 +374,16 @@
},
{
"cell_type": "code",
- "execution_count": 19,
+ "execution_count": 12,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
- "3.132436174322294e-09"
+ "3.132436174322299e-09"
]
},
- "execution_count": 19,
+ "execution_count": 12,
"metadata": {},
"output_type": "execute_result"
}
@@ -545,25 +396,25 @@
"cell_type": "markdown",
"metadata": {},
"source": [
- "The probability that **A** sweeps a game of **War** is a little over 3 in a billion. \n",
+ "A little over 3 in a billion. \n",
"\n",
- "(By the way, this computation took less than 1/10 second, a big improvement over millions of years.) \n",
+ "(By the way, this computation took about 0.2 seconds, a big improvement over millions of years.) \n",
"\n",
"That's the answer to *my* question about the probability of **A** sweeping, but 538 Riddler actually posed the question somewhat differently: \"*How many games of War would you expect to play until you had a game that lasted just 26 turns with no wars, like Duane’s friend’s granddaughter?*\" That is, they want to know the inverse of the probability, and they are allowing for either **A** or **B** to sweep. So that would be:"
]
},
{
"cell_type": "code",
- "execution_count": 20,
+ "execution_count": 13,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
- "159620171.7049113"
+ "159620171.70491105"
]
},
- "execution_count": 20,
+ "execution_count": 13,
"metadata": {},
"output_type": "execute_result"
}
@@ -585,99 +436,21 @@
"source": [
"# Working through the algorithm\n",
"\n",
- "Let's work through how `p_sweep`, `play_turn`, `select1` and `select2` work. We'll start by reminding ourself what the starting `deck` is, and then `select1` card from it:"
+ "Let's work through how `p_sweep` and `play_turn` work. We'll start with `P` being the initial probability distribution, and then play four turns, each time displaying how `P` is updated:"
]
},
{
"cell_type": "code",
- "execution_count": 21,
+ "execution_count": 14,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
- "(0, 0, 0, 0, 13)"
+ "PDist({(4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4): 1.0})"
]
},
- "execution_count": 21,
- "metadata": {},
- "output_type": "execute_result"
- }
- ],
- "source": [
- "deck"
- ]
- },
- {
- "cell_type": "code",
- "execution_count": 22,
- "metadata": {},
- "outputs": [
- {
- "data": {
- "text/plain": [
- "[((0, 0, 0, 1, 12), 1.0, 4)]"
- ]
- },
- "execution_count": 22,
- "metadata": {},
- "output_type": "execute_result"
- }
- ],
- "source": [
- "select1(deck)"
- ]
- },
- {
- "cell_type": "markdown",
- "metadata": {},
- "source": [
- "This is saying with probability 1.0 the result of selecting one card is a deck where there are 12 ranks with four-of-a-kind and 1 rank with three-of-a-kind. The selected card came from index 4 in the deck (it was one of a four-of-a-kind).\n",
- "\n",
- "Now we'll try selecting two cards:"
- ]
- },
- {
- "cell_type": "code",
- "execution_count": 23,
- "metadata": {},
- "outputs": [
- {
- "data": {
- "text/plain": [
- "[((0, 0, 0, 2, 11), 0.47058823529411764)]"
- ]
- },
- "execution_count": 23,
- "metadata": {},
- "output_type": "execute_result"
- }
- ],
- "source": [
- "select2(deck)"
- ]
- },
- {
- "cell_type": "markdown",
- "metadata": {},
- "source": [
- "This says that the only result of the first turn in which **A** wins has 11 ranks with four-of-a-kind and 2 ranks with 3-of-a-kind. The probability of this outcome is 0.47058823529411764, which, as we calculated earlier, is 24/51. The rest of the probability goes to an equiprobable result in which **B** wins, and to `(0, 0, 1, 0, 12)`, which indicates a tie on the first turn. Since these outcomes don't result in **A** winning, they do not appear in the result from `select2`.\n",
- "\n",
- "Now let's work through how `p_sweep` repeatedly calls `play_turn`, updating the `PDist` `P`, which is initially:"
- ]
- },
- {
- "cell_type": "code",
- "execution_count": 24,
- "metadata": {},
- "outputs": [
- {
- "data": {
- "text/plain": [
- "PDist({(0, 0, 0, 0, 13): 1.0})"
- ]
- },
- "execution_count": 24,
+ "execution_count": 14,
"metadata": {},
"output_type": "execute_result"
}
@@ -687,31 +460,24 @@
"P"
]
},
- {
- "cell_type": "markdown",
- "metadata": {},
- "source": [
- "Here's the outcome of playing the first turn:"
- ]
- },
{
"cell_type": "code",
- "execution_count": 25,
+ "execution_count": 15,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
- "PDist({(0, 0, 0, 2, 11): 0.47058823529411764})"
+ "PDist({(3, 3, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4): 0.47058823529411703})"
]
},
- "execution_count": 25,
+ "execution_count": 15,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
- "P = play_turn(P)\n",
+ "P = play_turn(P) # Turn 1\n",
"P"
]
},
@@ -719,29 +485,91 @@
"cell_type": "markdown",
"metadata": {},
"source": [
- "Now for the second turn:"
+ "That makes sense: the only result where **A** still has a chance to sweep is when the two players select different ranks, giving us two ranks with three-of-a-kind remaining (and 11 with four-of-a-kind remaining). The probability of **A** sweeping one turn is 24/51 = 0.47058823529411703."
]
},
{
"cell_type": "code",
- "execution_count": 26,
+ "execution_count": 16,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
- "PDist({(0, 0, 2, 0, 11): 0.0017286914765906362,\n",
- " (0, 0, 1, 2, 10): 0.05070828331332534,\n",
- " (0, 0, 0, 4, 9): 0.16902761104441777})"
+ "PDist({(2, 2, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4): 0.001728691476590634,\n",
+ " (2, 3, 3, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4): 0.050708283313325254,\n",
+ " (3, 3, 3, 3, 4, 4, 4, 4, 4, 4, 4, 4, 4): 0.16902761104441777})"
]
},
- "execution_count": 26,
+ "execution_count": 16,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
- "P = play_turn(P)\n",
+ "P = play_turn(P) # Turn 2\n",
+ "P"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 17,
+ "metadata": {},
+ "outputs": [
+ {
+ "data": {
+ "text/plain": [
+ "PDist({(1, 1, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4): 3.065055809557862e-06,\n",
+ " (1, 2, 3, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4): 0.00040458736686163773,\n",
+ " (2, 2, 3, 3, 4, 4, 4, 4, 4, 4, 4, 4, 4): 0.010114684171540949,\n",
+ " (1, 3, 3, 3, 4, 4, 4, 4, 4, 4, 4, 4, 4): 0.0017981660749406122,\n",
+ " (2, 2, 2, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4): 0.00020229368343081884,\n",
+ " (2, 3, 3, 3, 3, 4, 4, 4, 4, 4, 4, 4, 4): 0.04855048402339655,\n",
+ " (3, 3, 3, 3, 3, 3, 4, 4, 4, 4, 4, 4, 4): 0.043155985798574735})"
+ ]
+ },
+ "execution_count": 17,
+ "metadata": {},
+ "output_type": "execute_result"
+ }
+ ],
+ "source": [
+ "P = play_turn(P) # Turn 3\n",
+ "P"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 18,
+ "metadata": {},
+ "outputs": [
+ {
+ "data": {
+ "text/plain": [
+ "PDist({(4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4): 1.4807032896414793e-09,\n",
+ " (1, 3, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4): 5.212075579538007e-07,\n",
+ " (1, 1, 3, 3, 4, 4, 4, 4, 4, 4, 4, 4, 4): 3.648452905676605e-05,\n",
+ " (2, 2, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4): 5.863585026980257e-07,\n",
+ " (2, 3, 3, 4, 4, 4, 4, 4, 4, 4, 4, 4): 1.563622673861402e-05,\n",
+ " (1, 1, 2, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4): 2.3454340107921024e-06,\n",
+ " (1, 2, 2, 3, 4, 4, 4, 4, 4, 4, 4, 4, 4): 0.00018763472086336856,\n",
+ " (1, 2, 3, 3, 3, 4, 4, 4, 4, 4, 4, 4, 4): 0.0016887124877703135,\n",
+ " (2, 2, 2, 2, 4, 4, 4, 4, 4, 4, 4, 4, 4): 4.397688770235195e-05,\n",
+ " (2, 2, 2, 3, 3, 4, 4, 4, 4, 4, 4, 4, 4): 0.0023923426910079496,\n",
+ " (2, 2, 3, 3, 3, 3, 4, 4, 4, 4, 4, 4, 4): 0.014635508227342713,\n",
+ " (3, 3, 3, 3, 4, 4, 4, 4, 4, 4, 4, 4): 3.127245347722804e-05,\n",
+ " (1, 3, 3, 3, 3, 3, 4, 4, 4, 4, 4, 4, 4): 0.002001437022542595,\n",
+ " (2, 3, 3, 3, 3, 3, 3, 4, 4, 4, 4, 4, 4): 0.02101508873669728,\n",
+ " (3, 3, 3, 3, 3, 3, 3, 3, 4, 4, 4, 4, 4): 0.007005029578899085})"
+ ]
+ },
+ "execution_count": 18,
+ "metadata": {},
+ "output_type": "execute_result"
+ }
+ ],
+ "source": [
+ "P = play_turn(P) # Turn 4\n",
"P"
]
},
@@ -749,68 +577,7 @@
"cell_type": "markdown",
"metadata": {},
"source": [
- "And the third turn:"
- ]
- },
- {
- "cell_type": "code",
- "execution_count": 27,
- "metadata": {},
- "outputs": [
- {
- "data": {
- "text/plain": [
- "PDist({(0, 2, 0, 0, 11): 3.0650558095578654e-06,\n",
- " (0, 1, 1, 1, 10): 0.00040458736686163827,\n",
- " (0, 0, 2, 2, 9): 0.010114684171540957,\n",
- " (0, 1, 0, 3, 9): 0.0017981660749406148,\n",
- " (0, 0, 3, 0, 10): 0.00020229368343081916,\n",
- " (0, 0, 1, 4, 8): 0.048550484023396595,\n",
- " (0, 0, 0, 6, 7): 0.04315598579857475})"
- ]
- },
- "execution_count": 27,
- "metadata": {},
- "output_type": "execute_result"
- }
- ],
- "source": [
- "P = play_turn(P)\n",
- "P"
- ]
- },
- {
- "cell_type": "markdown",
- "metadata": {},
- "source": [
- "We'll leave it as an exercise for the reader to work through these, but one thing you can clearly see is that the total probability of winning is becoming smaller with every turn:"
- ]
- },
- {
- "cell_type": "code",
- "execution_count": 28,
- "metadata": {},
- "outputs": [
- {
- "data": {
- "text/plain": [
- "0.10422926617455494"
- ]
- },
- "execution_count": 28,
- "metadata": {},
- "output_type": "execute_result"
- }
- ],
- "source": [
- "sum(P.values())"
- ]
- },
- {
- "cell_type": "markdown",
- "metadata": {},
- "source": [
- "There is only a 10% chance of winning three turns in a row. This will steadily fall to 3-in-a-billion after 26 turns. "
+ "We'll leave it as an exercise for the reader to verify that these turns are correct."
]
},
{
@@ -819,21 +586,21 @@
"source": [
"# Gaining Confidence in the Answer\n",
"\n",
- "The answer should be somewhere close to the arithmetic calculation of $(24/51)^{26}$; in fact we see that it differs by less than 2%:"
+ "One way to gain confidence: The answer should be somewhere close to the arithmetic calculation of $(24/51)^{26}$."
]
},
{
"cell_type": "code",
- "execution_count": 29,
+ "execution_count": 19,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
- "1.0167508045532467"
+ "1.0167508045532485"
]
},
- "execution_count": 29,
+ "execution_count": 19,
"metadata": {},
"output_type": "execute_result"
}
@@ -846,83 +613,43 @@
"cell_type": "markdown",
"metadata": {},
"source": [
- "Although the brute force algorithm can't handle a 52-card deck, we can use it to verify that there is no difference between `brute_force_p_sweep(deck)` and `p_sweep(deck)` for small decks; this gives us more confidence that both are correct (or maybe both have similar errors):"
+ "We see that it is indeed close, differing by less than 2%.\n",
+ "\n",
+ "Another way to gain confidence: We can verify that there is no difference (∆ = 0) between the brute force `p_sweeps` and the abstract incremental `p_sweep` for small decks. This gives us some confidence that both are correct (or maybe both have similar errors). `p_sweeps` takes less than 50 milliseconds to solve decks of 8 cards or less, but takes about 4 seconds to do the 10! (almost 4 million) permutations of ten card decks."
]
},
{
"cell_type": "code",
- "execution_count": 30,
+ "execution_count": 20,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
- "deck(2, 1): p = 0.50000; diff = 0.000000000\n",
- "deck(4, 1): p = 0.25000; diff = 0.000000000\n",
- "deck(6, 1): p = 0.12500; diff = 0.000000000\n",
- "deck(8, 1): p = 0.06250; diff = 0.000000000\n",
- "deck(2, 2): p = 0.16667; diff = 0.000000000\n",
- "deck(2, 3): p = 0.05000; diff = 0.000000000\n",
- "deck(3, 2): p = 0.06667; diff = 0.000000000\n",
- "deck(2, 4): p = 0.01429; diff = 0.000000000\n",
- "deck(4, 2): p = 0.03571; diff = 0.000000000\n"
+ "deck(2, 1): p = 0.5000; ∆ = 0.000000000\n",
+ "deck(4, 1): p = 0.2500; ∆ = 0.000000000\n",
+ "deck(2, 2): p = 0.1667; ∆ = 0.000000000\n",
+ "deck(6, 1): p = 0.1250; ∆ = 0.000000000\n",
+ "deck(3, 2): p = 0.0667; ∆ = 0.000000000\n",
+ "deck(2, 3): p = 0.0500; ∆ = 0.000000000\n",
+ "deck(8, 1): p = 0.0625; ∆ = 0.000000000\n",
+ "deck(4, 2): p = 0.0357; ∆ = 0.000000000\n",
+ "deck(2, 4): p = 0.0143; ∆ = 0.000000000\n",
+ "deck(5, 2): p = 0.0180; ∆ = 0.000000000\n",
+ "deck(2, 5): p = 0.0040; ∆ = 0.000000000\n"
]
}
],
"source": [
- "for ranks, suits in [(2, 1), (4, 1), (6, 1), (8, 1),\n",
- " (2, 2), (2, 3), (3, 2), (2, 4), (4, 2)]:\n",
- " p1 = brute_force_p_sweep(make_deck(ranks, suits))\n",
- " p2 = p_sweep(make_abdeck(ranks, suits))\n",
- " print(f'deck({ranks}, {suits}): p = {p1:.5f}; diff = {abs(p1 - p2):.9f}')"
- ]
- },
- {
- "cell_type": "markdown",
- "metadata": {},
- "source": [
- "And we can do some unit tests on components (we should do more):"
- ]
- },
- {
- "cell_type": "code",
- "execution_count": 31,
- "metadata": {},
- "outputs": [
- {
- "data": {
- "text/plain": [
- "True"
- ]
- },
- "execution_count": 31,
- "metadata": {},
- "output_type": "execute_result"
- }
- ],
- "source": [
- "def test() -> bool:\n",
- " # If there are `2N` cards in a deck, all with different ranks, \n",
- " # then the probability that A sweeps is `(1/2) ** N`.\n",
- " for N in range(200):\n",
- " assert p_sweep((0, 2 * N)) == (1/2) ** N\n",
- " # Remove a single card from abstract deck\n",
- " assert remove(4, (0, 0, 0, 0, 13)) == (0, 0, 0, 1, 12)\n",
- " assert remove(4, (0, 0, 0, 1, 12)) == (0, 0, 0, 2, 11)\n",
- " assert remove(3, (0, 0, 0, 1, 12)) == (0, 0, 1, 0, 12)\n",
- " assert remove(1, (0, 2)) == (0, 1)\n",
- " assert remove(1, (0, 1)) == (0, 0)\n",
- " # Count cards in abstract deck\n",
- " assert [ncards(d) for d in (deck, tie1, dif1, end4)] == [52, 50, 50, 4]\n",
- " # Select cards\n",
- " assert select1((0, 0, 0, 0, 13)) == [((0, 0, 0, 1, 12), 1.0, 4)]\n",
- " assert select2((0, 0, 0, 0, 13)) == [((0, 0, 0, 2, 11), 0.47058823529411764)]\n",
- " assert select1((0, 2)) == [((0, 1), 1.0, 1)]\n",
- " assert select2((0, 2)) == [((0, 0), 0.5)]\n",
- " return True\n",
- " \n",
- "test()"
+ "for r, s in [(2, 1), \n",
+ " (4, 1), (2, 2), \n",
+ " (6, 1), (3, 2), (2, 3),\n",
+ " (8, 1), (4, 2), (2, 4), \n",
+ " (5, 2), (2, 5)]:\n",
+ " p1 = p_sweeps(permutations(make_deck(r, s)))\n",
+ " p2 = p_sweep(make_abdeck(r, s))\n",
+ " print(f'deck({r}, {s}): p = {p1:.4f}; ∆ = {abs(p1 - p2):.9f}')"
]
}
],