{ "cells": [ { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "# Week 2-3: Discrete random variables\n", "\n", "\n", " #### [Back to main page](https://petrosyan.page/fall2020math3215)\n" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "___\n", "### Example 1\n", "\n", "Consider the double coin flip experiment $S=\\{HH, HT, TH, TT\\}$, and define a random variable which counts the number of heads in each outcome. Then \n", "$$X(HH)=2,\\;X(HT)=1,\\;X(TH)=1,\\; X(TT)=0.$$\n", " \n", "Observe that Range$(X)=\\{0,1,3\\}$ and the pmf of the random variable is \n", " \n", "$$f(0)=0.25,\\; f(1)=0.5,\\; f(2)=0.15.$$\n", "\n", "Below we draw the line graph and the histogram of the pmf of this random variable. Red circles represent the range of $X$." ] }, { "cell_type": "code", "execution_count": 1, "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAskAAAE/CAYAAAC0Fl50AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjMsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+AADFEAAAUgElEQVR4nO3de5DdZ33f8c/XlgEFyyi2pcb4JjWVnZDAENAYQ4fEvThoCAxMbmNMaN0xTZMQT8uEdkKbZmxSN5M0Q4fY0MQJkEBcA0MZBhrHoZfYcVKbWAanjgFf6svYsYNkOWvZjWRk9PSPc1TWj9bSSrs6v91zXq+ZHe05+ztnn33m+Ku3fvvbdbXWAgAAfMtxQy8AAABWGpEMAAAdkQwAAB2RDAAAHZEMAAAdkQwAAB2RzIpSVa+vqruHXsdSVdUFVfXI0OsAmJSququqLhh6HbBcRDKDqKoHq+of9ve31m5urZ07xJoAeH4Lze2quqSq/iRJWmvf01q78TDPsamqWlWtOYZLhWUhkmEBBjjA6mN2s5xEMitKf5nC+MzFe6rqf1fVk1X1yap60byPv6mq7qiquar6X1X1ikM89w9W1d3j5/lQVd1UVe8cf+ySqvrTqvqPVfVEksur6jur6n9W1a6qeryqrq2q9d3a3ltVX6mqv66qj85f2/iYn6uqHVX1WFX9k2XdLIAVZP6Z5qo6r6q2V9Xuqvp6Vb1/fNgfj/+cq6qnq+q1VXVcVf1CVT00npcfq6qXzHvefzT+2K6q+rfd57m8qj5dVb9XVbuTXDL+3LeM/154rKqurqoXzHu+VlU/U1X3VtVTVfVL43l/y3i9n5p/PLNLJLMa/HiSbUk2J3lFkkuSpKpeleQjSf5ZklOS/GaSz1XVC/snqKpTk3w6yXvHx96d5HXdYa9Jcn+SjUmuTFJJfjnJS5N8d5Izk1zePebtSd6Q5DuTnJPkF+Z97DuSvCTJ6UkuTfLBqvr2I/vSAValDyT5QGvtpIzm46fG93//+M/1rbUTW2u3ZDTTL0ny95L87SQnJrk6SarqZUk+lNGsPS3fmqnzvSWj+b4+ybVJvpnk3UlOTfLaJP8gyc90j9mW5NVJzk/yr5JcM/4cZyb53iRvW8LXzpQQyawGv95ae7S19kSSzyd55fj+f5rkN1trX2ytfbO19rtJnslo6PXemOSu1tpnWmvPJvn1JH/VHfNoa+2q1tqzrbU9rbX7Wmv/rbX2TGttZ5L3J/mB7jFXt9YeHq/tyjx3sO5L8r7W2r7W2vVJnk7iemtgNfvs+AztXFXNZRSwC9mX5O9U1amttadba7ce4jnfnuT9rbX7W2tPZ3Qy46LxpRM/muTzrbU/aa19I8kvJmnd429prX22tbZ/PLtvb63dOp7lD2Z0AqWf3b/SWtvdWrsryV8k+cL48z+Z5A+SfN/it4RpJZJZDebH7N9kdJYhSc5O8nPdwD4zozO/vZcmefjAjdZaS9L/9omH59+oqo1V9Ymq+svxt/F+L6MzE8/3mIe6z71rHOQLrR1gNXpra239gbccfIb2gEsz+u7a16rqtqp60yGe86UZzc8DHkqyJsnfysGz+2+S7Ooe38/uc6rqv1bVX41n97/PwbP76/Pe37PAbbMakcyq9nCSK+cP7Nbat7XWrlvg2MeSnHHgRlXV/Ntj/dmJXx7f94rxtwx/IqNLMOY7c977ZyV59Ci+DoCp0lq7t7X2towuX/uVJJ+uqhfn4DmbjObm2fNun5Xk2YzCtZ/dazO6ZO45n667/Z+SfC3JlvHs/tc5eHbDYYlkhnRCVb1o3tuR/lTybyX5qap6TY28uKp+qKrWLXDs7yd5eVW9dfx53pXRNcOHsi6jSyTmqur0JP9ygWPeVVVnVNXJGQ3iTx7h1wAwdarqJ6pqQ2ttf5K58d3fTLIzyf6Mrj0+4Lok766qzVV1YkZnfj85/k7cp5O8uapeN/5huity+OBdl2R3kqer6ruS/PSyfWHMFJHMkK7P6NtaB94uP5IHt9a2Z3Rd8tVJ/jrJfRn/UN8Cxz6e5MeS/GpG36p7WZLtGV3D/HyuSPKqJE9mFNmfWeCY/5zkCxn9wN/9Sf7dkXwNAFNqW5K7qurpjH6I76LW2t7x5RJXJvnT8WVy52f0A9gfz+g3XzyQZG+Sy5JkfM3wZUk+kdFZ5aeS7MihZ/d7klw8Pva34uQFR6lGl2bCbKmq4zK6JvntrbU/OsrneDDJO1tr/3051wbAwsZnmucyupTigaHXw3RzJpmZUVVvqKr1418Rd+AatUP9xDUAA6uqN1fVt42vaf61JHcmeXDYVTELRDKz5LVJ/k+Sx5O8OaOf0t4z7JIAOIy3ZPTDfY8m2ZLRpRu+Dc4x53ILAADoOJMMAAAdkQwAAJ0j/b20E3Hqqae2TZs2Db0MgCN2++23P95a2zD0OibJzAZWq0PN7BUZyZs2bcr27duHXgbAEauqhw5/1HQxs4HV6lAz2+UWAADQEckAANARyQAA0BHJAADQEckAANARyQAA0BHJAADQEckAANARyQAA0BHJAADQEckAANARyQAA0BHJAADQEckAANARyQAA0BHJAADQEckAANARyQAA0BHJAADQEckAANARyQAA0BHJAADQEckAANARyQAA0BHJAADQEckAANARyQAA0BHJAADQEckAANARyQAA0BHJAADQEckAANARyQAA0BHJAADQEckAANARyQAA0FlUJFfVtqq6u6ruq6qfX+DjF1TVk1V1x/jtFxf7WACWl5kNsHRrDndAVR2f5INJLkzySJLbqupzrbWvdIfe3Fp701E+FoBlYGYDLI/FnEk+L8l9rbX7W2vfSPKJJG9Z5PMv5bEAHDkzG2AZHPZMcpLTkzw87/YjSV6zwHGvrao/T/Jokve01u46gscCsDzMbBZlx45kbm7oVcyW9euTjRuHXgWLtZhIrgXua93tLyU5u7X2dFW9Mclnk2xZ5GNHn6TqJ5P8ZJKcddZZi1gWLN2uXclVVyWXXZaccsrQq4FlYWZzWDt2JFu2JLt3D72S2XLSScm99wrl1WIxkfxIkjPn3T4jozMP/19rbfe896+vqg9V1amLeey8x12T5Jok2bp164JDGZbbrl3JFVckF18skpkaZjaHNTc3CuQbbkg2bx56NbPhgQeSbdtGey+SV4fFRPJtSbZU1eYkf5nkoiQXzz+gqr4jyddba62qzsvoWuddSeYO91gAlpWZzaJt3pycc87Qq4CV6bCR3Fp7tqp+NskfJjk+yUdaa3dV1U+NP/4bSX40yU9X1bNJ9iS5qLXWkiz42GP0tQDMPDMbYHks5kxyWmvXJ7m+u+835r1/dZKrF/tYAI4dMxtg6fwf9wAAoCOSAQCgI5IBAKAjkgEAoCOSAQCgI5IBAKAjkgEAoCOSAQCgI5IBAKAjkgEAoCOSAQCgI5IBAKAjkgEAoCOSAQCgI5IBAKAjkgEAoCOSAQCgI5IBAKAjkgEAoCOSAQCgI5IBAKAjkgEAoCOSAQCgI5IBAKAjkgEAoCOSAQCgI5IBAKAjkgEAoCOSAQCgI5IBAKAjkgEAoCOSAQCgI5IBAKAjkgEAoCOSAQCgI5IBAKAjkgEAoCOSAQCgI5IBAKAjkgEAoCOSAQCgI5IBAKAjkgEAoCOSAQCgI5IBAKAjkgEAoCOSAQCgI5IBAKAjkgEAoCOSAQCgI5IBAKAjkgEAoCOSAQCgI5IBAKAjkgEAoCOSAQCgI5IBAKAjkgEAoCOSAQCgI5IBAKAjkgEAoCOSAQCgI5IBAKAjkgEAoCOSAQCgI5IBAKAjkgEAoCOSAQCgI5IBAKAjkgEAoCOSAQCgI5IBAKAjkgEAoCOSAQCgI5IBAKAjkgEAoCOSAQCgI5IBAKAjkgEAoCOSAQCgI5IBAKAjkgEAoCOSAQCgI5IBAKAjkgEAoCOSAQCgI5IBAKAjkgEAoCOSAQCgI5IBAKAjkgEAoCOSAQCgI5IBAKAjkgEAoCOSAQCgI5IBAKCzqEiuqn9TVd8Yv/3BAh//YFXtGb/trqofm/exZ6tq7/hj/3c5Fw/AwcxsgKU7bCRX1QlJLk9yYZJvT3JBVb25O+zPk3x3a21tkv+Q5KPdx1/RWlvbWnvx0pcMwPMxswGWx5pFHHNJkidbazclSVXdlORdST5/4IDW2jXzjv94RgMagMm7JKt0Zu/YkczNDb2K2fDAA0OvYHbZ+8lZvz7ZuPHoH7+YSD43yePzbj+Y5HWHOP6qJPfOu92S3FFVLclnWmvvWOhBVfXxJD+cJCeffPIiljWddu1Krroqueyy5JRThl4NLD+v8WNuVc7sHTuSLVuS3buX/FQcgV27hl7B7Diw19u2DbuOWXLSScm99x59KC8mkmuB+9qCB1a9O8kPJvmeeXef11r7clW9LMntVfXF1trVBz3haBC/I0m2bt264PPPgl27kiuuSC6+WEAwnbzGj7lVObPn5kaBfMMNyebNS302Dufmm5N3vtM/SibpwF7/9m8nr3/9sGuZBQ88MPoHydzcsY3kryX5x/Nub0ryWH9QVf1Ikl9N8kOttfsO3N9a+/L4z69U1ReTvCHJQQMXgGWxqmf25s3JOedM6rPNLt/yH84ZZ3iNrxaL+e0WH0vykqp6fVW9OMkPJPnQ/AOq6vwk1yV5V2vtC/Pu31BVpx14P8mrkty6XIsH4CBmNsAyOGwkt9aeSfJLSf5HkrkkN7fWPldV11bVtePDfiejs9If6H5t0MuS3F9Ve5I8nOTW1tqVy/1FADBiZgMsj8VcbpHW2vuSvK+77+3z3v+u53ncTUnWLmWBABwZMxtg6fwf9wAAoCOSAQCgI5IBAKAjkgEAoCOSAQCgI5IBAKAjkgEAoCOSAQCgI5IBAKAjkgEAoCOSAQCgI5IBAKAjkgEAoCOSAQCgI5IBAKAjkgEAoCOSAQCgI5IBAKAjkgEAoCOSAQCgI5IBAKAjkgEAoCOSAQCgI5IBAKAjkgEAoCOSAQCgI5IBAKAjkgEAoCOSAQCgI5IBAKAjkgEAoCOSAQCgI5IBAKAjkgEAoCOSAQCgI5IBAKAjkgEAoCOSAQCgI5IBAKAjkgEAoCOSAQCgI5IBAKAjkgEAoCOSAQCgI5IBAKAjkgEAoCOSAQCgI5IBAKAjkgEAoCOSAQCgI5IBAKAjkgEAoCOSAQCgI5IBAKAjkgEAoCOSAQCgI5IBAKAjkgEAoCOSAQCgI5IBAKAjkgEAoCOSAQCgI5IBAKAjkgEAoCOSAQCgI5IBAKAjkgEAoCOSAQCgI5IBAKAjkgEAoCOSAQCgI5IBAKAjkgEAoCOSAQCgI5IBAKAjkgEAoCOSAQCgI5IBAKAjkgEAoCOSAQCgI5IBAKAjkgEAoCOSAQCgI5IBAKAjkgEAoCOSAQCgI5IBAKAjkgEAoCOSAQCgI5IBAKAjkgEAoCOSAQCgI5IBAKAjkgEAoCOSAQCgI5IBAKAjkleaJ5547p8cW/Z78uw502TnzuS220Z/Mhn2fLJmeL9F8kpy3XXJBReM3r/ggtFtjh37PXn2nGly3XXJ2WcnF144+tPr+diz55M14/tdrbWh13CQdevWtVe/+tVDL2Oy9u1Lbr01X9tf+XpOzLo8leOyPznppKRq6NVNn9aS3buzP8flqayz35PwnD0/Lefliaw97pnk/POTE04YenXL5qabbrq9tbZ16HVM0nLM7D17kj/7s+TlL0/Wrl2mhR1L+55N7rgjafu/dV8dl7zylckJa4Zb1yLNzSX33JOcc06yfv3Qq1kkez5Zq3y/9+xJ7rwzOe+8Q8+UQ83slf9Vzoq9e5OqHJf9eWH2pjJ+Ue7fnxx//LBrm0b7R/tb9nty5u35uuzOmuwb/YNk796pimSOzpo1o//07rxz6JUs1pok3d+rLcmXh1jL0bvnnqFXcCTs+WSt/v0+/vjRbDlaKzKSzz333Nx4441DL2Oydu4cfStjz54kz4zuW7s2uf/+ZMOGQZc2lez35D1nz58a3feCtckNN0zVntcMfidiuWb2jh2js22rwhNPjC4Zembvt+574YuSG29MTj55qFUdkX37Vtm/T+35ZE3Bfq9fn2zceOhjDjWzV2Qkz6QNG5IPfzi59NLRf0H79o1uT1E8rCj2e/LsOYexcePh/0JbOU5OPvreg1/P56+OeFid7Plk2e8VeU3y1q1b2/bt24dexjB27kwefDDZtEk8TIL9nrwp3/Oqmrlrks3sB6f29bwi2fPJmvL9PtTMdiZ5pdmwYSpfhCuW/Z48e8408XqePHs+WTO8334FHAAAdEQyAAB0RDIAAHREMgAAdEQyAAB0RDIAAHREMgAAdEQyAAB0RDIAAHREMgAAdEQyAAB0RDIAAHREMgAAdEQyAAB0RDIAAHREMgAAdKq1NvQaDlJVO5M8NPQ6BnRqkseHXsQMsd+TN817fnZrbcPQi5gkM3uqX88rlT2frGne7+ed2SsykmddVW1vrW0deh2zwn5Pnj1nmng9T549n6xZ3W+XWwAAQEckAwBARySvTNcMvYAZY78nz54zTbyeJ8+eT9ZM7rdrkgEAoONMMgAAdETyClNV26rq7qq6r6p+fuj1TLOq+khV7aiqvxh6LbOgqs6sqj+qqq9W1V1V9c+HXhMslZk9Web25JjZLrdYUarq+CT3JLkwySNJbkvyttbaVwZd2JSqqu9P8nSSj7XWvnfo9Uy7qjotyWmttS9V1boktyd5q9c3q5WZPXnm9uSY2c4krzTnJbmvtXZ/a+0bST6R5C0Dr2lqtdb+OMkTQ69jVrTWHmutfWn8/lNJvprk9GFXBUtiZk+YuT05ZrZIXmlOT/LwvNuPZMZekMyGqtqU5PuSfHHYlcCSmNnMhFmd2SJ5ZakF7nM9DFOlqk5M8l+S/IvW2u6h1wNLYGYz9WZ5ZovkleWRJGfOu31GkkcHWgssu6o6IaNhe21r7TNDrweWyMxmqs36zBbJK8ttSbZU1eaqekGSi5J8buA1wbKoqkry4SRfba29f+j1wDIws5laZrZIXlFaa88m+dkkf5jRBfKfaq3dNeyqpldVXZfkliTnVtUjVXXp0Guacn83yTuS/P2qumP89sahFwVHy8yePHN7omZ+ZvsVcAAA0HEmGQAAOiIZAAA6IhkAADoiGQAAOiIZAAA6IhkAADoiGQAAOiIZAAA6/w8dZIEaE4YnvAAAAABJRU5ErkJggg==\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "# nbi:hide_in\n", "# library\n", "import matplotlib.pyplot as plt\n", "import numpy as np\n", "plt.rcParams[\"figure.figsize\"] = (12, 5)\n", "\n", "range_x = np.array([0, 1, 2])\n", "pmf_values = np.array([1/4, 1/2, 1/4])\n", "\n", "fig, [ax1, ax2] = plt.subplots(1,2, num=1, clear=True)\n", "\n", "ax1.set_ylim(-0.05,0.6) \n", "ax1.set_xlim(-0.7, 2.6)\n", "ax1.axhline(y=0, color='k')\n", "ax1.set_xticks(range_x)\n", "ax1.set_yticks(pmf_values)\n", "\n", "ax2.set_ylim(-0.05, 0.6) \n", "ax2.set_xlim(-0.7, 2.6)\n", "ax2.axhline(y=0, color='k')\n", "ax2.set_xticks(range_x)\n", "ax2.set_yticks(pmf_values)\n", "\n", "\n", "# Plotting line graphs using plt.stem with stems removed\n", "ax1.scatter(range_x,np.zeros(range_x.shape), color =\"red\", s=20)\n", "markers,stems,base = ax1.stem(range_x, pmf_values, markerfmt=' ', linefmt=\"blue\", basefmt=\"black\", use_line_collection=True)\n", "stems.set_linewidth(1.3)\n", "ax1.set_title(\"Line graph\")\n", "\n", "\n", "# PLotting with plt.bar instead of plt.hist works better when f(x) are knowwn\n", "ax2.scatter(range_x,np.zeros(range_x.shape), color =\"red\", s=20)\n", "ax2.bar(range_x, pmf_values, width=1, color=(1, 1, 1, 0), edgecolor='blue', linewidth=1.3, label=\"Histogran\")\n", "ax2.set_title(\"Histogram\")\n", "\n", "plt.show();\n", " " ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "___\n", "### Example 2\n", "\n", "If Range$(X)=\\{x_1,\\dots,x_k\\}$ and $P(x_1)=\\cdots=P(x_k)=\\frac{1}{k}$ then we say that $X$ has uniform distribution.\n", "\n", "in the example below, the range of $X$ is \n", "\n", "$$ \\text{Range}(X)=\\{-3, -6, 5, 9, 2, 1, 3, 8, 7, 10\\}$$\n", "\n", "with $f(x)=0.1$ for every $x\\in \\text{Range}(X).$" ] }, { "cell_type": "code", "execution_count": 2, "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAsIAAAE/CAYAAABM9qWDAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjMsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+AADFEAAAfJ0lEQVR4nO3dfbRdd13n8feHPECAUh5anpLctmp8yHJapNdQR0QUqilVKz6MrQjCULPqsj7MLJHAOI4ZhpnUh1nqap1MxlbUUTKdRdGIGRvAseAokhZbaIplMiW0oZW24KSAwSbtd/64p8zp7U3uvvece3dOfu/XWmfd/fDbZ3/P3jnffLLPvjmpKiRJkqTWPKnvAiRJkqQ+GIQlSZLUJIOwJEmSmmQQliRJUpMMwpIkSWqSQViSJElNMgirSZnx20n+PsmH+65HkgRJzk5SSVYO5p+X5ANJPp/kV/uuT6eelX0XIPXkpcCFwLqq+mLfxUiS5rQFeBB4RvnFB1oCXhFWq84CDhqCJemkdhZwhyFYS8UgrImS5GCStyS5Y3Bbw28neUqSlyc5lOTnktyf5L4k35vkVUk+keRzSd46eI43Ar8FfFOSLyTZ1u+rkqRTX5L1SW5I8kCSzya5OsmKJL+S5MEkdwEXD41/B/CjwM8NevUr+6pdpy5vjdAkeg3wncAXgT8Gfh54H/B84CnAWuD1wH8B3gucD0wBtyTZVVXXJnkEuLyqXrr85UtSW5KsAN4D/BnwWuARYBr4MeC7gG9gpqe/67Ftqur1SQAOVdXPL3fNaoNXhDWJrq6qe6rqc8DbgcsGy48Cb6+qo8Au4Azg16vq81W1H9gPnNtLxZLUtk3AC4E3VdUXq+pLVfUXwD8Dfm2op/+HXqtUcwzCmkT3DE1/ipnmCvDZqnpkMH1k8PMzQ2OPAE9f4tokSU+0HvhUVR2btfyFPLGnS8vGIKxJtH5oegq4t69CJEmd3ANMPfbfog25jyf2dGnZGIQ1iX4iybokzwbeCvy3vguSJJ3Qh5kJvduTPG3wS87fDFwP/NSgpz8L2NprlWqOQViT6A+AvcBdg8e/67ccSdKJDG5b+27gq4C7gUPADzHzS803ArcBHwFu6KtGtSn+13yaJEkOMvO/Pbyv71okSdJk84qwJEmSmmQQliRJUpO8NUKSJElN8oqwJEmSmmQQliRJUpNm/8fWy+aMM86os88+u6/dS9JIbrnllger6sy+61gu9mxJk+x4Pbu3IHz22Wdz880397V7SRpJkqa+CtaeLWmSHa9ne2uEJEmSmmQQliRJUpMMwpIkSWqSQViSJElNMghLkiSpSQZhSZIkNckgLEmSpCYZhCVJktQkg7AkSZKaZBCWJElSkwzCkiRJapJBWJIkSU0yCEuSJKlJBmFJkiQ1ySAsSZKkJhmEJUmS1CSDsCRJkpo0bxBO8okkjyb50nHWJ8mtSR5OciTJD4+/TElSV/ZtSeqmyxXh3wB+5ATr/zWwFngycCWwcwx1SZIWz74tSR3MG4Sr6mrg7hMMuRTYVTOuBVYnOW9cBUqSFsa+LUndjOMe4ecAdwzNfwE4dwzPK0laGvZtSQJWjuE5MseyR+ccmPwe8H0Az372s8ewa0nSInTq2/ZsSeN03ra9HD5ytO8yHmccQfhBYOPQ/NOB2+caWFWvBV4LMD09XWPYtyRp4Tr1bXu2pHE6fOQoB7df3Mu+c9Xcy8dxa8T1wKWD30J+I/BwVd02hueVJC0N+7Yk0eGKcJJPAeuAJyU5BvwesBqgql4DbAO+H3gYeAT4sSWrVpI0L/u2JHUzbxCuqrPmWV/APxlbRZKkkdi3Jakbv1lOkiRJTTIIS5IkqUkGYUmSJDXJICxJkqQmGYQlSZLUJIOwJEmSmmQQliRJUpMMwpIkSWqSQViSJElNMghLkiSpSQZhSZIkNckgLEmSpCYZhCVJktQkg7AkSZKaZBCWJElSkzoF4SSbk9yZ5ECSrXOsf1aSdyf5aJIPJ/n68ZcqSerCni1J3cwbhJOsAK4BLgI2Apcl2Thr2FuBW6vqXOB1wK+Pu1BJ0vzs2ZLUXZcrwpuAA1V1V1U9DOwCLpk1ZiPwfoCq+lvg7CTPG2ulkqQu7NmS1NHKDmPWAvcMzR8CXjJrzG3A9wF/kWQTcBawDvjM8KAkW4AtAFNTU4ssWZJ0AvZsSUvmvG17OXzk6KK2PX3NqjFXM7ouQThzLKtZ89uBX09yK/Ax4G+AY0/YqGonsBNgenp69nNIkkZnz5a0ZA4fOcrB7Rf3XcbYdAnCh4D1Q/PrgHuHB1TVQ8AbAJIE+OTgIUlaXvZsSeqoyz3C+4ANSc5Jshq4FNg9PCDJMwfrAC4HPjBotJKk5WXPlqSO5r0iXFXHklwJ3AisAK6rqv1Jrhis3wF8HfC7SR4B7gDeuIQ1S5KOw54tSd11uTWCqtoD7Jm1bMfQ9F8BG8ZbmiRpMezZktSN3ywnSZKkJhmEJUmS1CSDsCRJkppkEJYkSVKTDMKSJElqkkFYkiRJTTIIS5IkqUkGYUmSJDXJICxJkqQmGYQlSZLUJIOwJEmSmmQQliRJUpMMwpIkSWqSQViSJElN6hSEk2xOcmeSA0m2zrH+9CR/nOS2JPuTvGH8pUqSurBnS1I38wbhJCuAa4CLgI3AZUk2zhr2E8AdVXUe8HLgV5OsHnOtkqR52LMlqbsuV4Q3AQeq6q6qehjYBVwya0wBpyUJ8HTgc8CxsVYqSerCni1JHXUJwmuBe4bmDw2WDbsa+DrgXuBjwE9X1aNjqVCStBD2bEnqaGWHMZljWc2a/07gVuDbga8E3pvkg1X10OOeKNkCbAGYmppaeLWSpPnYsyWd0Hnb9nL4yNFFbXv6mlVjrqZfXYLwIWD90Pw6Zq4iDHsDsL2qCjiQ5JPA1wIfHh5UVTuBnQDT09OzG7MkaXT2bEkndPjIUQ5uv7jvMk4KXW6N2AdsSHLO4JcpLgV2zxpzN/AKgCTPA74GuGuchUqSOrFnS1JH814RrqpjSa4EbgRWANdV1f4kVwzW7wDeBrwjyceY+VjuzVX14BLWLUmagz1bkrrrcmsEVbUH2DNr2Y6h6XuB7xhvaZKkxbBnS1I3frOcJEmSmmQQliRJUpMMwpIkSWqSQViSJElNMghLkiSpSQZhSZIkNckgLEmSpCYZhCVJktQkg7AkSZKaZBCWJElSkwzCkiRJapJBWJIkSU0yCEuSJKlJBmFJkiQ1qVMQTrI5yZ1JDiTZOsf6NyW5dfC4PckjSZ49/nIlSfOxZ0tSN/MG4SQrgGuAi4CNwGVJNg6PqapfrqoXVdWLgLcAN1XV55aiYEnS8dmzJam7LleENwEHququqnoY2AVccoLxlwHvHEdxkqQFs2dLUkddgvBa4J6h+UODZU+Q5KnAZuBdo5cmSVoEe7YkdbSyw5jMsayOM/a7gf91vI/YkmwBtgBMTU11KlCStCD2bKkB523by+EjRxe17elrVo25msnVJQgfAtYPza8D7j3O2Es5wUdsVbUT2AkwPT19vMYsSVo8e7bUgMNHjnJw+8V9lzHxutwasQ/YkOScJKuZaZy7Zw9KcjrwrcAfjbdESdIC2LMlqaN5rwhX1bEkVwI3AiuA66pqf5IrBut3DIa+GthbVV9csmolSSdkz5ak7rrcGkFV7QH2zFq2Y9b8O4B3jKswSdLi2LMlqRu/WU6SJElNMghLkiSpSQZhSZIkNckgLEmSpCYZhCVJktQkg7AkSZKaZBCWJElSkwzCkiRJapJBWJIkSU0yCEuSJKlJBmFJkiQ1ySAsSZKkJhmEJUmS1CSDsCRJkprUKQgn2ZzkziQHkmw9zpiXJ7k1yf4kN423TElSV/ZsSepm5XwDkqwArgEuBA4B+5Lsrqo7hsY8E/hNYHNV3Z3kuUtVsCTp+OzZktRdlyvCm4ADVXVXVT0M7AIumTXmh4EbqupugKq6f7xlSpI6smdLUkddgvBa4J6h+UODZcO+GnhWkj9PckuS142rQEnSgtizJamjeW+NADLHsprjec4HXgGsAf4qyYeq6hOPe6JkC7AFYGpqauHVSpLmY8+Wlsl52/Zy+MjRXvZ9+ppVvez3VNMlCB8C1g/NrwPunWPMg1X1ReCLST4AnAc8rqlW1U5gJ8D09PTsxixJGp09W1omh48c5eD2i/suQyPocmvEPmBDknOSrAYuBXbPGvNHwLckWZnkqcBLgI+Pt1RJUgf2bEnqaN4rwlV1LMmVwI3ACuC6qtqf5IrB+h1V9fEkfwp8FHgU+K2qun0pC5ckPZE9W5K663JrBFW1B9gza9mOWfO/DPzy+EqTJC2GPVuSuvGb5SRJktQkg7AkSZKaZBCWJElSkwzCkiRJapJBWJIkSU0yCEuSJKlJBmFJkiQ1ySAsSZKkJhmEJUmS1CSDsCRJkppkEJYkSVKTDMKSJElqkkFYkiRJTTIIS5IkqUkGYUmSJDWpUxBOsjnJnUkOJNk6x/qXJzmc5NbB4xfGX6okqQt7tiR1s3K+AUlWANcAFwKHgH1JdlfVHbOGfrCqvmsJapQkdWTPlqTuulwR3gQcqKq7quphYBdwydKWJUlaJHu2JHU07xVhYC1wz9D8IeAlc4z7piS3AfcCP1tV+2cPSLIF2AIwNTW18GolSfOxZ6s5523by+EjR5d9v6evWbXs+9R4dQnCmWNZzZr/CHBWVX0hyauAPwQ2PGGjqp3AToDp6enZzyFJGp09W805fOQoB7df3HcZmkBdbo04BKwfml/HzBWEL6uqh6rqC4PpPcCqJGeMrUpJUlf2bEnqqEsQ3gdsSHJOktXApcDu4QFJnp8kg+lNg+f97LiLlSTNy54tSR3Ne2tEVR1LciVwI7ACuK6q9ie5YrB+B/ADwI8nOQYcAS6tKj9Gk6RlZs+WpO663CP82Edne2Yt2zE0fTVw9XhLkyQthj1bkrrxm+UkSZLUJIOwJEmSmmQQliRJUpMMwpIkSWqSQViSJElNMghLkiSpSQZhSZIkNckgLEmSpCYZhCVJktQkg7AkSZKaZBCWJElSkwzCkiRJapJBWJIkSU0yCEuSJKlJnYJwks1J7kxyIMnWE4z7xiSPJPmB8ZUoSVoIe7YkdTNvEE6yArgGuAjYCFyWZONxxl0F3DjuIiVJ3dizJam7LleENwEHququqnoY2AVcMse4nwTeBdw/xvokSQtjz5akjroE4bXAPUPzhwbLvizJWuDVwI7xlSZJWgR7tiR1tLLDmMyxrGbN/xrw5qp6JJlr+OCJki3AFoCpqamuNUqSurNnayKdt20vh48cXdS2p69ZNeZq1IouQfgQsH5ofh1w76wx08CuQUM9A3hVkmNV9YfDg6pqJ7ATYHp6enZjliSNzp6tiXT4yFEObr+47zLUmC5BeB+wIck5wKeBS4EfHh5QVec8Np3kHcB7ZjdUSdKysGdLUkfzBuGqOpbkSmZ+s3gFcF1V7U9yxWC995hJ0knCni1J3XW5IkxV7QH2zFo2ZzOtqtePXpYkabHs2ZLUjd8sJ0mSpCYZhCVJktQkg7AkSZKaZBCWJElSkwzCkiRJapJBWJIkSU0yCEuSJKlJBmFJkiQ1ySAsSZKkJhmEJUmS1CSDsCRJkppkEJYkSVKTDMKSJElqkkFYkiRJTeoUhJNsTnJnkgNJts6x/pIkH01ya5Kbk7x0/KVKkrqwZ0tSNyvnG5BkBXANcCFwCNiXZHdV3TE07P3A7qqqJOcC1wNfuxQFS5KOz54tSd11uSK8CThQVXdV1cPALuCS4QFV9YWqqsHs04BCktQHe7YkddQlCK8F7hmaPzRY9jhJXp3kb4E/Af75eMqTJC2QPVuSOpr31gggcyx7wtWDqno38O4kLwPeBrzyCU+UbAG2AExNTS2sUklSF/Zs9eK8bXs5fOToorc/fc2qMVYjddMlCB8C1g/NrwPuPd7gqvpAkq9MckZVPThr3U5gJ8D09LQfxUnS+Nmz1YvDR45ycPvFfZchLUiXWyP2ARuSnJNkNXApsHt4QJKvSpLB9IuB1cBnx12sJGle9mxJ6mjeK8JVdSzJlcCNwArguqran+SKwfodwPcDr0tyFDgC/NDQL2JIkpaJPVuSuutyawRVtQfYM2vZjqHpq4CrxluaJGkx7NmS1I3fLCdJkqQmGYQlSZLUJIOwJEmSmmQQliRJUpMMwpIkSWqSQViSJElNMghLkiSpSQZhSZIkNckgLEmSpCYZhCVJktQkg7AkSZKaZBCWJElSkwzCkiRJapJBWJIkSU3qFISTbE5yZ5IDSbbOsf41ST46ePxlkvPGX6okqQt7tiR1M28QTrICuAa4CNgIXJZk46xhnwS+tarOBd4G7Bx3oZKk+dmzJam7LleENwEHququqnoY2AVcMjygqv6yqv5+MPshYN14y5QkdWTPlqSOugThtcA9Q/OHBsuO543A/xilKEnSotmzJamjlR3GZI5lNefA5NuYaaovPc76LcAWgKmpqY4lSpIWwJ59kjhv214OHznadxnL5vQ1q/ouQVqwLkH4ELB+aH4dcO/sQUnOBX4LuKiqPjvXE1XVTgb3ok1PT8/ZmCVJI7FnnyQOHznKwe0X912GpBPocmvEPmBDknOSrAYuBXYPD0gyBdwAvLaqPjH+MiVJHdmzJamjea8IV9WxJFcCNwIrgOuqan+SKwbrdwC/ADwH+M0kAMeqanrpypYkzcWeLUnddbk1gqraA+yZtWzH0PTlwOXjLU2StBj2bEnqxm+WkyRJUpMMwpIkSWqSQViSJElNMghLkiSpSQZhSZIkNckgLEmSpCYZhCVJktQkg7AkSZKaZBCWJElSkwzCkiRJapJBWJIkSU0yCEuSJKlJBmFJkiQ1ySAsSZKkJhmEJUmS1KROQTjJ5iR3JjmQZOsc6782yV8l+cckPzv+MiVJXdmzJamblfMNSLICuAa4EDgE7Euyu6ruGBr2OeCngO9dkiolSZ3YsyWpuy5XhDcBB6rqrqp6GNgFXDI8oKrur6p9wNElqFGS1J09W5I6mveKMLAWuGdo/hDwksXsLMkWYAvA1NTUYp5iLH7xF2cek2QSax7VKK951OPV1/Ge1P32da5afF90cMr17PO27eXwkcnL7KevWdV3CZLmkao68YDkB4HvrKrLB/OvBTZV1U/OMfYXgS9U1a/Mt+Pp6em6+eabF1X0qBKY52WfdCax5lGN8ppHPV59He9J3W9f56rP90WSW6pqup+9H9+p2LPP3vonHNx+cS/7lnRqOF7P7nJrxCFg/dD8OuDecRUmSRore7YkddQlCO8DNiQ5J8lq4FJg99KWJUlaJHu2JHU07z3CVXUsyZXAjcAK4Lqq2p/kisH6HUmeD9wMPAN4NMnPABur6qElrF2SNIs9W5K66/LLclTVHmDPrGU7hqb/jpmP3yRJPbNnS1I3frOcJEmSmmQQliRJUpMMwpIkSWqSQViSJElNMghLkiSpSQZhSZIkNckgLEmSpCYZhCVJktQkg7AkSZKaZBCWJElSkwzCkiRJapJBWJIkSU0yCEuSJKlJBmFJkiQ1qVMQTrI5yZ1JDiTZOsf6JPmNwfqPJnnx+EuVJHVhz5akbuYNwklWANcAFwEbgcuSbJw17CJgw+CxBfhPY65zfB544PE/J8Ek1jyqUV7zqMerr+M9qfvt61y1+L7owJ59knjgAdi3b/LqXqxRX+8o2/d1rPs8x30dr0k8T/PockV4E3Cgqu6qqoeBXcAls8ZcAvxuzfgQ8MwkLxhzraN75zvhrLNmps86a2b+ZDeJNY9qlNc86vHq63hP6n77Olctvi+6s2f37bG6L7xwsuperFFf7yjb93Ws+zzHfR2vSTxPHaSqTjwg+QFgc1VdPph/LfCSqrpyaMx7gO1V9ReD+fcDb66qm4/3vKeddlqdf/75Y3gJHR09Ch/6EAcefZRPczqnc3hm+TOeAcny1bEQVfDQQwAcnpSaRzXKax71ePV1vCd1v32dq8dtu5Zv5dPwpCfBBRfAqlXz1z0mN9100y1VNb1sO+zoVOvZPPooB9/8Hs6+6rt6Oc8LNlT3l01C3Ys16usdZfu+jnWf57iv4zWJ52mW4/XslR22netvpdnpucsYkmxh5mM4nvzkJ3fY9Rh96Utf/gv2yXzp/y9/9FFYsWJ5a+lq6A/NxNQ8qlFe86jHq6/jPan77etcDW17GjOBmGTmPX4qBo2FW5Kevfr5X8XBC940enUL8S0zP5505POPFXTyn+ehv2u+bBLqXqxRX+8o2/d1rPs8x30dr0k8T11V1QkfwDcBNw7NvwV4y6wx/xm4bGj+TuAFJ3re888/v5bV/fdXrVlTNXM9aeaxZs3M8pPVJNY8qlFe86jHq6/jPan77etcnSTvC+Dmmqd/9vGwZ/dsUuterEntI6Po8xxPYt89Sd4Tx+vZXe4R3gdsSHJOktXApcDuWWN2A68b/CbyBcDhqrpvUcl8qZx5Jlx7LaxZM/Px65o1M/Nnntl3Zcc3iTWPapTXPOrx6ut4T+p++zpXLb4vFsae3adJrXuxJrWPjKLPczyJffckf0/Me48wQJJXAb8GrACuq6q3J7kCoKp2JAlwNbAZ+AfgDXWCe80Apqen6+abTzhkaTzwABw8CGeffdKchHlNYs2jGuU1j3q8+jrek7rfvs5Vz++LJCflPcJgzz4pTGrdizWpfWQUfZ7jSey7J2nP7hSEl0JvTVWSxuBkDsJLwZ4taZIdr2f7zXKSJElqkkFYkiRJTTIIS5IkqUkGYUmSJDXJICxJkqQmGYQlSZLUJIOwJEmSmmQQliRJUpMMwpIkSWqSQViSJElNMghLkiSpSamqfnacPAB8qpedwxnAgz3te7EmseZRjfKaRz1efR3vSd1vX+eqz/fFWVV1Zk/7Xnb27EWZ1LoXa1L7yCj6PMeT2HdPup7dWxDuU5Kbq2q67zoWYhJrHtUor3nU49XX8Z7U/fZ1rlp8X7RoUs/zpNa9WJPaR0bR5zmexL57Mr4nvDVCkiRJTTIIS5IkqUmtBuGdfRewCJNY86hGec2jHq++jvek7revc9Xi+6JFk3qeJ7XuxZrUPjKKPs/xJPbdk+490eQ9wpIkSVKrV4QlSZLUuGaDcJKfTHJnkv1JfqnverpI8rYkH01ya5K9SV7Yd01LIcl1Se5Pcvsitl2f5H8m+fjg3P70ArZ9SpIPJ7ltsO22he5/FEkOJvnY4PzevEz7/JrB/h57PJTkZxaw/b8YHKvbk7wzyVMWsO1PD7bb32Wfc/25SPLsJO9N8r8HP5/Vdf+aPJPWt+3Znba1Zy9sn7317MH2nfv2xPTsqmruAXwb8D7gyYP55/ZdU8e6nzE0/VPAjr5rWqLX+TLgxcDti9j2BcCLB9OnAZ8ANnbcNsDTB9OrgL8GLljG130QOKPH474C+Dtm/q/FLuPXAp8E1gzmrwde33HbrwduB54KrBy8Hzcs9M8F8EvA1sH0VuCqvo6fj6V9TGLftmd32taevfj9L1vPHoxfUN+elJ7d6hXhHwe2V9U/AlTV/T3X00lVPTQ0+zTglLzBu6o+AHxukdveV1UfGUx/Hvg4M2/+LttWVX1hMLtq8Dglj/FxvAL4P1W1kC9NWAmsSbKSmeZ4b8ftvg74UFX9Q1UdA24CXn2iDY7z5+IS4HcG078DfG/XwjVxJq5v27M7bWvPXrzl7NmwwL49KT271SD81cC3JPnrJDcl+ca+C+oqyduT3AO8BviFvus5mSU5G/gGZq4SdN1mRZJbgfuB91ZV523HoIC9SW5JsmUZ9/uYS4F3dh1cVZ8GfgW4G7gPOFxVeztufjvwsiTPSfJU4FXA+gXWC/C8qrpvUM99wHMX8RyaDBPZt+3Z3dmzF2w5ezaMp2+fdD37lA3CSd43uI9l9uMSZv5F9CzgAuBNwPVJ0mvBA/PUTVX9q6paD/w+cGW/1Z68kjwdeBfwM7OuypxQVT1SVS8C1gGbknz9UtU4h2+uqhcDFwE/keRly7XjJKuB7wH++wK2eRYz/7o/B3gh8LQkP9Jl26r6OHAV8F7gT4HbgGMLLFunmEns2/bs8bBnL8xy92w4dfv2yr4LWCpV9crjrUvy48ANNXOTyoeTPMrM918/sFz1Hc+J6p7lD4A/Af7NEpYzkZKsYqah/n5V3bCY56iq/5vkz4HNzPwreMlV1b2Dn/cneTewCfjAcuybmUb+kar6zAK2eSXwyap6ACDJDcA/Bf5rl42r6lrg2sG2/x44tKCKZ3wmyQuq6r4kL2DmqpAm1CT2bXv26OzZi7LsPRvG0rdPup59yl4RnscfAt8OkOSrgdXAg71W1EGSDUOz3wP8bV+1nKwGV4iuBT5eVf9xgduemeSZg+k1zDSNZTnGSZ6W5LTHpoHvYJma+cBlLOAjtoG7gQuSPHVw3F/BzP19nSR57uDnFPB9i9g/wG7gRwfTPwr80SKeQ5Nh4vq2PXt+9uxFW/aeDWPp2yddz27yCzUGHylcB7wIeBj42ar6s36rml+SdwFfAzwKfAq4YnDPzyklyTuBlzNzteczwL8Z/Cu0y7YvBT4IfIyZ4wTw1qra02Hbc5m5eX8FM/9IvL6q/u2CX8AiJPkK4N2D2ZXAH1TV25dp308F7gG+oqoOL3DbbcAPMfPx2N8Alz/2y0wdtv0g8BzgKPAvq+r984x/wp8LZsLR9cAUM03+B6tqUb+0o5PbJPZte3anbe3ZC993Lz17sH3nvj0pPbvJICxJkiS1emuEJEmSGmcQliRJUpMMwpIkSWqSQViSJElNMghLkiSpSQZhSZIkNckgLEmSpCYZhCVJktSk/we00Qw1DghtxAAAAABJRU5ErkJggg==\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "# nbi:hide_in\n", "# library\n", "import matplotlib.pyplot as plt\n", "import numpy as np\n", "plt.rcParams[\"figure.figsize\"] = (12, 5)\n", "\n", "K=10\n", "range_x=np.array([-3, -6, 5, 9, 2, 1, 3, 8, 7, 10])\n", "pmf_values=np.ones(range_x.size)/range_x.size\n", "\n", "# Sort x and y according to order in x\n", "sortargs = range_x.argsort()\n", "range_x = range_x[sortargs]\n", "pmf_values = pmf_values[sortargs]\n", "\n", "#cdf values using cumsum function with padding\n", "cdf_values = np.cumsum(pmf_values)\n", "\n", "def padding(cdf_values, range_x):\n", " edge = (range_x[0]-2, range_x[-1]+2)\n", " range_x_padded = np.pad(range_x, (1,1), 'constant', constant_values=edge)\n", " cdf_values_padded = np.pad(cdf_values, (1,1), 'constant', constant_values=(0, 1))\n", " return cdf_values_padded, range_x_padded, edge\n", "\n", "cdf_values_padded, range_x_padded, edge = padding(cdf_values, range_x)\n", "# setting up the figure\n", "fig, [ax1, ax2] = plt.subplots(1,2, num=1, clear=True)\n", "\n", "ax1.set_ylim(-0.07,1.1)\n", "ax1.axhline(y=0, color='k')\n", "ax1.set_xlim(edge)\n", "ax1.set_xticks(range_x)\n", "ax1.set_yticks(cdf_values_padded)\n", "\n", "ax2.set_ylim(-0.07,1.1) \n", "ax2.axhline(y=0, color='k')\n", "ax2.set_xlim(edge)\n", "ax2.set_xticks(range_x)\n", "ax2.set_yticks(cdf_values_padded)\n", "\n", "\n", "# plot line grapghs\n", "ax1.scatter(range_x,np.zeros(range_x.shape), color =\"red\", s=20)\n", "markers,stems,base = ax1.stem(range_x, pmf_values, markerfmt=' ',\n", " linefmt=\"blue\", basefmt=\"black\", use_line_collection=True)\n", "stems.set_linewidth(1)\n", "ax1.set_title(\"pmf\")\n", "\n", "# plot cdf using step function\n", "ax2.scatter(range_x,np.zeros(range_x.shape), color =\"red\", s=20)\n", "ax2.step(range_x_padded, cdf_values_padded, linewidth=1, where='post')\n", "ax2.set_title(\"cdf\")\n", "\n", "plt.show();" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "___\n", "### Relative frequency histogram for pmf\n", "\n", "Suppose our data consists of values of the random variable on a sequence of trial outcomes. As mentioned earlier, probability of an event represents how frequent the experiment outcome terminates in the event, in a large number of repetitive trials. Hence, the pmf at $x\\in\\text{Range}(X)$ can be empirically estimated using the relative frequency \n", "\n", " $$f_{\\text{emp}}(x)= \\frac{\\text{number of elements in data = }x}{ \\text{ size of data}}.$$\n", " \n", "Resulting relative frequency histogram will approximate the pmf histogram.\n", "\n", "In the example below, we consider the following experiment: a dice is tossed twice and the random variable is the maximum of the two tosses: \n", "\n", "$$S=\\{(i,j):\\; 1\\leq i\\leq n, 1\\leq j\\leq 6\\}$$\n", "\n", "and for any $s=(i,j)$,\n", " \n", " $$X(i,j)=\\max\\{i,j\\}.$$\n", " \n", " It can be seen that in this case, $\\text{range}(X)=\\{1,\\dots, n\\}$ and, for any $x\\in \\text{range}(X)$,\n", " \n", " $$f(x)=\\frac{2x-1}{n^2}.$$\n", "\n", "The data consist of values of $X$ on 1000 random pairs of tosses." ] }, { "cell_type": "code", "execution_count": 287, "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAeoAAAEzCAYAAAD+XEDdAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjMsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+AADFEAAAgAElEQVR4nO3dfXgU5bk/8O9NBPEAxYqoCVFCPVAKeTNZo1y0EpEKHBJoLSpY3zhHVryq/lrU1teSnF6//tQjtrWlJYvV9rSW+FY9hOLLaStYWm3ZpUjxBaWW1GSppCgIVcrb/ftjs5udZGcygWSegef7ua69zMzO7Nx7O+z97Ow8zyOqCiIiIgqnfqYDICIiIncs1ERERCHGQk1ERBRiLNREREQhxkJNREQUYizUREREIearUIvINBHZLCJbRORWj+3OFpGDIjK7p/sSERFRV90WahHJA7AEwHQA4wDMFZFxLtvdA+C5nu5LREREufn5Rl0FYIuqvq2q+wA0ApiVY7sbADwJYPth7EtEREQ5+CnUIwC8k7Xc0r4uQ0RGAPg8gKU93ZeIiIjcHedjG8mxrvO4o98G8DVVPSji2NzPvqkNRaIAogAwaNCgyrFjx/oIjYiI6OiXSCT+rqrDcz3np1C3ADg9a7kQQLLTNhEAje1F+mQA/yYiB3zuCwBQ1RiAGABEIhGNx+M+QiMiIjr6iUiz23N+CvU6AKNFZBSAVgBzAFyWvYGqjso62I8ArFTVp0XkuO72JSIiInfdFmpVPSAi1yN1N3cegIdU9VURWdD+fOffpbvdt3dCJyIiOvZJGKe55KVvIiKyiYgkVDWS6zk/l76JiAK1f/9+tLS0YO/evaZDIepVAwcORGFhIfr37+97HxZqIgqdlpYWDBkyBEVFRejUk4ToqKWq2LFjB1paWjBq1Kjud2jHsb6JKHT27t2LYcOGsUjTMUVEMGzYsB5fKWKhJqJQYpGmY9HhnNcs1EREnezYsQPl5eUoLy/HaaedhhEjRmSW9+3b12vHefDBB/HlL38553NTp07F7t27Xfe9//77+Ru+JVioiSj0iooAkd57FBV5H2/YsGHYsGEDNmzYgAULFuArX/lKZnnAgAEAUr83Hjp0qM/e83PPPYchQ4a4Pt+bhfrAgQO98jrUN1ioiSj0mpsB1d57NLuOAeVty5YtKC4uxoIFC1BRUYF33nkHJ554Yub5xsZGXHPNNQCAd999FxdddBEikQiqqqrw8ssv53zNlpYWTJ06FaNHj8Ztt92WWV9YWIidO3di9+7dmD59OsrKylBcXIwnnngC3/rWt7B9+3Z85jOfwZQpUwAAP/3pT1FSUoLi4mLcfvvtmddpaGjAmDFjUF1djWuuuSbzDf7yyy/HTTfdhPPPPx+33347Xn75ZUyYMAFnnXUWJk6ciLfeegtA6lv/7Nmzc8ZIweBd30REPfDaa6/h4YcfxtKlSz2/id5444346le/inPPPRdbt25FTU0NNm3a1GW7V155BevXr8dxxx2HMWPG4IYbbkBBQUHm+VWrVqGoqAjPPPMMAGDXrl0YOnQoFi9ejN/85jc48cQT0dLSgjvvvBPxeBxDhw7FlClTsHLlSpSVleHuu+/G+vXrMWjQIFRXV6Oqqirz2n/+85/xq1/9Cv369cOuXbuwdu1a5OXl4dlnn8Wdd96JRx991FeM1LdYqImIeuDMM8/E2Wef3e12v/zlL7F58+bM8vvvv4+PPvoIJ5xwgmO7KVOmZC5xjx07Fn/9618dRbC0tBS33norbr31VtTW1mLixIldjvX73/8ekydPxsknnwwAuOyyy/Diiy9i7969mDx5Mj7+8Y8DAGbPno2//vWvmf0uvvhi9OuXurC6c+dOXHnllfjzn//c5fW7i5H6Fi99ExH1wKBBgzJ/9+vXD9mjO2b/Zqyq+MMf/pD5bbu1tbVLkQaA448/PvN3Xl5el2/pn/rUpxCPxzF+/Hjccsst+OY3v9nlNdxGmOxu5Mns93LHHXdg6tSp2LRpE55++mnHe+kuRupbLNRERIepX79++PjHP4633noLhw4dwlNPPZV5bsqUKViyZElmecOGDYd1jNbWVgwePBhXXHEFFi5ciPXr1wMAhgwZkrkr/Nxzz8ULL7yAHTt24MCBA2hsbMSkSZNwzjnn4IUXXsDOnTuxf/9+/PznP3c9zq5duzBixAgAwI9+9KPDipX6Bgs1EdERuOeeezBt2jRccMEFKCwszKxfsmQJfvvb36K0tBTjxo3DsmXLDuv1X3nlFZx99tkoLy/Hvffem7lRLBqNYsqUKZgyZQoKCwvxn//5n6iurkZ5eTnOPfdczJgxA2eccQZuueUWVFVV4cILL8T48eMxdOjQnMf52te+hltuuSXnpXUyi5NyEFHovP766/jUpz6VWS4qOvw7tXMZORLYurX3Xi/M9uzZg8GDB2P//v2YNWsWrrvuOtTW1poOy2qdz2+Ak3IQ0VHOlqLaF+666y6sXr0ae/fuxbRp01BTU2M6JOohFmoiomPYt771LdMh0BHib9REREQhxkJNREQUYizUREREIcZCTUREFGIs1EREOeTl5aG8vBzFxcWora3Fzp07u91n8ODBns/v3LkT3//+9zPLyWQSs2fPPuJYAeCWW27JjF52rLj66qvxxBNPdFnfXd465/lox0JNREeFaFMUUi+ZR3J3Ek2bmxzrYokYADjW1S5P9RmuXV7rWN+dE044ARs2bMCmTZtw0kknOUYZO1ydC0hBQUHOQnQ4GhoasH79evzXf/2XY/2xONxnd3nr7UJ98ODBXnutw6KqoXtUVlYqEdnrtddecyxXNFQEHsOgQYMyf//gBz/Q6667LrN87733aiQS0ZKSEv3617/eZZ/du3fr5MmT9ayzztLi4mJ9+umnVVX10ksv1YEDB2pZWZnefPPN+pe//EXHjx+vqqpVVVW6adOmzGtNmjRJ4/G47tmzR+fNm6eRSETLy8szr5WttrZW+/Xrp2VlZdrY2KhXXXWVfuUrX9Hq6mpduHCh62t8+OGHeumll2pJSYlecsklWlVVpevWrevy/h9//HG96qqrVFV1+/btetFFF2kkEtFIJKJr165VVdVFixbpvHnzdNKkSTpq1Cj9zne+k9n/xz/+sZaUlGhpaalefvnl+sEHH2hRUZHu27dPVVV37dqlI0eOzCynXXXVVXrDDTfohAkTdNSoUfr444+rqjrytmnTJj377LO1rKxMS0pK9M033+yS50OHDunNN9+s48eP1+LiYm1sbFRV1YMHD+p1112n48aN0xkzZuj06dMzxxg5cqTW19frxIkTdfny5RqLxTQSiWhpaaledNFF+o9//CMT44IFC7S6ulpHjRqlq1ev1nnz5unYsWMzOeus8/mtqgogri410XhRzvVgoSayW+cPMtShV1+/Id7Q7TbpQnXgwAGdPXu2PvPMM6qq+txzz+n8+fP10KFDevDgQZ0xY4auWbPGsc/+/ft1165dqqra1tamZ555ph46dMhRYFSdBef+++/PFP1kMqmjR49WVdXbbrtNf/KTn6iq6vvvv6+jR4/WPXv2uMarmioeM2bM0AMHDni+xuLFi3XevHmqqvrKK69oXl5et4V67ty5+pvf/EZVVZubm3Xs2LGqmirUEyZM0L1792pbW5uedNJJum/fPt20aZOOGTNG29raVFV1x44dqqp69dVX61NPPZX6/9HQoAsXLuzynq666iqdPXu2Hjx4UF999VU988wzu+Tt+uuv15/+9KeqqvrPf/5TP/zwwy55fuKJJ3TKlCl64MAB/dvf/qann366JpNJffzxx3X69Ol68OBB3bZtm5544omOQn3PPfdkXuPvf/975u877rhDH3jggUyMl156qR46dEiffvppHTJkiG7cuFEPHjyoFRUV+sc//rHL++ppoealbyKyzrUrr+12m48++gjl5eUYNmwY3nvvPXz2s58FADz//PN4/vnncdZZZ6GiogJvvPEG3nrrLce+qorbb78dpaWlmDJlClpbW/Huu+96Hu+SSy7B448/DgB47LHHcPHFF2eOd/fdd6O8vBzV1dXYu3evY6pKNxdffDHy8vI8X+PFF1/E5ZdfDiA1nWZpaWm3r/vLX/4S119/PcrLyzFz5kx88MEHmclBZsyYgeOPPx4nn3wyTjnlFLz77rv49a9/jdmzZ2em4DzppJMAANdccw0efvhhAMDDDz+MefPm5Tze5z73OfTr1w/jxo3LmcMJEybgm9/8Ju655x40NzfnnKFs7dq1mDt3LvLy8nDqqadi0qRJWLduHdauXZuZ6vO0007D+eef79jv0ksvzfy9adMmfOYzn0FJSQkeeeQRvPrqq5nnamtrISIoKSnBqaeeipKSEvTr1w/jx4/H1l4YVs/XyGQiMg3AdwDkAXhQVe/u9PwsAN8AcAjAAQBfVtW17c9tBbAbwEEAB9RlLFMiIjf5g/MDP2b6N+pdu3ahpqYGS5YswY033ghVxW233YZrr3Uv9o888gja2tqQSCTQv39/FBUVOaaNzGXEiBEYNmwYNm7ciEcffRQNDQ0AUkX/ySefxCc/+ckexZ89haXXa4jk/r0+e3127IcOHcJLL73ke8pOVc15jIkTJ2Lr1q1Ys2YNDh48iOLi4pxxZL+m5pib4rLLLsM555yDX/ziF5g6dSoefPBBfOITn3Bsk2s/r/Vp2Tm8+uqr8fTTT6OsrAw/+tGPsHr16i4x9uvXzxFvv379euUegW6/UYtIHoAlAKYDGAdgroiM67TZrwCUqWo5gH8H8GCn589X1XIWaSI6HMmbksaOPXToUDzwwAO47777sH//fkydOhUPPfQQ9uzZAyA1DeX27dsd++zatQunnHIK+vfvjxdeeAHN7TOKZE9NmcucOXNw7733YteuXSgpKQEATJ06Fd/97nczReWPf/xjj9+D22ucd955eOSRRwCkvjFu3Lgxs8+pp56K119/vcv0nRdeeCG+973vZZa7m77zggsuwGOPPYYdO3YAAN57773Mc1deeSXmzp3r+m3aj7fffhuf+MQncOONN2LmzJnYuHFjlzyfd955ePTRR3Hw4EG0tbXhxRdfRFVVFT796U/jySefxKFDh/Duu+86im9nu3fvRn5+Pvbv35/JWVD8XPquArBFVd9W1X0AGgHMyt5AVfdoR9NkEIDwTclFREetutV1vfp6K+as6NH2Z511FsrKytDY2IgLL7wQl112GSZMmICSkhLMnj27S/H94he/iHg8jkgkgkceeQRjx44FAAwbNgwTJ05EcXFxzm5Us2fPRmNjIy655JLMurvuugv79+9HaWkpiouLcdddd/X4/bq9xnXXXYc9e/agtLQU9957L6qqqjL73H333aipqcHkyZORn99xReOBBx5APB7PTN+5dOlSz2OPHz8ed9xxByZNmoSysjIsXLjQkaf3338fc+fO7fF7Snv00UdRXFyM8vJyvPHGG7jyyiu75Pnzn/88SktLUVZWhsmTJ+Pee+/Faaedhi984QsoLCxEcXExrr32Wpxzzjmu04B+4xvfwDnnnIPPfvazmf+fQel2mksRmQ1gmqpe0758BYBzVPX6Ttt9HsD/A3AKgBmq+lL7+r8AeB+p4t2gqjGX40QBRAHgjDPOqGzuzTntiOio0nkawM7dqeLzU9PgRpZ1XKRbNGkR6qrrULC4ANv2bAMAVORXIBFNINoUxbL1HfNBty5sRcGQgr58C0el6upq3HfffYhEgrn4+cQTT+B//ud/8JOf/CSQ4+WSngZ0x44dqKqqwm9/+1ucdtppfXrMvpjmMtcPGF2qu6o+BeApETkPqd+rp7Q/NVFVkyJyCoD/FZE3VPXFHPvHAMSA1HzUPuIiIkvoIpffGHOsz3WZPFYbQ6w253cEMuSGG27AM888g1WrVhmNo6amBjt37sS+fftw11139XmRPhx+CnULgNOzlgsBuP5gpKovisiZInKyqv5dVZPt67eLyFNIXUrvUqiJiMgsr99oe9t3v/vdwI7lJcj3fLj8/Ea9DsBoERklIgMAzAHg+IFHRP5V2m/rE5EKAAMA7BCRQSIypH39IAAXAtjUm2+AiIjoWNbtN2pVPSAi1wN4DqnuWQ+p6qsisqD9+aUAvgDgShHZD+AjAJeqqorIqUhdDk8f62eq+mwfvRciOoa4deshOpp1d19YLr76UavqKgCrOq1bmvX3PQDuybHf2wDKehwVEVlt4MCB2LFjB4YNG8ZiTccMVcWOHTswcODAHu3nq1ATEQWpsLAQLS0taGtrMx0KUa8aOHAgCgsLe7QPCzURhU7//v0xatQo02EQhQLH+iYiIgoxFmoiIqIQY6EmIiIKMRZqIiKiEGOhJiIiCjEWaiIiohBjoSYiIgoxFmoiIqIQY6EmIiIKMRZqIiKiEGOhJiIiCjEWaiIiohBjoSYiIgoxFmoiIqIQY6EmIiIKMRZqIiKiEGOhJiIiCjEWaiIiohBjoSYiIgoxFmoiIqIQY6EmIiIKMV+FWkSmichmEdkiIrfmeH6WiGwUkQ0iEheRT/vdl4iIwq+oCBDhI/0oKgou98d1t4GI5AFYAuCzAFoArBORFar6WtZmvwKwQlVVREoBPAZgrM99iYgo5JqbAVXTUYSHSHDH8vONugrAFlV9W1X3AWgEMCt7A1Xdo5r5XzgIgPrdl4iIiNz5KdQjALyTtdzSvs5BRD4vIm8A+AWAf+/JvkREdHQpWFwAqRdIvaAyVgkAiDZFM+ukXpDcnUTT5ibHulgiBgCOdbXLawEAtctrHesBIJaIOdY1bW5CcnfSsS7aFAUAVMYqM+sKFhcAAOpW1zm2TSQTSCQTjnV1q+u6fU9F3y4KKrVdiHZzLUNELgYwVVWvaV++AkCVqt7gsv15AL6uqlN6sq+IRAFEAeCMM86obG5uPoK3RUREvUmk49J33eo61FXXGY0naFIv0EUd9TI7H73y+iIJVY3kes7PN+oWAKdnLRcCSLptrKovAjhTRE7uyb6qGlPViKpGhg8f7iMsIiIyoX5NvekQAte6sNXYsf0U6nUARovIKBEZAGAOgBXZG4jIv4qkfloXkQoAAwDs8LMvERFR2CWSCWPH7vaub1U9ICLXA3gOQB6Ah1T1VRFZ0P78UgBfAHCliOwH8BGAS9tvLsu5bx+9FyIioj4xs3Gm49J3kLot1ACgqqsArOq0bmnW3/cAuMfvvkREdPSKz4+bDsEqHJmMiIgoxFioiYioRyLLct6cfExrqGkwdmwWaiIiom5EK6PGjs1CTURE1I30ACwmsFATEVGPLJq0yHQIVmGhJiKiHrFtVDLTWKiJiKhH0uNo26RmTI2xY7NQExFRj2zbs810CIFrmttk7Ngs1ERERN1Iz/BlAgs1ERH1SEV+hekQArfyzZXGjs1CTUREPZKImpugwkYs1ERE1CPRJnODf9iIhZqIiHpk2fplpkMInKmZswAWaiIiom7FEjFjx2ahJiIi6sa1K681dmwWaiIi6pHWha2mQ7AKCzUREfVIIsm7voPEQk1ERD0ys3Gm6RACt2LOCmPHZqEmIiLqRmVBpbFjs1ATERF1Y8T9I4wd+zhjRyYiOsrEEjHH3b8r5qxAZUGl40N8fsV8xGpjqIxVYv229QCA/MH5SN6URN3qOtSvqc9sG58fBwBElkUy6xZNWoS66joULC7ITH5RkV+BRDSBaFPU0Ye5dWErEsmE41J0Q00DopVRSL1k1tWMqUHT3CbULq91DIWpi9TzPY0cOhJbv7y1Sx4aahr8J42OmKia68TtJhKJaDweNx0GEVGG1IvRQS9MyH7PIkAIy0VgOv//7+18iEhCVSO5nuOlbyIiyml+xfyc67O/rdvCLRdB8FWoRWSaiGwWkS0icmuO578oIhvbH78TkbKs57aKyJ9EZIOI8GsyEdFRIlZrbjSusDGZi24LtYjkAVgCYDqAcQDmisi4Tpv9BcAkVS0F8A0And/R+apa7va1nogo7GrG1JgOIXCVMXN3OoeNyVz4+UZdBWCLqr6tqvsANAKYlb2Bqv5OVd9vX3wZQGHvhklEZFbT3CbTIQQufTNcZzY2WtxyEQQ/hXoEgHeyllva17n5DwDPZC0rgOdFJCEinBuNiI5KtctrTYcQGjY2WkzyU6hz3TWQ8143ETkfqUL9tazVE1W1AqlL518SkfNc9o2KSFxE4m1tbT7CIiLqO0VFqTt704+Vb650LNvwwO78jr+z2NhoyR+cb+zYfgp1C4DTs5YLASQ7byQipQAeBDBLVXek16tqsv2/2wE8hdSl9C5UNaaqEVWNDB8+3P87ICLqA83Nqe436QfgXLbicV/S8f7Tsvti2yJ5U5eyFxg/hXodgNEiMkpEBgCYA8Ax6KmInAHg5wCuUNU3s9YPEpEh6b8BXAhgU28FT0REfadudZ3pEELDZC66LdSqegDA9QCeA/A6gMdU9VURWSAiC9o3+zqAYQC+36kb1qkA1orIKwD+AOAXqvpsr78LIqI+ZttgJwAco6jZzmQufA0hqqqrAKzqtG5p1t/XALgmx35vAyjrvJ6I6GgTS8QQreT9sICdjRaTODIZEZEP2eNh2y6W4EAoQWKhJiKinNKThnRmY6PFLRdBYKEmIiIKMRZqIiIfVsxZ0f1Gx5js6TdtZzIXLNRERD5UFnDc6zQbGy0msVATEfkw4n6vkZPtwkZLsFioiYgop0WTFuVcb2OjxS0XQWChJiKinOqq60yHEBomc8FCTUTkw/yK+aZDCFzB4gLTIYSGyVywUBMR+RCrtW+Qj217tuVcb2OjxS0XQWChJiLyoTLGG6jSbGy0mMRCTUTkw/pt602HELiK/Iqc621stLjlIggs1ERElFMimsi53sZGi1sugsBCTUTkQ/7gfNMhBC7axNnC0kzmgoWaiMiH5E1J0yEEbtn6ZTnX29hocctFEFioiYh8qFtdZzqE0LCx0WISCzURkQ/1a+pNhxAabLQEi4WaiIhyal3YmnO9jY0Wt1wEgYWaiIhySiTN3ekcNiZzwUJNRORDfH7cdAiBm9k403QIoWEyFyzURETUIzY2WkxioSYi8iGyLGI6BLIUCzUREeXUUNOQc72NjRa3XASBhZqIiHKKVnJksjSTufBVqEVkmohsFpEtInJrjue/KCIb2x+/E5Eyv/sSER0NFk1aZDqEwEm9mA4hNEzmottCLSJ5AJYAmA5gHIC5IjKu02Z/ATBJVUsBfANArAf7EhGFXl11nekQQsPGRotJfr5RVwHYoqpvq+o+AI0AZmVvoKq/U9X32xdfBlDod18ioqNBweIC0yGEBhstwfJTqEcAeCdruaV9nZv/APBMT/cVkaiIxEUk3tbW5iMsIqLgbNuzzXQIgasZU5NzvY2NFrdcBMFPoc51YV5zbihyPlKF+ms93VdVY6oaUdXI8OHDfYRFRER9qWluU871NjZa3HIRBD+FugXA6VnLhQC6TJ0iIqUAHgQwS1V39GRfIqKwq8ivMB1C4GqX15oOITRM5sJPoV4HYLSIjBKRAQDmAFiRvYGInAHg5wCuUNU3e7IvEdHRIBG1b9zrlW+uzLnexkaLWy6C0G2hVtUDAK4H8ByA1wE8pqqvisgCEVnQvtnXAQwD8H0R2SAica99++B9EFEvKCoCRPiQHD/aRZvYpzjNxkaLScf52UhVVwFY1Wnd0qy/rwFwjd99iSicmpsBzXkXiX06F+tl65chVhszE0zIRJuizEWAODIZERHlpItyt9qWrV8WcCTmueUiCCzURESUUyzBb81pJnPBQk1E5EPrwlbTIQTu2pXXmg4hNEzmgoWaiMiHRJI3UKXZ2GgxiYWaiMiHmY0zTYcQGmy0BIuFmoiIcloxJ/ewFzY2WtxyEQQWaiIiyqmyoNJ0CKFhMhcs1EREPjTUNJgOIXAj7veaf8kuJnPBQk1E5EO0kiOTpdnYaDGJhZqIyAepzzUZoJ3YaAkWCzUREeU0v2J+zvU2NlrcchEEFmoiIsqJ43l3MJkLFmoiIh9qxtSYDiFwlTHe9Z1mMhcs1EREPjTNbTIdQuDWb1ufc72NjRa3XASBhZqIyIfa5bWmQwgNGxstJvmaj5qI7JN9w9D8ivmI1cZQGavMfLPIH5yP5E1J1K2uQ/2a+sy28flxAEBkWSSzbtGkRairrkPB4gJs27MNAFCRX4FENIFoU9QxbWLrwlYkkgnH6FcNNQ2IVkYdMdWMqUHT3CbULq/FyjdXZtbrIkUsEXNMorBizgpUFlQ6+sJ6vadcso9hi/zB+TnX1y6vta5Yu+UiCKIhnCU+EoloPB43HQaRdUSA9EdCcncSBUMKzAZkUHYugFTDxeScxKZl58P2XABdz48jfz1JqGok13O89E1EOdk48ULd6jrTIYQK89HBZC5YqIkoJxsnXsi+hN+Zjd8gvfJhG5O5YKEmIvIhlmCf4jQbGy0msVATZSn6dhEKFqd+l61bXQepl8wjkUwgkUw41qUvhxUsLsisS/e3jDZFHdsmdyfRtLnJsS794Z+9Ln13ce3yWsd6IFUsstc1bW5CcnfSsS7alBresTJWmVnn9Z7In+yb02zHRkuweDMZURbbbpJJJBOO6fuyb5CJJWLWjemcnQ/eTOaeD9tzAfBmMiJjKvIrTIcQqOwuVJ3ZVqSJwspXoRaRaSKyWUS2iMitOZ4fKyIvicg/ReTmTs9tFZE/icgGEeHXZAq1RJSXgtNsnHjBq+GyYs6KACMJB6982MZkLrot1CKSB2AJgOkAxgGYKyLjOm32HoAbAdzn8jLnq2q529d6orBI/75L1Fn2ZU/b2dhoMcnPN+oqAFtU9W1V3QegEcCs7A1UdbuqrgOwvw9iJApM9ghZNlg0aZHpEI4a2aOa2Y6NlmD5KdQjALyTtdzSvs4vBfC8iCREhF9XiEKkrrrO9TkbJ15gw8XJLR82NlpMnht+CnWuH6p6cq/bRFWtQOrS+ZdE5LycBxGJikhcROJtbW09eHkiOlzpblu52DaWM+DdcLER89HBZC78FOoWAKdnLRcCyD1qfQ6qmmz/73YATyF1KT3XdjFVjahqZPjw4X5fnqhXtS5sNR1CoNITZORi42xRXg2X+RXzA4wkHGrdlZgAABQvSURBVLzyYRuTufBTqNcBGC0io0RkAIA5AHzdSSAig0RkSPpvABcC2HS4wRL1NQ4A0sHG2aK8Gi6xWvsG+XDLh42NFq9zo691W6hV9QCA6wE8B+B1AI+p6qsiskBEFgCAiJwmIi0AFgK4U0RaRORjAE4FsFZEXgHwBwC/UNVn++rNEB0p28a3tq3f+JFIjzhHdjZaTPI1H7WqrgKwqtO6pVl//w2pS+KdfQCg7EgCJKK+w37jTl4Nl/Sc1TZxy0dlrNK6c8dko5YjkxFZzKvfuG1DRAJsuHTmlg8bGy0mzw0WaqIsDTUNpkMIlFe/cRsnXvBquOQPzg8wknDgAEAdTOaChZqsVlSUGlw//bg2EnUsH+sPoOtymo2zRXk1XJI3+e7scsxwy4eNjRaTgyGxUJPVmptTM+CkH6gTx/Kx/gC6LlNu6SlNyc5Gi0ks1EQWs63f+JGoX1NvOoTQYKMlWCzURBbz6jdu48QLbLg4ueXDxkaLyXODhZooi23jW3v1G7dx4gUOeOPEfHQwmQsWaqIsNo5v7cbGiRe8Gi7x+fEAIwkH2wYA8mIyFyzURFlsHN+aqKdsbLSYxEJNlMW28a1t6zd+JCLLIqZDIEuxUBNZLFrpPoiDjRMvsOHi5JYPGxstJs8NFmoii0l9runmU2yceMGr4WIj5qODyVywUBNlsXF8azc2zhbl1XBZNGlRgJGEg1c+bGMyFyzURFlsHN/ajY0TL3ipq64zHUJo2NhoMYmFmiiLbeNb29Zv/EgULC4wHUJosNESLBZqIot59Ru3ceIFr4bLtj3bAowkHNzyYWOjxWSjloWayGJe/cZtnHiBA944ueXDxkaLyXODhZooi23jW3v1G7dx4gWvhktFfkWAkYQDBwDqYDIXLNREWWwc39qNjRMveDVcElH7xr12y4eNjRaTgyGxUBNlsXF8a/In2sQ+xWk2NlpMYqEmshj7jfu3bP0y0yGEBhstwWKhJrKYV79xGydeYMPFyS0fNjZaTJ4bLNREWWwb39q2fuPd4YA3TsxHB5O58FWoRWSaiGwWkS0icmuO58eKyEsi8k8Rubkn+xKFiY3jW7uxceIFr4ZL68LWACMJBzbkOpjMRbeFWkTyACwBMB3AOABzRWRcp83eA3AjgPsOY1+i0LBxfGvyJ5HkDVRpNjZaTPLzjboKwBZVfVtV9wFoBDArewNV3a6q6wDs7+m+RGFi2/jWtvUbPxIzG2eaDiE02GgJlp9CPQLAO1nLLe3r/DiSfYmoj3n1G7dx4gU2XJzc8mFjo8XkueGnUOea28vv7W++9xWRqIjERSTe1tbm8+WJepdt41t79Ru3ceIFDnjjxHx0MJkLP4W6BcDpWcuFAPwOAux7X1WNqWpEVSPDhw/3+fJEvcvG8a3d2DjxglfDpaGmIcBIwoEDAHUwmQs/hXodgNEiMkpEBgCYA8DvNYAj2ZcocDaOb+3GxokXvEQrOchHmo2NFpO6LdSqegDA9QCeA/A6gMdU9VURWSAiCwBARE4TkRYACwHcKSItIvIxt3376s0QHSnbxre2rd/4kZD6XL/k2YmNlmAd52cjVV0FYFWndUuz/v4bUpe1fe1LROHg1W/cxokX2HBxcsuH1It1o7iZPDc4MhmRxbz6jds48QIHvHFiPjqYzAULNVEW28a39uo3buPEC14Nl5oxNQFGEg4cAKiDyVywUBNRTjZOvODVcGma2xRgJOHglg8bGy0mB0NioSbKYtv41rb1Gz8StctrTYcQGjY2WkxioSayGPuNO3k1XFa+uTLASMLBLR82NlpMNmpZqIks5tVv3MaJF9hwcXLLh42NFpPnBgu1hYqKABE+JEe3WNvGt/bqN27jxAsc8MaJ+ehgMhcs1BZqbgZU+dAc3UBtHN/ajY0TL3g1XGzrNwzYNwCQF5O5YKEmymLj+NbkTyzBPsVpNjZaTGKhtlxlrBJSL5B6yRSputV1mXVSL0gkE0gkE4516ctABYsLMuvS/QyjTVHHtsndSTRtbnKsS3/oZa9L36BSu7zWsR5IfUhmr2va3ITk7qRjXbrfr5/35Ma28a1t6zd+JK5dea3pEEKDjZZgiea6/mdYJBLReJwfIH1FJPdlX1skkonMlHWdc2Hb0IjZuQCc+YglYtaN6cxzw8ktH7bnAuj9z1ERSahqzv6h/EZtORtHn/Ji2/jWXv3GbSvSRGHFQm05G0ef8ipONo5v7cbG2aK8zo0Vc+ybode2AYC8mMwFCzVRFl5hIDfZlz1tZ2OjxSQWaqIstl1hsK3f+JEYcf8I0yGEBhstwWKhtpyNo0+xOHXw6jdu48QLPDec3PJhY6PF5LnBQm05K0ef4qAmGV79xm2ceIHnhhPz0cFkLlioLWfj6FNexcm2Kwxe/cZtnHjB69yYXzE/wEjCgQMAdTCZCxZqso5XcbLxCoMbGyde8Do3YrX2DfLhlg8bGy0mB0NioSbKYtsVBtv6jR+J9Mh7ZGejxSQWass11DSYDiFwLE4d2G/cyevcWL9tfYCRhINbPmxstJj83GChtpyNo0+xOHXw6jdu2xCRAM+NztzyYWOjxeS5wUJtORtHn/IqTrZdYfDqN27jxAte50b+4PwAIwkHDgDUwWQuWKjJOl7FycYrDG5snC3K69xI3pQMMJJwcMuHjY0Wk4Mh+SrUIjJNRDaLyBYRuTXH8yIiD7Q/v1FEKrKe2yoifxKRDSLCKbEo1Gy8wkD+pKd2JTsbLSZ1W6hFJA/AEgDTAYwDMFdExnXabDqA0e2PKIAfdHr+fFUtd5vCi8yxcfQp6mBbv/EjUb+m3nQIocFGS7D8fKOuArBFVd9W1X0AGgHM6rTNLAD/rSkvAzhRROy7NnIUsnH0KRanDl79xm2ceIHnhpNbPmxstJg8N/wU6hEA3slabmlf53cbBfC8iCRExPUHQBGJikhcROJtbW0+wqLeYOPoU17FybYrDF79xm2ceIED3jgxHx1M5sJPoc71o13nfhte20xU1QqkLo9/SUTOy3UQVY2pakRVI8OHD/cRFvUGG0ef8ipONl5hcGPjxAte50Z8vn232Ng2AJAXk7nwU6hbAJyetVwIoPOdBK7bqGr6v9sBPIXUpXSiULLxCgNRT9nYaDHJT6FeB2C0iIwSkQEA5gDo/OPVCgBXtt/9fS6AXaq6TUQGicgQABCRQQAuBLCpF+Mn6lW2XWGwrd/4kYgs472wZMZx3W2gqgdE5HoAzwHIA/CQqr4qIgvan18KYBWAfwOwBcCHAOa1734qgKdEJH2sn6nqs73+Luiw2Tj6FItTB69+4zZOvMBzw8ktH5FlEes+O0yeG6IavmRHIhGNx3lppa+IAOn/7bFEzOpBPrJzAaT6Udv0AdT5/XbOh81sPzc6y86H7bkAev/fiogk3Lowc2Qyy9k4+pTXoCa2f/hks3HiBa9zY9GkRQFGEg4cAKiDyVywUBNlsXF8azc2Trzgpa66znQIoWFjo8UkFmqiLLZdYbCt3/iRKFhcYDqE0GCjJVgs1JazcfQpFqcOXv3GbZx4wevc2LZnW4CRhINbPmxstJj83GChtpyNo09xUJMOXv3GbZx4geeGk1s+bGy0mDw3WKgtZ+PoU17FybYrDF79xm2ceMHr3KjIr3B97ljFAYA6mMwFCzVZx6s42XiFwY2NEy94nRuJqH3jXrvlw8ZGi8nBkKws1FIvmUe6lVS7vNaxHkjdAZy9rmlzE5K7k4510aZUH+TKWGVmXfr3m7rVdY5tE8kEEsmEY136W0vB4oLMunS3mGhT1LFtcncSTZubHOvSdyl7vaeibxcFldqjno1XGMif9L91srPRYpSqhu5RWVmpfaUh3tBnrx1WqINzOWtx/or5AUdjXnY+APfnbJSdDxtzwXPDyS0ftn9uqHY9P4749YG4utREK75RFxWlRpERSXW/Sf9tywN16lzOEqu1r98wBzXp4NVv3MaJF3huOLnlY9n6ZQFHYp7Jc8OKQt3cnBrqLT3cW/pvWx4N8ZhjOZuNo095FSfbxre2rd94dzjgjRPz0cFkLqwo1Lbz+jC2cfQpr3zYeIXBjY2zRXmdG60LWwOMJBzYkOtgMhfWFWrbut9Qz9h4hYH8SSR5A1WajY0Wk6wr1Ox+42Tj6FNebLvCwIarfzMbZ5oOITTYaAmWdYXaxu43Xh/GNo4+xeLUwavhauPECzw3nNzyYWOjxeS5YV2htpHXh7GNo0955cO2KwxeDVcbJ17gFTcn5qODyVywUFvA68PYxtGnvPJh4xUGNzZOvOB1bjTUNAQYSTjYeAXSjclcWFeobet+Qz1j4xUGNzZOvOAlWsmRydJsbLSYZF2hZvcb8mLbFQY2XP1LDy1MbLQEzbpCbWP3G68PYxtHn2Jx6uDVcLVx4gWeG05u+bCx0WLy3LCuUNvW/QbgVYTOmI8OXg1XGyde4LnhxHx0MJkL6wq1jbw+jG0cfcorH7ZdYfBquNo4W5TXuVEzpibASMLBxiuQbkzmwrpCbVv3G8DOqwhemA9/bJx4wevcaJrbFGAk4eCWDxsbLSY/N3wVahGZJiKbRWSLiNya43kRkQfan98oIhV+9w0au9+QF9uuMNjYcD1c6Xneyc5Gi0ndFmoRyQOwBMB0AOMAzBWRcZ02mw5gdPsjCuAHPdg3UMa637S1AevWpf4bMK8PYyOjTxnMBRDC4mQwH6FsuIb038rKN1cGGEm7kP5bMdZoCem50efcJqpOPwBMAPBc1vJtAG7rtE0DgLlZy5sB5PvZN9ejsrKylyfkzvrbxOTvP/uZ6gknqA4dmvrvz34WfAxZenvC8x4JeS4CPz8M52PRC4scy9n5aP2gNdBYVNV4PrLZfm50xs9Rp97+HAUQV5eaeJyPWj4CwDtZyy0AzvGxzQif+3axefNmVFdX+wjNv8zLbQWqX+jd1/a0fz/w8svAoUPARx+l1l1+OfC97wH9+wcSwtadW1F0YpFjXTofL73zEiacPiGQOMKQC6BrPrJPtZE7RwZ3foQgH2u2rsHqotWOdel87PhwB4b9y7BA4gAQinx4nRuBfnaEIBeARz628nMU6HR+9CE/hTpXhzn1uY2ffVMvIBJF6rI5jj/+eB9h+Xf88cCaNe0LI1MfThnb2n9Oz8+6UWDnSGBXEVD4EpC3L7Vu32BgWyUw7E1gcNaITS0TgAG7gVM2dazbMQbYkw+MbD/O6QA+OgnYXgKc8ifghPeA1t+lnmuelHq9YW927L+9GNg3JHX8tD35qdfNTwAD9qTWHRyQOv7QrcCJze7v6cDxaH6lyJGTjnzsCzAfv0vlAln5OPlPHbkwlI81WW8fQ4HmEy3Kx4Hjne8/Ox8jNzmfsPDfivPcGMl/K/wczTzdy2XKk59C3YKOUwYACgF0/mHLbZsBPvYFAKhqDEAMACKRiK5evdpHaEeBtjZg5Ejgo/cArAG2AzjhBKC5GRg+3HR0wcrk4iMA7fnYfQKwabt9uQCYj874b6UDzw0nC84NEfdBZPzc9b0OwGgRGSUiAwDMAdB5vq8VAK5sv/v7XAC7VHWbz32PbcOHAz/8Yeqk+tjHUv/94Q+PmZOrR5gLJ+bDifnowFw4WZ4PSf2G3c1GIv8G4NsA8gA8pKr/V0QWAICqLpVUU+B7AKYB+BDAPFWNu+3b3fEikYjG48fYwBNtbcDWrUBRkTUnlyvmwon5cGI+OjAXTsdwPkQkoao5+4f6KtRBOyYLNRERkQuvQm3dyGRERERHExZqIiKiEGOhJiIiCjEWaiIiohBjoSYiIgoxFmoiIqIQY6EmIiIKMRZqIiKiEAvlgCci0gagudsNjz4nA/i76SBCgrlwYj6cmI8OzIXTsZqPkaqac7i1UBbqY5WIxN1GnrENc+HEfDgxHx2YCycb88FL30RERCHGQk1ERBRiLNTBipkOIESYCyfmw4n56MBcOFmXD/5GTUREFGL8Rk1ERBRiLNQBEJGHRGS7iGwyHYtpInK6iLwgIq+LyKsi8n9Mx2SSiAwUkT+IyCvt+ag3HZNpIpInIn8UkZWmYzFNRLaKyJ9EZIOIxE3HY5KInCgiT4jIG+2fHxNMxxQUXvoOgIicB2APgP9W1WLT8ZgkIvkA8lV1vYgMAZAA8DlVfc1waEaIiAAYpKp7RKQ/gLUA/o+qvmw4NGNEZCGACICPqWqN6XhMEpGtACKqeiz2G+4REfkxgN+o6oMiMgDAv6jqTtNxBYHfqAOgqi8CeM90HGGgqttUdX3737sBvA5ghNmozNGUPe2L/dsf1raeRaQQwAwAD5qOhcJDRD4G4DwAPwQAVd1nS5EGWKjJIBEpAnAWgN+bjcSs9ku9GwBsB/C/qmpzPr4N4KsADpkOJCQUwPMikhCRqOlgDPoEgDYAD7f/LPKgiAwyHVRQWKjJCBEZDOBJAF9W1Q9Mx2OSqh5U1XIAhQCqRMTKn0dEpAbAdlVNmI4lRCaqagWA6QC+1P4zmo2OA1AB4AeqehaAfwC41WxIwWGhpsC1/xb7JIBHVPXnpuMJi/ZLeasBTDMciikTAcxs/122EcBkEfmp2ZDMUtVk+3+3A3gKQJXZiIxpAdCSdbXpCaQKtxVYqClQ7TdP/RDA66p6v+l4TBOR4SJyYvvfJwCYAuANs1GZoaq3qWqhqhYBmAPg16p6ueGwjBGRQe03XKL9Mu+FAKzsOaKqfwPwjoh8sn3VBQCsuQH1ONMB2EBElgOoBnCyiLQAWKSqPzQblTETAVwB4E/tv8sCwO2quspgTCblA/ixiOQh1XB+TFWt75ZEAIBTATyVatviOAA/U9VnzYZk1A0AHmm/4/ttAPMMxxMYds8iIiIKMV76JiIiCjEWaiIiohBjoSYiIgoxFmoiIqIQY6EmIiIKMRZqIiKiEGOhJiIiCjEWaiIiohD7/xwTNV/bQ5tGAAAAAElFTkSuQmCC\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "# nbi:hide_in\n", "# library\n", "import matplotlib.pyplot as plt\n", "import numpy as np\n", "plt.rcParams[\"figure.figsize\"] = (8, 5)\n", "\n", "n = 6\n", "num_samples=1000\n", "\n", "range_x = np.arange(1,n+1)\n", "pmf_values = np.array([(2*i-1)/n**2 for i in range(1,n+1)])\n", "\n", "# generate data\n", "toss = np.random.randint(1,n+1,(2,num_samples))\n", "data = np.amax(toss, axis=0).squeeze()\n", "\n", "# compute empirical pmf\n", "def epmf(data):\n", " erange_x, counts = np.unique(data, return_counts=True)\n", " epmf_values = counts/data.size\n", " return epmf_values, erange_x\n", "\n", "epmf_values, erange_x = epmf(data)\n", "\n", "# plot \n", "plt.ylim(-0.02,0.4) \n", "plt.axhline(y=0, color='k')\n", "plt.xticks(range_x)\n", "\n", "plt.scatter(range_x,np.zeros(range_x.shape), color =\"red\", s=20)\n", "plt.bar(range_x, pmf_values, width=1, color=(1, 1, 1, 0), edgecolor='blue', linewidth=1, label=\"True histogran\")\n", "plt.bar(erange_x, epmf_values, width=0.9, color=(1, 1, 1, 0), edgecolor='green', linewidth=1,linestyle=\"--\", label=\"Relative frequency histogram\")\n", "plt.legend()\n", "plt.show();" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "___\n", "### Empirical cdf\n", "\n", "Similar to empirical pmf, the empirical cfd is computed as\n", "\n", "$$F_{\\text{emp}}(x)=\\frac{\\text{number of outomes }\\leq x}{\\text{size of data}}.$$\n", "\n", "We consider the double die toss experiment. Let the random variable $X$ be the maximum of two tosses. Again, $\\text{Range}(X)=\\{1,\\dots, 6\\}$ and it can be checked that\n", " $$F(x)=\n", "\\begin{cases}\n", "0 & x<1\\\\\n", "\\sum\\limits_{i=1}^k\\frac{2i-1}{6^2}=\\frac{k^2}{6^2} & k\\leq x" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "# nbi:hide_in\n", "# library\n", "import matplotlib.pyplot as plt\n", "import numpy as np\n", "plt.rcParams[\"figure.figsize\"] = (8, 5)\n", "\n", "n = 6\n", "num_samples=100\n", "range_x = np.arange(1,n+1)\n", "\n", "# generate data\n", "toss = np.random.randint(1,n+1,(2,num_samples))\n", "data = np.amax(toss, axis=0).squeeze()\n", "\n", "# compute true cdf values\n", "cdf_values = np.array([i**2/n**2 for i in range(1,n+1)])\n", "\n", "# compute empirical cdf values\n", "def ecdf(data):\n", " erange_x, counts = np.unique(data, return_counts=True)\n", " cdf_emp = np.cumsum(counts)/data.size\n", " return cdf_emp, erange_x\n", "\n", "ecdf_values, erange_x = ecdf(data)\n", "\n", "# add padding \n", "def padding(cdf_values, range_x):\n", " edge = (range_x[0]-2, range_x[-1]+2)\n", " range_x_padded = np.pad(range_x, (1,1), 'constant', constant_values=edge)\n", " cdf_values_padded = np.pad(cdf_values, (1,1), 'constant', constant_values=(0, 1))\n", " return cdf_values_padded, range_x_padded\n", "\n", "cdf_values_padded, range_x_padded = padding(cdf_values, range_x)\n", "ecdf_values_padded, erange_x_padded = padding(ecdf_values, erange_x)\n", "\n", "# plot setup\n", "plt.ylim(-0.07,1.1)\n", "plt.xlim(-1.2,8.2)\n", "plt.axhline(y=0, color='k')\n", "plt.xticks(range_x)\n", "\n", "# plot cdf using step function\n", "plt.scatter(range_x,np.zeros(range_x.shape), color =\"red\", s=20)\n", "plt.step(range_x_padded, cdf_values_padded, where='post', color=\"blue\", linewidth=1, label=\"True cdf\")\n", "plt.step(erange_x_padded, ecdf_values_padded, where='post', color=\"green\", linewidth=1, linestyle=\"--\", label='Empirical cdf')\n", "plt.legend()\n", "\n", "plt.show();\n" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "___\n", "### Hypergeometric distribution" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Suppose $N$ balls in an urn, $K$ of which are red, the rest are blue. $n$ balls are selected without order and without replacement. \n", " Let $S$ be the set of all such selections. Consider the following random variable:for $s\\in S$, define $X(s)$ to be the number of red balls in $s$. \n", " \n", "The histogram of this random variable is shown below. Red circles represent the range of $X$. You can change the values of the paramaters by playing with the slider. " ] }, { "cell_type": "code", "execution_count": 1, "metadata": { "scrolled": false }, "outputs": [ { "data": { "application/vnd.jupyter.widget-view+json": { "model_id": "f8451098e2b847759614fae57561fa0c", "version_major": 2, "version_minor": 0 }, "text/plain": [ "interactive(children=(FloatSlider(value=80.0, description='N', min=1.0, readout_format='', step=1.0), FloatSliā€¦" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "# nbi:hide_in\n", "# library\n", "import matplotlib.pyplot as plt\n", "import numpy as np\n", "from scipy.special import comb\n", "from ipywidgets import interact, FloatSlider\n", "\n", "\n", "\n", "def hypergeometric_pmf(N, K, n):\n", " # the range of the hyper-geometric distribution\n", " range_x = np.arange(max(0, n-(N-K)), min(n, K)+1, 1)\n", " \n", " # compute the pmf at value i\n", " def hyper_pmf(N,K,n,i):\n", " pmf_val = comb(K, i, exact=True) * comb(N-K, n-i, exact=True) / comb(N, n, exact=True)\n", " return pmf_val\n", " pmf_values = np.array([hyper_pmf(N,K,n,i) for i in range_x])\n", "\n", " # plot setup\n", " plt.figure(figsize=(14,7)) \n", " plt.axhline(y=0, color='k')\n", " plt.xlim(-2,N+2)\n", " plt.xticks(np.arange(0, N+1, 5))\n", "\n", " # plotting with plt.bar instead of plt.hist works better when f(x) are knowwn\n", " plt.scatter(range_x,np.zeros(range_x.shape), color =\"red\", s=20)\n", " plt.bar(range_x, pmf_values, width=1, color=(1, 1, 1, 0), edgecolor='blue', linewidth=1.3, label=\"Histogran\")\n", " plt.plot();\n", " \n", "# create interactive variables\n", "N = FloatSlider(min=1.0, max=100.0, step=1.0, value=80, readout_format='')\n", "K = FloatSlider(min=1.0, max=N.value, step=1.0, value=40, readout_format='')\n", "n = FloatSlider(min=1.0, max=N.value, step=1.0, value=30, readout_format='')\n", "\n", "# enforce K<=N and n<=N\n", "def update_range(*args):\n", " K.max = N.value\n", " n.max = N.value\n", "N.observe(update_range, 'value')\n", "\n", "# display the interactive plot\n", "interact(hypergeometric_pmf, N=N, K=K, n=n);\n", "\n", "\n" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "\n", "**Problem (1.2-15 in the textbook):**\n", " \n", "
\n", "\n", "\n", "Five cards are selected at random without replacement from a standard, thoroughly shuffled 52-card deck of playing cards. Let X equal the number of face cards (kings, queens, jacks) in the hand. Forty observations of $X$ yielded the following data:\n", " \n", "$$\\text{Data}=\\{2, 1, 2, 1, 0, 0, 1, 0, 1, 1, 0, 2, 0, 2, 3, 0, 1, 1, 0, 3,\n", "\t\t\t\t\t\t\t1, 2, 0, 2, 0, 2, 0, 1, 0, 1, 1, 2, 1, 0, 1, 1, 2, 1, 1, 0.\\}$$\n", " \n", "1. Determine the pmf of $X$.
\n", "2. Draw a probability histogram for this distribution.
\n", "3. Determine the relative frequencies of 0, 1, 2, 3, and superimpose the relative frequency histogram on your probability histogram. \n", " \n", "
\n", " " ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "**Solution:** \n", "\n", "1. $X$ has Hypergeometric distribution with parameters $(N=52, K=12, n=5)$. $\\text{Range}(X)=\\{0,\\dots, 5\\}$ and \n", "\n", "$$f(x)=\\frac{{12\\choose x}{40 \\choose 5-x}}{{52 \\choose 5}}.$$\n", "\n", "2. Values of pmf are computed below.\n", "\n", "3. We compute the frequency pmf values and display the histogram. $\\blacksquare$\n" ] }, { "cell_type": "code", "execution_count": 2, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "True pmf = [2.53181273e-01 4.21968788e-01 2.50900360e-01 6.60264106e-02\n", " 7.61843199e-03 3.04737280e-04]\n" ] } ], "source": [ "# nbi:hide_in\n", "import numpy as np\n", "from scipy.special import comb\n", "\n", "N = 52\n", "K = 12\n", "n = 5\n", "\n", "# compute the range of X\n", "range_x = np.arange(max(0, n-(N-K)), min(n, K)+1, 1)\n", "\n", "def hyper_pmf(N,K,n,i):\n", " pmf_val = comb(K, i, exact=True) * comb(N-K, n-i, exact=True) / comb(N, n, exact=True)\n", " return pmf_val\n", "\n", "pmf_values = np.array([hyper_pmf(N,K,n,i) for i in range_x])\n", "print(\"True pmf = \",pmf_values)" ] }, { "cell_type": "code", "execution_count": 3, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "array([0.325, 0.4 , 0.225, 0.05 ])" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "# nbi:hide_in\n", "def epmf(data):\n", " erange_x, counts = np.unique(data, return_counts=True)\n", " epmf_values = counts/data.size\n", " return epmf_values, erange_x\n", "\n", "data = np.array([2, 1, 2, 1, 0, 0, 1, 0, 1, 1, 0, 2, 0, 2, 3, 0, 1, 1, 0, 3,\n", " 1, 2, 0, 2, 0, 2, 0, 1, 0, 1, 1, 2, 1, 0, 1, 1, 2, 1, 1, 0])\n", "\n", "epmf_values, erange_x = epmf(data)\n", "display(epmf_values)" ] }, { "cell_type": "code", "execution_count": 4, "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAlMAAAEvCAYAAABhSUTPAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjMsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+AADFEAAAgAElEQVR4nO3dfXRUhZ3/8c+XiOIRihWfErCE+pMqhIDJiHJoC7UouibYWnzAWh9+KwHPT90Wa+tDLWF7tkfdYru2thJs7W7rAZ+qS6jWblvB2q2tCaKLz9aChlBNUZBspTx9f38kjEkIc+cyN9zcue/XOTlm7kxmvvkS4Zvv3PmMubsAAACwbwbEXQAAAECSMUwBAAAUgGEKAACgAAxTAAAABWCYAgAAKADDFAAAQAEOiOuBDz/8cC8vL4/r4QEAAPLW3Nz8V3c/orfrYhumysvL1dTUFNfDAwAA5M3M1u3tOp7mAwAAKADDFAAAQAEYpgAAAAoQ2zlTAIDk2r59u1paWrR169a4SwEiNWjQII0YMUIDBw7M+2sYpgAAobW0tGjIkCEqLy+XmcVdDhAJd9fGjRvV0tKiUaNG5f11PM0HAAht69atGjZsGIMUioqZadiwYaE3rgxTAIB9wiCFYrQvP9cMUwCAxNm4caMmTJigCRMm6Oijj9bw4cOzl7dt2xbZ49x111364he/2Ot106dP15YtW/b6tbfddhvnlKUEwxQAoGDl5ZJZdB9Bb5AxbNgwrV69WqtXr9bcuXP1pS99KXv5wAMPlNRx/suuXbv67Ht+7LHHNGTIkL1eH+UwtWPHjkjuB32DYQoAULB16yT36D7W7TVrOrfXXntNFRUVmjt3rqqqqvTmm2/q0EMPzV6/dOlSXX755ZKkt956S+ecc44ymYwmTpyop556qtf7bGlp0fTp03Xcccfp+uuvzx4fMWKENm3apC1btujMM8/U+PHjVVFRoQceeEDf/va39fbbb+sTn/iEpk2bJkn66U9/qnHjxqmiokI33HBD9n4WLVqk0aNHa+rUqbr88suzm7CLLrpI11xzjT71qU/phhtu0FNPPaVJkybpxBNP1OTJk/Xqq69K6tiezZw5s9casX/waj6EUl6+73/JoXcjR0pr18ZdBVA8XnjhBd1999268847c250rr76an3lK1/RKaecorVr16qmpkZr1qzZ43bPPvusVq1apQMOOECjR4/WVVddpbKysuz1jzzyiMrLy/Xoo49KkjZv3qyhQ4dq4cKF+u1vf6tDDz1ULS0t+trXvqampiYNHTpU06ZN0/LlyzV+/HjdfPPNWrVqlQ455BBNnTpVEydOzN73n/70J/3617/WgAEDtHnzZj355JMqKSnRL37xC33ta1/Tvffem1eN6FsMUwhl92+fiA7n8ALROvbYY3XSSScF3u5Xv/qVXn755ezld999V++//74OPvjgbrebNm1a9um8448/Xm+88Ua3QaWyslLXXXedrrvuOtXW1mry5Ml7PNYf/vAHnXrqqTr88MMlSRdeeKGeeOIJbd26Vaeeeqo+/OEPS5JmzpypN954I/t15557rgYM6HgSadOmTbr44ov1pz/9aY/7D6oRfYun+QAAReWQQw7Jfj5gwAB5l98Au57D5O764x//mD3Xav369XsMUpJ00EEHZT8vKSnZY9t1wgknqKmpSWPHjtW1116rb37zm3vch+/lt9C9He/te7nxxhs1ffp0rVmzRg8//HC37yWoRvQthikAQNEaMGCAPvzhD+vVV1/Vrl279NBDD2WvmzZtmu64447s5dWrV+/TY6xfv16DBw/WF77wBc2bN0+rVq2SJA0ZMiT7ar9TTjlFjz/+uDZu3KgdO3Zo6dKlmjJlik4++WQ9/vjj2rRpk7Zv366f/exne32czZs3a/jw4ZKkH//4x/tUK/oGwxQAoKjdcsstOuOMM/TpT39aI0aMyB6/44479Lvf/U6VlZUaM2aMFi9evE/3/+yzz+qkk07ShAkTdOutt2ZPLq+rq9O0adM0bdo0jRgxQv/8z/+sqVOnasKECTrllFN01lln6SMf+YiuvfZaTZw4UaeffrrGjh2roUOH9vo4X/3qV3Xttdf2+jQi4mVBK8a+kslkvKmpKZbHxr4z45ypqNFTJNGLL76oE044IXs56henpOmFGe3t7Ro8eLC2b9+us88+W1dccYVqa2vjLivVev58S5KZNbt7prfbcwI6AKBgaRl8+sJNN92kFStWaOvWrTrjjDNUU1MTd0kIiWEKAIAYffvb3467BBSIc6YAAAAKkNcwZWZnmNnLZvaamV2X43YnmdlOM5sZXYkAAAD9V+AwZWYlku6QdKakMZJmmdmYvdzuFkmPRV0kAABAf5XPZmqipNfc/XV33yZpqaSze7ndVZIelPR2hPUBAAD0a/kMU8MlvdnlckvnsSwzGy7ps5LujK40AAD2rqSkRBMmTFBFRYVqa2u1adOmwK8ZPHhwzus3bdqk73//+9nLra2tmjkzmjNXrr322mxKerG49NJL9cADD+xxPKhvPfucdPkMU729c1jPVJzvSPqqu+/MeUdmdWbWZGZNbW1t+daIfqxsYcd7P9WvqJctsOxHc2uzmlubux2rX1Gf/Zrdx6obqiVJdY113W7buqVVjS83djvW0NwgSd2O1S7pyGKpXVLb7bgkNTQ3dDvW+HKjWre0djtW11gnSapuqM4e6+vvCShWffH/cS4HH3ywVq9erTVr1uiwww7rlma+r3r+I19WVtbrsLAvFi1apFWrVulf//Vfux0vxrd+Cepb1MPUzp05x4++5+45PyRNkvRYl8vXS7q+x23+LGlt50e7Op7q+0yu+62urnYkj9Tjcr16vyF61Vu/evYUSIIXXnih2+WqRVX7vYZDDjkk+/kPfvADv+KKK7KXb731Vs9kMj5u3Dj/+te/vsfXbNmyxU899VQ/8cQTvaKiwh9++GF3dz///PN90KBBPn78eP/yl7/sf/7zn33s2LHu7j5x4kRfs2ZN9r6mTJniTU1N3t7e7pdddplnMhmfMGFC9r66qq2t9QEDBvj48eN96dKlfskll/iXvvQlnzp1qs+bN2+v9/G3v/3Nzz//fB83bpyfd955PnHiRH/66af3+P7vv/9+v+SSS9zd/e233/ZzzjnHM5mMZzIZf/LJJ93dff78+X7ZZZf5lClTfNSoUf5v//Zv2a//93//dx83bpxXVlb6RRdd5O+9956Xl5f7tm3b3N198+bNPnLkyOzl3S655BK/6qqrfNKkST5q1Ci///773d279W3NmjV+0kkn+fjx433cuHH+yiuv7NHnXbt2+Ze//GUfO3asV1RU+NKlS93dfefOnX7FFVf4mDFj/KyzzvIzzzwz+xgjR470BQsW+OTJk33JkiXe0NDgmUzGKysr/ZxzzvH//d//zdY4d+5cnzp1qo8aNcpXrFjhl112mR9//PHZnvXU8+fb3V1Sk+9tVtrbFf7BoHSApNcljZJ0oKRnJY3NcfsfS5oZdL8MU8nEMFUYhikUi57/2ET9d8GipkWBt9k9TOzYscNnzpzpjz76qLu7P/bYYz579mzftWuX79y508866yxfuXJlt6/Zvn27b9682d3d29ra/Nhjj/Vdu3Z1GwLcuw8Ft912W3Ywa21t9eOOO87d3a+//nr/yU9+4u7u7777rh933HHe3t6+13rdO/6BP+uss3zHjh0572PhwoV+2WWXubv7s88+6yUlJYHD1KxZs/y3v/2tu7uvW7fOjz/+eHfvGKYmTZrkW7du9ba2Nj/ssMN827ZtvmbNGh89erS3tbW5u/vGjRvd3f3SSy/1hx56qOPPY9Einzdv3h7f0yWXXOIzZ870nTt3+vPPP+/HHnvsHn278sor/ac//am7u//973/3v/3tb3v0+YEHHvBp06b5jh07/C9/+Ysfc8wx3tra6vfff7+feeaZvnPnTt+wYYMfeuih3YapW265JXsff/3rX7Of33jjjX777bdnazz//PN9165d/vDDD/uQIUP8ueee8507d3pVVZU/88wze3xfYYepwKf53H2HpCvV8Sq9FyXd5+7Pm9lcM5sbxXYMyVVVWhV3CQCK0JzlcwJv8/7772vChAkaNmyY3nnnHZ122mmSpF/+8pf65S9/qRNPPFFVVVV66aWX9Oqrr3b7WnfXDTfcoMrKSk2bNk3r16/XW2+9lfPxzjvvPN1///2SpPvuu0/nnntu9vFuvvlmTZgwQVOnTtXWrVv1xhtvBNZ/7rnnqqSkJOd9PPHEE7roooskSZWVlaqsrAy831/96le68sorNWHCBM2YMUPvvfde9g2XzzrrLB100EE6/PDDdeSRR+qtt97Sb37zG82cOVOHH364JOmwww6TJF1++eW6++67JUl33323Lrvssl4f7zOf+YwGDBigMWPG9NrDSZMm6Zvf/KZuueUWrVu3TgcffPAet3nyySc1a9YslZSU6KijjtKUKVP09NNP68knn9S5556rAQMG6Oijj9anPvWpbl93/vnnZz9fs2aNPvGJT2jcuHG655579Pzzz2evq62tlZlp3LhxOuqoozRu3DgNGDBAY8eO1doI4vvzSkB390ckPdLjWK8nm7v7pQVXhcRormuOu4REWVSzKO4SgD5ROrh0vz/m7nOmNm/erJqaGt1xxx26+uqr5e66/vrrNWfO3geye+65R21tbWpubtbAgQNVXl6urVu35ny84cOHa9iwYXruued07733atGijv+f3V0PPvigPvaxj4Wq/5BDDsl+nus+zHo/f6zr8a6179q1S7///e97HVoOOuig7OclJSXasWOH3L3Xx5g8ebLWrl2rlStXaufOnaqoqOi1jq736b280eiFF16ok08+WT//+c81ffp03XXXXfroRz/a7Ta9fV2u47t17eGll16qhx9+WOPHj9ePf/xjrVixYo8aBwwY0K3eAQMGRHLOGgnoKMjuE7iRn7pq+oXi1HpNa2yPPXToUN1+++361re+pe3bt2v69On60Y9+pPb2dknS+vXr9fbb3VN7Nm/erCOPPFIDBw7U448/rnWd79I8ZMiQ7BanNxdccIFuvfVWbd68WePGjZMkTZ8+Xd/97nez//A/88wzob+Hvd3HJz/5Sd1zzz2SOjYvzz33XPZrjjrqKL344ovatWuXHnrooezx008/Xd/73veyl1evXp3zsT/96U/rvvvu08aNGyVJ77zzTva6iy++WLNmzdrrViofr7/+uj760Y/q6quv1owZM/Tcc8/t0edPfvKTuvfee7Vz5061tbXpiSee0MSJE/Xxj39cDz74oHbt2qW33nqr24DU05YtW1RaWqrt27dne7a/MEyhIItXLY67hETJ5xVKQBLtfmVrVJZdsCzU7U888USNHz9eS5cu1emnn64LL7xQkyZN0rhx4zRz5sw9BqTPf/7zampqUiaT0T333KPjjz9ekjRs2DBNnjxZFRUVvUYYzJw5U0uXLtV5552XPXbTTTdp+/btqqysVEVFhW666abQ3+/e7uOKK65Qe3u7Kisrdeutt2rixInZr7n55ptVU1OjU089VaWlH2wGb7/9djU1NamyslJjxozRnXfmTi0aO3asbrzxRk2ZMkXjx4/XvHnzuvXp3Xff1axZs0J/T7vde++9qqio0IQJE/TSSy/p4osv3qPPn/3sZ1VZWanx48fr1FNP1a233qqjjz5an/vc5zRixAhVVFRozpw5OvnkkzV06NBeH+cb3/iGTj75ZJ122mnZP8/9xYJWaH0lk8l4U1NTLI+NfWcmdf2RsQUmnx/Pz1AS9davnj0FkuDFF1/UCSeckL3c8xeFptkdf79nFmeyx+ZPma/6qfUqW1imDe0bJHWcd9lc16y6xrpuv5ytn7deZUPK+vJbSKSpU6fqW9/6ljKZTPCNI/DAAw/oP//zP/WTn/xkvzxeb9rb2zV48GBt3LhREydO1O9+9zsdffTRffqYPX++JcnMmt2918bndc4UAAC57O2Xqt6O9/aUYENtgxpqGyKvC/vuqquu0qOPPqpHHnkk+MZ9qKamRps2bdK2bdt000039fkgtS8YplCQ9fPWx11CotSMrom7BAAJluucoah997vf3W+Plcv+/J73FedMoSDNrbyaL4zGWSSgA0CxYZhCQWYsnRF3CYmy+20zgGIQ1zm3QF/al59rhilgP1r+yvK4SwAiMWjQIG3cuJGBCkXF3bVx40YNGjQo1NdxzhQAILQRI0aopaVFvGk9is2gQYM0YsSIUF/DMIWCkOgNpNPAgQM1atSouMsA+gWe5kNBSPQOh0wuACg+DFMoCIne4TQ0k6MDAMWGYQrYj+Ys3/sbrwIAkolhCgAAoAAMUygIid4AgLRjmEJBSPQOZ9kFy+IuAQAQMYYpFIRE73Cqy6rjLgEAEDGGKRSERO9wht82PO4SAAARY5gCAAAoAMMUAABAARimUBASvcOZXTU77hIAABFjmEJBSPQOp6GWfgFAsWGYQkFI9A6nuoFX8wFAsWGYAvajVRtWxV0CACBiDFMAAAAFYJhCQUj0Dqd0cGncJQAAIsYwhYKQ6B1O6zWtcZcAAIgYwxQKQqJ3OPUr6uMuAQAQMYYpYD9asHJB3CUAACLGMAUAAFAAhikUhERvAEDaMUyhICR6h9M0uynuEgAAEWOYQkFI9AYApB3DFApConc4mcWZuEsAAESMYQoAAKAADFMoCIneAIC0Y5hCQUj0Dmf+lPlxlwAAiBjDFApConc49VPr4y4BABCx1A1TjS83qnVLq2yBZT/qGuskdbwybfexsoVlkjqGha63bW5tVnNrc7djuweKsoVl2WO7X+VW11jX7batW1rV+HJjt2MNzR3xAl2P1S6plSTVLqntdlySGpobuh3r6+8pFxK9wwnqJwAgeczdY3ngTCbjTU37P3PHFph8fjzfcxL17JeZ1PVHhn6G01u/evYUAND/mFmzu/f6kuzUbaYAAACixDCFnKpKq3JeT6J3OEH9BAAkT+qGqUU1i+IuIVGa65rjLqGo0E8AKD6pG6bqquviLiFRdp/IvjckeocT1E8AQPKkbpja/Yo45GfxqsVxl1BU6CcAFJ/UDVMAAABRYphCQUj0BgCkXeqGqZrRNXGXkCjr563PeT2J3uEE9RMAkDypG6YaZzXGXUKiNLfmfvUZid7hBPUTAJA8eQ1TZnaGmb1sZq+Z2XW9XH+2mT1nZqvNrMnMPh59qdHY/TYtyM+MpTNyXr+hfcN+qqQ4BPUTAJA8BwTdwMxKJN0h6TRJLZKeNrNl7v5Cl5v9WtIyd3czq5R0n6Tj+6LgQi1/ZXncJQAAgCKSz2ZqoqTX3P11d98maamks7vewN3b/YM3+TtEEu80lhIkegMA0i6fYWq4pDe7XG7pPNaNmX3WzF6S9HNJ/7e3OzKzus6nAZva2tr2pV7sZ0GJ8SR6h0MCPwAUn3yGqd5SLvfYPLn7Q+5+vKTPSPpGb3fk7g3unnH3zBFHHBGu0oj4fJZmYQQlxpPoHQ4J/ABQfPIZplokHdPl8ghJrXu7sbs/IelYMzu8wNr6RENzQ9wlJEpQYjyJ3uGQwA8AxSefYeppSceZ2SgzO1DSBZKWdb2Bmf0fM7POz6skHShpY9TFRmHO8jlxlwAAAIpI4Kv53H2HmV0p6TFJJZJ+5O7Pm9nczuvvlPQ5SReb2XZJ70s6v8sJ6QAAAEUrcJiSJHd/RNIjPY7d2eXzWyTdEm1p6A+CEuNJ9A6HBH4AKD6pS0BfdsGy4BshKygxnkTvcEjgB4Dik7phqrqsOu4SEiUoMZ5E73BI4AeA4pO6YWr4bXtEZCEHEuOjRT8BoPikbpgCAACIEsMUCkKiNwAg7VI3TM2umh13CYkSlBhPonc4JPADQPFJ3TDVUEsCehhBifEkeodDAj8AFJ/UDVPVDbyaLwwS46NFPwGg+KRumFq1YVXcJQAAgCKSumEK0SLRGwCQdqkbpkoHl8ZdQqIEJcaT6B0OCfwAUHxSN0y1XtMadwmJEpQYT6J3OCTwA0DxSd0wVb+iPu4SEiUoMZ5E73BI4AeA4pO6YWrBygVxlwAAAIpI6oYpAACAKDFMIaegxHgSvcMhgR8Aik/qhqmm2U1xl5AoQYnxJHqHQwI/ABSf1A1TCCcoMZ5E73BI4AeA4pO6YSqzOBN3CYlCYny06CcAFJ/UDVMAAABRYphCTkGJ8SR6h0MCPwAUn6IfpsrLJbMPPrRifrfLfOT+2PDl1u7964FE73BI4AeA4lP0w9S6dZJ7l4/H67tf5iPnx/we/eqJRO9wSOAHgOJT9MNUT2ULy+IuIVFIjI/W3voZ9waymD7Ky/fvnykAHBB3AfvbhvYNcZcA7KG3rR/2TW9PRwNAX0rdZgrRItEbAJB2qRumqkqr4i4hUYIS40n0DocEfgAoPqkbpprrmuMuoaiQ6A0ASLvUDVN1jXVxl5AoQYnxJHqHQwI/ABSf1A1Ti1ctjrsEAABQRFI3TCFaJHoDANKOYQo5zZ8yP+f1JHqHE9RPAEDypG6YWj9vfdwlJEr91Prc15PoHUpQPwEAyZO6Yaq5lVfzhRGUGE9Cejgk8ANA8UndMDVj6Yy4S0gUEuOjRT8BoPikbpgCAACIEsMUcgpKjCfROxwS+AGg+KTujY4X1SyKu4REITE+Wvn0s/w75Tr92NPVUNug6obqbDBq6eBStV7TqvoV9d3OVds90HYNBJ0/Zb7qp9arbGFZ9qnFqtIqNdc1q66xrlve2vp569Xc2tztKfBFNYtUV10nW/DBuwbXjK5R46xG1S6p1fJXlmeP+3xXQ3OD5iyfkz227IJlqi6r1vDbhmePza6a3WffEwDEyTymt6vPZDLe1NT3Ww0zKaZvsSjUNdZ1e/+9nv20BSafT4Pz1bOfEj0tVNDPKABEwcya3b3Xt7FI3dN8XX/TRjAS46NFP6NHTwHELXXDFNDfkYUGAMnCMIWCkOgdPbLQACBZUjdM1YyuibuERAnakpDoHU4+Wyey0MJhkwcgbqkbphpnNcZdQqIEbUlI9A6HrVP06CmAuKVumKpdUht3CYkStCUh0Tsctk7Ro6cA4pa6YaprPg7QH5GFBgDJkrphCtEi0Tt6ddV1cZcAAAiBYQo5BW1JSJ8OJ5+tE1lo4bDJAxC3vIYpMzvDzF42s9fM7Lperv+8mT3X+fHfZjY++lKjQbJ0OEFbkrpGtihhsHWKHj0FELfAYcrMSiTdIelMSWMkzTKzMT1u9mdJU9y9UtI3JDWon2po7rel9UtBWxLSp8Nh6xQ9egogbvlspiZKes3dX3f3bZKWSjq76w3c/b/d/d3Oi09JGhFtmdHp+masQH9EFhoAJEs+w9RwSW92udzSeWxv/lHSo4UUBaQZWWgAkCz5DFO97dB7PfHIzD6ljmHqq3u5vs7Mmsysqa2tLf8qEZugLQnp0+Hks3UiCy0cNnkA4pbPMNUi6Zgul0dIau15IzOrlHSXpLPdfWNvd+TuDe6ecffMEUccsS/1FmzZBctiedykCtqSkD4dTj5bJ7LQwmGTByBu+QxTT0s6zsxGmdmBki6Q1G0iMbOPSPqZpC+4+yvRlxmd6rLquEtIlKAtCenT4bB1ih49BRC3wGHK3XdIulLSY5JelHSfuz9vZnPNbG7nzb4uaZik75vZajNr6rOKCzT8tlyne6EntiTRop/Ro6cA4nZAPjdy90ckPdLj2J1dPr9c0uXRlgakE1loAJAsJKCjIKRPR48sNABIltQNU7OrZsddQqIEbUlInw4nn60TWWjhsMkDELfUDVMNtfzWH0bQloT06XDYOkWPngKIW+qGqeoGXs0XBluSaNHP6NFTAHFL3TC1asOquEsAciILDQCSJXXDFKJF+nT0yEIDgGRJ3TBVOrg07hISJWhLQvp0OPlsnchCC4dNHoC4pW6Yar1mj3fCQQ5BWxLSp8Nh6xQ9egogbqkbpupX1MddQqIEbUlInw6HrVP06CmAuKVumFqwckHcJQA5kYUGAMmSumEK6O/IQgOAZGGYQk5BWxLSp8PJZ+tEFlo4bPIAxC11w1TT7Ka4S0iUoC0J6dPh5LN1IgstHDZ5AOKWumEK4QRtSUifDoetU/ToKYC4pW6YyizOxF1CorAliVY+/SQLLRx+RgHELXXDFNDfkYUGAMnCMIWcgrYkpE+Hk8/WiSy0cNjkAYhb6oap+VPmx11CogRtSUifDiefrRNZaOGwyQMQt9QNU/VT6+MuIVGCtiSkT4fD1il69BRA3FI3TJUtLIu7hERhSxIt+hk9egogbqkbpja0b4i7BCAnstAAIFlSN0whWqRPAwDSLnXDVFVpVdwlJErQloT06XDy2TqRhRYOmzwAcUvdMNVc1xx3CUWF9GkAQNqlbpiqa6yLu4RECdqSkD4dDlun6NFTAHFL3TC1eNXiuEsAciILDQCSJXXDFKJF+nT0yEIDgGRhmEJOQVsS0qfDyWfrRBZaOGzyAMTN3D2WB85kMt7U1PevwjGTun6LrVtaVTaEf6z2Vc9+1q+oZ5NSoJ49tQUmnx/P/5fFoGc/ASAKZtbs7r2epJm6zVRzK6/mCyNoS0L6dDhsnaJHTwHELXXD1IylM+IuIVFIjI9WPv0kCy0cfkYBxC11wxTQ35GFBgDJwjCFnIK2JKRPh5PP1okstHDY5AGIW+qGqUU1i+IuIVHYkkQrn36ShRYOP6MA4pa6Yaqumt/6wwjakpA+HQ5bp+jRUwBxS90wZQss7hIShS1JtOhn9OgpgLilbpgC+rv189bHXQIAIASGKRSE9OnokYUGAMmSumGqZnRN3CUkStCWhPTzcPLZOpGFFg6bPABxS90w1TirMe4SEiVoS0L6dDhsnaJHTwHELXXDVO2S2rhLSJSgLQnp0+GwdYoePQUQt9QNU8tfWR53CUBOZKEBQLKkbphCtEifjh5ZaACQLAxTyCloS0L6dDj5bJ3IQguHTR6AuKVumPL5HncJiRK0JSF9Ohy2TtGjpwDilrphqqG5Ie4SEiVoS0L6dDhsnaJHTwHELXXD1Jzlc+IuAciJLDQASJbUDVNAf0cWGgAkS17DlJmdYWYvm9lrZnZdL9cfb2a/N7O/m7oWAw0AAA2OSURBVNmXoy8TcQnakpA+HU4+Wyey0MJhkwcgboHDlJmVSLpD0pmSxkiaZWZjetzsHUlXS/pW5BVGbNkFy+IuIVGCtiSkT4eTz9aJLLRw2OQBiFs+m6mJkl5z99fdfZukpZLO7noDd3/b3Z+WtL0PaoxUdVl13CUkStCWhPTpcNg6RY+eAohbPsPUcElvdrnc0nkskYbfltjSY8GWJFr0M3r0FEDc8hmmenvd8T6FNZlZnZk1mVlTW1vbvtwFUPTIQgOAZMlnmGqRdEyXyyMkte7Lg7l7g7tn3D1zxBFH7MtdoJ8hfTp6ZKEBQLLkM0w9Lek4MxtlZgdKukBSYs/inl01O+4SEiVoS0L6dDj5bJ3IQguHTR6AuAUOU+6+Q9KVkh6T9KKk+9z9eTOba2ZzJcnMjjazFknzJH3NzFrM7EN9Wfi+aqjlt/4wgrYkpE+Hw9YpevQUQNzyyply90fcfbS7H+vu/9J57E53v7Pz87+4+wh3/5C7H9r5+Xt9Wfi+qm7g1XxhsCWJFv2MHj0FELfUJaCv2rAq7hKAnMhCA4BkSd0whWiRPh09stAAIFlSN0yVDi6Nu4RECdqSkD4dTj5bJ7LQwmGTByBuqRumWq/Zp1SH1ArakpA+HQ5bp+jRUwBxS90wVb+iPu4SEiVoS0L6dDhsnaJHTwHELXXD1IKVC+IuAciJLDQASJbUDVNAf0cWGgAkC8MUcgrakpA+HU4+Wyey0MJhkwcgbqkbpppmN8VdQqIEbUlInw4nn60TWWjhsMkDELfUDVMIJ2hLQvp0OGydokdPAcQtdcNUZnEm7hIShS1JtPLpJ1lo4fAzCiBuqRumgP6OLDQASBaGKeQUtCUhfTqcfLZOZKGFwyYPQNxSN0zNnzI/7hISJWhLQvp0OPlsnchCC4dNHoC4pW6Yqp9aH3cJiRK0JSF9Ohy2TtGjpwDilrphqmxhWdwlJApbkmjRz+j11lMzPqL8KC/f/3+uQJIcEHcB+9uG9g1xlwDkRBZa4Zws2UiZxV0B0L+lbjOFaJE+DQBIu9QNU1WlVXGXkChBWxLSp8PJZ+tEFlo4bPIAxC11w1RzXXPcJRQV0qcBAGmXumGqrrEu7hISJWhLQvp0OGydokdPAcQtdcPU4lWL4y4ByIksNABIltQNU4gW6dPRIwsNAJKFYQo5BW1JSJ8OJ5+tE1lo4bDJAxA385gCWTKZjDc19f2rcMy6Z860bmlV2RD+sdpXPftZv6KeTUqBevbUFph8PkFJ+6pnP8u/U651m9dJ6ngvyeqy6m7J/bOrZquhtkHVDdXZcwBLB5eq9ZpW1a+o7xYKuvuVg13P05o/Zb7qp9arbGFZNseuqrRKzXXNqmus63Zqwfp569Xc2qwZS2dkjy2qWaS66jrZgg/CnGpG16hxVqNql9Rq+SvLs8d9vquhuUFzls/JHuvL72lvPQXSyMya3b3XkzRTN0w1vtyo2o/V9vnjFouyhWXdtk/8w1+Ynv2U6Gmhgn5GG5obVFfNC08KwTAF5B6mUvc0X9ffCBGMxPho5dNPstDCCeopg1R4PNUMhJO6YQro78hCi1bXp8+QH36JAsJhmEJOQVsS0qfDyWfrRBZaOGzyAMQtdcPUoppFcZeQKGxJopVPP8lCC4ef0egxoALhpG6Y4vyJcIK2JKRPh8PWKXpBPa0ZXbOfKikeDKhAOKkbpjh/Ihy2JNGin9EL6mnjrMb9VEnxYOgHwkndMAX0d+vnrY+7hKJSu4QolLAY+oFwGKZQENKno9fcylMsUeoaegkAfSF1wxTnT4QTtCUh/TycfLZOZKGFwyYPQNxSN0xx/kQ4QVsSwv3CYesUPXoaPQZUIJzUDVOcPxFO0JaEcL9w2DpFL6invDVPeAyoQDipG6Y4fwL9HVlo0Wpoboi7hMRh6AfCSd0whWgR7hc9stCiNWf5nLhLAFDkGKaQU9CWhHC/cPLZOpGFFg6bPABxS90wxfkT4QRtSQj3C4etU/ToafQYUIFwUjdMcf5EOEFbEsL9wmHrFL2gni67YNl+qqR4MKAC4aRumOL8CfR3ZKFFq7qsOu4SEoehHwgndcMU0N+RhRat4bcNj7sEAEWOYQo5BW1JCPcLJ5+tE1lo4bDJAxC31A1TnD8RTtCWhHC/cPLZOpGFFg6bvOgxoALhpG6Y4vyJcIK2JIT7hcPWKXpBPZ1dNXs/VVI8GFCBcFI3THH+RDhsSaJFP6MX1NOGWl7BGxZDPxBO6oYpoL8jCy1a1Q1so8PqbUA14yOqj/Ly/f9nir6V1zBlZmeY2ctm9pqZXdfL9WZmt3de/5yZ8R4jKUG4X/TIQovWqg2r4i6hKLjzEdXHunVx/2kiaoHDlJmVSLpD0pmSxkiaZWZjetzsTEnHdX7USfpBxHVGpl+eP9HWJj39dMd/+5mgLUm/DPdLcD+lfpqFlvCe9jv9uJ+JRU+jR0/zls9maqKk19z9dXffJmmppLN73OZsSf/hHZ6SdKiZlUZcayT63fkTS5ZII0dKp53W8d8lS+KuqJugLUm/C/dLeD/7pYT3tHRwP/urqJ/3U0rggJqAniYOPQ3F3HP/T2NmMyWd4e6Xd17+gqST3f3KLrdZLulmd3+y8/KvJX3V3Zv2dr9Dhgzx6uq+P5dh5UppypQPLjdvaFZ1aT85h2L7dumpp6Rduz44NmCAdMop0sCB8dXVxcq1KzWl/IMG9uxnz+tjlcB+SvS0UEE/o/1KAvopSRu2bFDpkA+GUHoarX7dTymRPd0fVq5c2ezumd6uOyCPr+9t9dBzAsvnNjKzOnU8DaiDDjooj4cu3EEHdfzgZo1s18q1XQ5s6Dy9q7TLeRWbRkqby6URv5dKtnUc2zZY2lAtDXtFGrzhg9u2TJIO3CIdueaDYxtHS+2l0sguj/P+YdLb46Qj/0c6+J0Pju+a0nF/w17ZfUB65j1p25COx9+tvbTjfkubpQPbO47tPLDj8YeulQ7t8iR8lN/TjoO69089+yn6WWA/pR49LR1MT6P8Ge3Lx96tmPopSTsO0iurum/0evu57R8GSvpE90O7JP13HLXkp6SkP/dTSmJPR46M98T+fIapFknHdLk8QlLrPtxG7t4gqUGSMpmMr1ixIkytxaetreMnoP19qfPvSh18sLTuN9IRR8RaWiLRz+jR02jRz+jt7un7739w7OCDO87ypqf7hp72ymzvp7Xkc87U05KOM7NRZnagpAsk9YwRXybp4s5X9Z0iabO7b+h5R+jhiCOkH/6w44f0Qx/q+O8Pf5jqH9aC0M/o0dNo0c/o0dPo0dPQAs+ZkiQz+wdJ35FUIulH7v4vZjZXktz9TusY174n6QxJf5N0Wa7zpaSOzVRTU86bpEdbm7R2bceOkh/WwtHP6NHTaNHP6NHT6NHTbsxsr+dM5TVM9QWGKQAAkBS5hikS0AEAAArAMAUAAFAAhikAAIACMEwBAAAUgGEKAACgAAxTAAAABWCYAgAAKADDFAAAQAEYpgAAAArAMAUAAFAAhikAAIACMEwBAAAUILY3OjazNknrYnnw/ulwSX+Nu4giQj+jR0+jRT+jR0+jR08/MNLdj+jtitiGKXRnZk17ezdqhEc/o0dPo0U/o0dPo0dP88PTfAAAAAVgmAIAACgAw1T/0RB3AUWGfkaPnkaLfkaPnkaPnuaBc6YAAAAKwGYKAACgAAxTMTOzM8zsZTN7zcyui7uepDOzH5nZ22a2Ju5aioGZHWNmj5vZi2b2vJn9U9w1JZ2ZDTKzP5rZs509XRB3TcXAzErM7BkzWx53LcXAzNaa2f+Y2Woza4q7nv6Op/liZGYlkl6RdJqkFklPS5rl7i/EWliCmdknJbVL+g93r4i7nqQzs1JJpe6+ysyGSGqW9Bl+RvedmZmkQ9y93cwGSnpS0j+5+1Mxl5ZoZjZPUkbSh9y9Ju56ks7M1krKuDsZU3lgMxWviZJec/fX3X2bpKWSzo65pkRz9yckvRN3HcXC3Te4+6rOz7dIelHS8HirSjbv0N55cWDnB7/VFsDMRkg6S9JdcdeCdGKYitdwSW92udwi/qFCP2Vm5ZJOlPSHeCtJvs6npFZLelvSf7k7PS3MdyR9RdKuuAspIi7pl2bWbGZ1cRfT3zFMxct6OcZvqOh3zGywpAclfdHd34u7nqRz953uPkHSCEkTzYynpPeRmdVIetvdm+OupchMdvcqSWdK+n+dp1BgLxim4tUi6Zgul0dIao2pFqBXnef1PCjpHnf/Wdz1FBN33yRphaQzYi4lySZLmtF5js9SSaea2U/jLSn53L21879vS3pIHaelYC8YpuL1tKTjzGyUmR0o6QJJy2KuCcjqPFn6h5JedPfb4q6nGJjZEWZ2aOfnB0uaJumleKtKLne/3t1HuHu5Ov4O/Y27XxRzWYlmZod0vuBEZnaIpNMl8QrpHBimYuTuOyRdKekxdZzYe5+7Px9vVclmZksk/V7Sx8ysxcz+Me6aEm6ypC+o47f91Z0f/xB3UQlXKulxM3tOHb9Q/Ze783J+9CdHSXrSzJ6V9EdJP3f3X8RcU79GNAIAAEAB2EwBAAAUgGEKAACgAAxTAAAABWCYAgAAKADDFAAAQAEYpgAAAArAMAUAAFAAhikAAIAC/H8ewtS4jw682AAAAABJRU5ErkJggg==\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "# nbi:hide_in\n", "# library\n", "import matplotlib.pyplot as plt\n", "from scipy.special import comb\n", "plt.rcParams[\"figure.figsize\"] = (10, 5)\n", "\n", "# plot setup\n", "plt.axhline(y=0, color='k')\n", "plt.xticks(range_x)\n", "\n", "plt.scatter(range_x,np.zeros(range_x.shape), color =\"red\", s=20)\n", "plt.bar(range_x, pmf_values, width=1, color=(1, 1, 1, 0), edgecolor='blue', linewidth=1, label=\"True histogran\")\n", "plt.bar(erange_x, epmf_values, width=0.9, color=(1, 1, 1, 0), edgecolor='green', linewidth=1,linestyle=\"--\", label=\"Relative frequency histogram\")\n", "plt.legend()\n", "plt.show();" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [] } ], "metadata": { "kernelspec": { "display_name": "Python 3", "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.7.7" } }, "nbformat": 4, "nbformat_minor": 4 }