probability_exercise.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Exercise \n",
    "\n",
    "We have a deck of 40 cards. All cards are dealt among the players, players can have 14 or 13 cards. If you receive 13 cards, what is the probability that 3 of them are the same number? (3 fives for example)\n",
    "\n",
    "Note 1: There are 10 unique numbers (ranks), and each number exists in 4 copies.\n",
    "\n",
    "Note 2: We are looking for hands in which _at least_ one number that occurs _at least_ 3 times."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "from math import comb as C  # named C to be consistent with notation"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The total number of (unordered) hands we could get is $C(40, 13)$.\n",
    "\n",
    "Of these, all valid hands contain at least one triplet or quadruplet of at least one number (rank). This means that given a hand consisting of $w$ many ranks that occur once, $x$ many ranks that occur twice, $y$ many ranks that occur three times, and $z$ many ranks that occur four times a valid hand satisfies $w+x+y+z \\leq 10$ (total ranks), $w+2x+3y+4z = 13$ (total cards), and $y>0 \\vee z>0$ (1+ triplet or 1+ quadruplet). This defines a system of inequalities over the domain of integers which we can brute force, because the problem itself is small.\n",
    "\n",
    "For each solution $(w, x, y, z)$ we can then compute the number of ways we can build this hand as the number of (unordered) ways to choose ranks $\\left(C(10, w) \\cdot C(10-w, x) \\cdot C(10-w-x, y) \\cdot C(10-w-x-y, z)\\right)$ times the number of (unordered) ways to choose within a rank $\\left(C(4, 1)^w \\cdot C(4, 2)^x \\cdot C(4, 3)^y \\cdot C(4, 4)^z \\right)$. We then sum all these combinations for each solution to obtain the total number of valid hands, which we can use to compute the desired probability."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "There is an approximately 70.85% chance to have at least one number at least three times.\n"
     ]
    }
   ],
   "source": [
    "cards_to_draw = 13\n",
    "total_ranks = 10\n",
    "\n",
    "integer_solutions = list()\n",
    "for w in range(cards_to_draw + 1):\n",
    "    for x in range(cards_to_draw + 1):\n",
    "        for y in range(cards_to_draw + 1):\n",
    "            for z in range(cards_to_draw + 1):\n",
    "                num_ranks = w + x + y + z\n",
    "                num_cards = w + 2 * x + 3 * y + 4 * z\n",
    "\n",
    "                if (\n",
    "                    num_ranks <= total_ranks\n",
    "                    and num_cards == cards_to_draw\n",
    "                    and (y > 0 or z > 0)\n",
    "                ):\n",
    "                    integer_solutions.append((w, x, y, z))\n",
    "\n",
    "total_combinations = C(40, 13)\n",
    "valid_combinations = sum(\n",
    "    (C(10, w) * C(10 - w, x) * C(10 - w - x, y) * C(10 - w - x - y, z))\n",
    "    * (C(4, 1) ** w * C(4, 2) ** x * C(4, 3) ** y * C(4, 4) ** z)\n",
    "    for w, x, y, z in integer_solutions\n",
    ")\n",
    "\n",
    "odds = valid_combinations / total_combinations\n",
    "print(\n",
    "    f\"There is an approximately {odds:.2%} chance to have at least one number at least three times.\"\n",
    ")\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Alternatively, we can compute the odds of the opposite event, i.e., the odds of getting a hand with exactly 0 triplets and 0 quadruplets. This means that the entire hand is made up of either singles or doubles. The approach to compute these is analogous to the one before, except that we only need to consider $w$ and $x$."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "There is an approximately 70.85% chance to have at least one number at least three times.\n"
     ]
    }
   ],
   "source": [
    "cards_to_draw = 13\n",
    "ranks_remaining = 10\n",
    "integer_solutions = list()\n",
    "for w in range(cards_to_draw + 1):\n",
    "    for x in range(cards_to_draw + 1):\n",
    "        num_ranks = w + x\n",
    "        num_cards = w + 2 * x\n",
    "        if num_ranks <= ranks_remaining and num_cards == cards_to_draw:\n",
    "            integer_solutions.append((w, x))\n",
    "\n",
    "total_combinations = C(40, 13)\n",
    "valid_combinations = sum(\n",
    "    C(ranks_remaining, w) * (C(4, 1) ** w) * C(ranks_remaining - w, x) * (C(4, 2) ** x)\n",
    "    for w, x in integer_solutions\n",
    ")\n",
    "\n",
    "# We counted the valid combinations of the opposite event, so the probability we are looking for\n",
    "# is one minus the computation above\n",
    "odds = 1 - (valid_combinations / total_combinations)\n",
    "print(\n",
    "    f\"There is an approximately {odds:.2%} chance to have at least one number at least three times.\"\n",
    ")\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "As a third solution, we can calculate the odds empirically by sampling hands and checking if it contains at least three of a kind. \n",
    "\n",
    "This will, obviously, require a large sample count, but is an easy way of checking if our analytical solution is correct. If it is, it should be in the same ballpark as our empirical estimate if our sample count (`num_trials`) is large enough."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Empirically, there is an approximately 70.81% chance to have at least one number at least three times.\n"
     ]
    }
   ],
   "source": [
    "import numpy as np\n",
    "\n",
    "rng = np.random.default_rng()\n",
    "deck = (np.arange(1, 11)[:, None] * np.ones((1, 4))).ravel().astype(int)\n",
    "\n",
    "successes = 0\n",
    "num_trials = int(1e6)\n",
    "for trial in range(num_trials):\n",
    "    hand = rng.choice(deck, 13, replace=False)\n",
    "    unique, counts = np.unique(hand, return_counts=True)\n",
    "    if np.any(counts >= 3):\n",
    "        successes += 1\n",
    "\n",
    "odds = successes / num_trials\n",
    "print(f\"Empirically, there is an approximately {odds:.2%} chance to have at least one number at least three times.\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3.8.8 64-bit",
   "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.8.8"
  },
  "orig_nbformat": 4,
  "vscode": {
   "interpreter": {
    "hash": "cd0338b0df120b07017d717bd7db869fbf170f43258ddd98004d1d6a646b6de8"
   }
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}