{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "# Bite Size Bayes\n", "\n", "Copyright 2020 Allen B. Downey\n", "\n", "License: [Attribution-NonCommercial-ShareAlike 4.0 International (CC BY-NC-SA 4.0)](https://creativecommons.org/licenses/by-nc-sa/4.0/)" ] }, { "cell_type": "code", "execution_count": 1, "metadata": {}, "outputs": [], "source": [ "import numpy as np\n", "import pandas as pd\n", "import matplotlib.pyplot as plt" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Estimating the difference between proportions\n", "\n", "From [the Reddit statistics forum](https://www.reddit.com/r/statistics/comments/fsieqf/q_distribution_of_difference_of_correlated_coin/)\n", "\n", "> You have two coins, and a sequence of flips from each. Your goal is to assess the distribution of the difference in coin flips. You could do Bayesian updating with a Beta prior to get a posterior for each coin, and assuming independence take the convolution. If you had the additional information that the average of the coin probabilities was p, how would you incorporate it? Is there any way to avoid the convolution and build a prior on the difference in coin flip probabilities directly? Thanks!\n", "\n", "I suggested:\n", "\n", "> Suppose the probability of heads for the two coins is p+x and p-x. If we know p, all we have to estimate is x.\n", ">\n", ">The likelihood function is the product of two binomials.\n", ">\n", ">I would use a grid method to estimate the posterior.\n", "\n", "Here's how:" ] }, { "cell_type": "code", "execution_count": 2, "metadata": {}, "outputs": [], "source": [ "xs = np.arange(101) / 100\n", "uniform = pd.Series(1, index=xs)\n", "uniform /= uniform.sum()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "We can use `binom.pmf` to compute the likelihood of the data for each possible value of $x$." ] }, { "cell_type": "code", "execution_count": 3, "metadata": {}, "outputs": [], "source": [ "from scipy.stats import binom\n", "\n", "def compute_likelihood(xs, p_avg, k, n):\n", " p_low = p_avg - xs\n", " p_high = p_avg + xs\n", "\n", " likelihood = (binom.pmf(k, n, p=p_low) + \n", " binom.pmf(k, n, p=p_high)) /2\n", " likelihood = np.nan_to_num(likelihood)\n", " \n", " return likelihood" ] }, { "cell_type": "code", "execution_count": 4, "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXoAAAEWCAYAAABollyxAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjMsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+AADFEAAAgAElEQVR4nO3de3Sc9X3n8fdX15FsSZZlWbIl+cLN4DgJoeYS0uakkLSQNHHaJiembUJ26aE5KWnTbbcl25Yl2ZzTptuzbC9sT0lISmkKpDTbOAnNpYFucwGCAZPYgMFg2bpYlizrZt0v3/1jnoHJIFlje0bPZT6vE4WZ53lm5js85qOff8/v+f3M3RERkeQqC7sAEREpLgW9iEjCKehFRBJOQS8iknAKehGRhFPQi4gknIJeYsvMfsbMDmY97zSzt5/F+9xuZv8QPN5kZqfMrDx4/u9m9uuFq3rJGj5sZt8r9udIaVLQS+QtFeDu/l1331bIz3L3o+6+2t3nC/m+58rM7jKzg2a2YGYfDrseiRcFvUg8PAN8FHgq7EIkfhT0Eltm9jYz615i38VmdtjMdgfPN5rZP5vZQLD9t5Z43RYzczOryNq82cy+b2ZjZvYtM1uXdfx7zOyAmQ0H3TyXZO27JNg2HBzznqx9TWa2x8xGzeyHwPmn+67ufqe7fweYyutfjkgWBb0kjpldBnwL+Ji7329mZcBXSbeK24BrgY+b2c/n+Za/AvwnYD1QBfxe8DkXAfcBHweagYeAr5pZlZlVBp/5reB1HwO+aGaZrqY7SYf2BuA/Bz8iRaGgl6T5GWAPcKO7fy3YdjnQ7O6fcvcZd38Z+CywO8/3/IK7v+Duk8CXgEuD7R8Avu7u33b3WeDPgRrgauAqYDXwp8FnPgx8DbghuND7y8Bt7j7u7vuBe871i4sspWL5Q0Ri5SPA/3P3R7K2bQY2mtlw1rZy4Lt5vmdf1uMJ0gEOsBE4ktnh7gtm1kX6bw1zQJe7L2S99kiwr5n0f3tdOftEikItekmajwCbzOyOrG1dwGF3X5P1U+fu7zzHz+ol/UsEADMzoAPoCfZ1BN1GGZuCfQOkfxF05OwTKQoFvcRFpZmlsn6W+tvoGHAd8FYz+9Ng2w+BUTP7AzOrMbNyM9thZpefY01fAt5lZtcGffK/C0wDPwAeB8aB3zezSjN7G/Bu4P5g6OaXgdvNrNbMtgM3nu6Dgn7/FGC8+u9C//1KXvQHReLiIWAy6+f2pQ5092HgHcD1ZvY/gmB9N+m+9cPACeBzQMO5FOTuB4FfA/4qeM93A+8O+uRngPcA1wf7/g/wIXd/Pnj5LaS7gPqAvwO+sMzHfYv0974auCt4/NZzqV9Kh2nhERGRZFOLXkQk4RT0IiIJp6AXEUk4Bb2ISMJF7oapdevW+ZYtW8IuQ0QkVp588skT7t682L7IBf2WLVvYu3dv2GWIiMSKmS15d7W6bkREEk5BLyKScAp6EZGEU9CLiCScgl5EJOEU9CIiCaegFxFJOAV9iVpYcP7l6R5GJmbDLkVEikxBX6Ie2NvFxx/Yx2179oddiogUmYK+BA1PzPBn33ieVGUZX9nXy5NHToZdkogUkYK+BP35tw4yOjXHF3/9KtbXVfPJrz7LwoIWoBFJKgV9idnfM8IXHz/KB6/azE9tbuTW6y/mR90jfPnpnrBLE5EiUdCXkIUF57av7KdpVRW/846LAHjvpW28sWMNf/aN5xmfngu5QhEpBgV9CflxzwhPHR3md95xEQ01lQCUlRn//d3b6R+b5p+f6g65QhEpBgV9CXn66BAA11y8/ie2X7apkea6avYdHQ6jLBEpMgV9CXmme4T1ddW01qdes+/1bQ3s7x0JoSoRKTYFfQnZ1zXMpR1rMLPX7NuxsZ5D/aeYmFE/vUjSKOhLxPDEDIdPjPPGjjWL7t/R1sCCw3PHxla4MhEpNgV9iXimO90t86bTBD2kh1+KSLIo6EvEvqPDmMHr2xsW3b+hIUXTqioFvUgCKehLxL6uIS5oXk1dqnLR/WbG69oa2N87usKViUixKehLgLvzTPcIly7RbZPx+rZ6Xjw+xtTs/ApVJiIrQUFfArpOTnJyfIZLNy0X9A3MLTgH+3RBViRJFPQlYF93+kaoN7afPuhftzHdf/9j9dOLJIqCvgTsOzpMqrKMba11pz2uvbGGNbWVHNCNUyKJoqAvAfu6htixsYHK8tOfbjNjx8YGtehFEkZBn3Cz8wvs7x1d9kJsxo62Bg72jTEzt1DkykRkpSjoE+6F4+nQXuqO2Fw72uqZnXdeOK4LsiJJoaBPuJcGxgG4sGV1Xse/vk0XZEWSRkGfcJ0n0kG/ee2qvI7vaKyluqKMlwdOFbMsEVlBeQW9mV1nZgfN7JCZ3brI/mozeyDY/7iZbcnZv8nMTpnZ7xWmbMlX54lxNjSkqKkqz+v4sjKjvbGGrpOTRa5MRFbKskFvZuXAncD1wHbgBjPbnnPYTcCQu18A3AF8Jmf/HcC/nnu5cqYOD46zpSm/1nxGx9pauoYmilSRiKy0fFr0VwCH3P1ld58B7gd25RyzC7gnePwgcK0Fk56b2XuBl4EDhSlZzkTniXG2rDuzoG9vrKF7SC16kaTIJ+jbgK6s593BtkWPcfc5YARoMrNVwB8Anzz3UuVMDU/MMDQxy3lnGPQdjbWMTM4yOjVbpMpEZCXlE/SvXY4IPM9jPgnc4e6nvbJnZjeb2V4z2zswMJBHSZKPw8GF2DNt0XesrQWg66S6b0SSIJ+g7wY6sp63A71LHWNmFUADcBK4EvgzM+sEPg78NzO7JfcD3P0ud9/p7jubm5vP+EvI4joH00G/dV3tGb2uozET9Oq+EUmCijyOeQK40My2Aj3AbuBXco7ZA9wIPAq8D3jY3R34mcwBZnY7cMrd/7oAdUseDp+YoMxebaHnq2NtDQDduiArkgjLBr27zwWt8G8C5cDn3f2AmX0K2Ovue4C7gXvN7BDplvzuYhYt+ek8Mc7GNTVUV+Q3tDKjoaaSuuoKdd2IJEQ+LXrc/SHgoZxtt2U9ngLev8x73H4W9ck56BwcZ+sZ9s9DenKz9rW1dGnkjUgi6M7YhHJ3Dp848zH0GR2NNWrRiySEgj6hTo7PMDY1d8YjbjI61tbSPTRJ+lKLiMSZgj6hMkMrz3TETUZHYw2Ts/OcODVTyLJEJAQK+oR6Nejzm7Uy1ytj6TXyRiT2FPQJ1Tk4TnkwQdnZ0E1TIsmhoE+ozhMTdDTWLLt84FIyvyA0541I/CnoE+rwWUxmlq22qoKmVVVq0YskgII+gdydzrOYnjhXu6YrFkkEBX0CDYxNMzEzf1Y3S2Xr0AIkIomgoE+gs521MlfH2lp6hyeZX9BYepE4U9AnUOYCasdZjrjJ6GisZW7BOTaiVr1InCnoE6hnOB3MG9ecY9AHs1iq+0Yk3hT0CdQ9NEFzXTWpyjObtTLXK/PS64KsSKwp6BOoZ3iStnNszUP6bwRmGksvEncK+gTqGZo86ztis1VVlLFudTV96qMXiTUFfcIsLDi9w1O0FSDoAVrrU/SNThfkvUQkHAr6hBk4Nc3M/ALtBei6AWhtSHF8ZKog7yUi4VDQJ0ymP7298eymJ86VbtEr6EXiTEGfMJkFvQvWddOQYmRylsmZ+YK8n4isPAV9wmTG0Bdi1A1AS30KQK16kRhT0CdM99AkjbWVrKrOa933ZbVmgl799CKxpaBPmJ6hyYJ12wC0NlQDcFwtepHYUtAnTM/wJO1rCnMhFtR1I5IECvoEcXe6hyYK2qKvS1WyqqpcXTciMaagT5CT4zNMzS4U5K7YbK0NKXXdiMSYgj5BCj3iJqO1QWPpReJMQZ8gmZulCtl1A+l+et0dKxJfCvoE6SnwXbEZrfUp+semWdBKUyKxpKBPkO6hCeqqK2ioqSzo+7Y2pJhbcE6Ma3IzkThS0CdIz3Bhx9BntOimKZFYU9AnSHeB5qHPpbtjReJNQZ8gPUOFWVkqV2tDOug1xFIknhT0CTEyOcvY9FzBL8QCrFtdTXmZaYilSEwp6BOip0hDKwHKy4z1ddX0jehirEgcKegTolg3S2W01OvuWJG4yivozew6MztoZofM7NZF9leb2QPB/sfNbEuw/Qoz2xf8PGNmv1jY8iWjp8ALjuTSSlMi8bVs0JtZOXAncD2wHbjBzLbnHHYTMOTuFwB3AJ8Jtu8Hdrr7pcB1wN+aWWEmSpef0DsyRVVFGU2rqory/lo7ViS+8mnRXwEccveX3X0GuB/YlXPMLuCe4PGDwLVmZu4+4e5zwfYUoFsri6RnOD3ixsyK8v4t9SnGpucYn55b/mARiZR8gr4N6Mp63h1sW/SYINhHgCYAM7vSzA4APwY+khX8rzCzm81sr5ntHRgYOPNvIfQOF2doZUZmARJ134jETz5Bv1gTMbdlvuQx7v64u78OuBz4hJmlXnOg+13uvtPddzY3N+dRkuTqGZpk45rX/KstGN0dKxJf+QR9N9CR9bwd6F3qmKAPvgE4mX2Auz8HjAM7zrZYWdz03Dz9Y9NsLGKLfkND+r0V9CLxk0/QPwFcaGZbzawK2A3syTlmD3Bj8Ph9wMPu7sFrKgDMbDOwDegsSOXyiuPB+PZiBn2rlhQUia1lR8C4+5yZ3QJ8EygHPu/uB8zsU8Bed98D3A3ca2aHSLfkdwcv/2ngVjObBRaAj7r7iWJ8kVKWGUPfXsSgr6kqpz5VobH0IjGU11BHd38IeChn221Zj6eA9y/yunuBe8+xRllGJuiL2aKHdD99/6jujhWJG90ZmwC9QdBnJh8rlpb6FMfH1KIXiRsFfQL0Dk+ybnU1qcryon7O+vpqtehFYkhBnwDFWnAkV0t9iv6xKS0pKBIzCvoESN8sVdxuG4CWumpm552hiZmif5aIFI6CPubcnZ7hSTY2FL9Fv74+swCJum9E4kRBH3NDE7NMzS4UfcQNQEt9ehoEXZAViRcFfcz1rtDQSoD1dekW/YBa9CKxoqCPuVdullqBi7HrMy163TQlEisK+pjLLCG4Ei366opyGmsr1XUjEjMK+pjrHZ4kVVlGY23linxeeklBdd2IxImCPuZ6RybZWMQFR3Ktr0/Rr64bkVhR0Mdcz/BUURccydVSV60WvUjMKOhjrmeouCtL5WqpTzFwapp53R0rEhsK+hibmp3nxKniLjiSq6W+mvkFZ3BcrXqRuFDQx1hmtaeVDPrM3bGa3EwkPhT0MZa5WWolu27W12ksvUjcKOhjrDuEoM8sEt4/pha9SFwo6GOse2iSMiv+giPZmtWiF4kdBX2MdQ9N0Fqfoqpi5U5jZXkZ61ZXaYilSIwo6GOse2iS9sbaFf/c9XW6aUokThT0MdYzNLkik5nlaqmv1nw3IjGioI+p2fkFjo2EFfSa70YkThT0MdU3MsWCE07XTX2KE6emmZtfWPHPFpEzp6CPqa6hCWBl5qHP1VJfjTucOKW1Y0XiQEEfU91DmQVHVr5F31KXWTtW/fQicaCgj6kwxtBnaKUpkXhR0MdUGGPoM3R3rEi8KOhjqntokrYQ+ucBmlZVUWZoLL1ITCjoY6onpJulACrKy1i3upo+Bb1ILCjoY2hufoG+0alQRtxkbGhI0aex9CKxoKCPoWMjU8wveKhB39qQom9kMrTPF5H8KehjKMyhlRkbGmo4NqKuG5E4UNDHUHeIN0tltDakGJua49T0XGg1iEh+FPQx1D00iVm6VR2WDcH4/T616kUiL6+gN7PrzOygmR0ys1sX2V9tZg8E+x83sy3B9neY2ZNm9uPgn9cUtvzS1D00GdoY+ozWegW9SFwsmxRmVg7cCVwPbAduMLPtOYfdBAy5+wXAHcBngu0ngHe7++uBG4F7C1V4Kesemgi12wZevSP3mC7IikRePk3CK4BD7v6yu88A9wO7co7ZBdwTPH4QuNbMzN2fdvfeYPsBIGVm1YUovJSFteBItha16EViI5+gbwO6sp53B9sWPcbd54ARoCnnmF8Gnnb31wy+NrObzWyvme0dGBjIt/aSFIUx9ACpynLWrqrimG6aEom8fILeFtnmZ3KMmb2OdHfObyz2Ae5+l7vvdPedzc3NeZRUuqIwhj6jtT7FcbXoRSIvn6DvBjqynrcDvUsdY2YVQANwMnjeDvxf4EPu/tK5FlzqojCGPmNDQ0pj6UViIJ+gfwK40My2mlkVsBvYk3PMHtIXWwHeBzzs7m5ma4CvA59w9+8XquhSlllwpG1NBFr0DSnNdyMSA8sGfdDnfgvwTeA54EvufsDMPmVm7wkOuxtoMrNDwH8BMkMwbwEuAP7YzPYFP+sL/i1KyJHBccrLLLSZK7NtaEhxcnyGqdn5sEsRkdOoyOcgd38IeChn221Zj6eA9y/yuk8Dnz7HGiVL52B6aGVlefj3urUGN2wdH51ic9OqkKsRkaWEnxZyRjpPjLMlIqG64ZWx9Oq+EYkyBX2MuDtHBifY0hT+hVh49aYpjaUXiTYFfYwMjs9wanqOLeui0aLPTIOgFr1ItCnoY6TzxDhAZLpuVlVXUJeq0Lz0IhGnoI+RzsH00MrNEem6gcxKU2rRi0SZgj5GOk+kh1ZG4WapjNaGGvXRi0Scgj5GOgfHaVtTE+r0xLk21OvuWJGoi05iyLKODE5E5kJsRmtDioFT08zOL4RdiogsQUEfE+4ejKGPTrcNpPvo3aF/7DWTkopIRCjoY+Lk+Axj03ORuwP11bH0GnkjElUK+pjoHEwPrdy6Lmot+vQ0COqnF4kuBX1MdJ7IDK2MaoteQS8SVQr6mDgyOE6ZQUeEhlYC1KcqqKksV4teJMIU9DFxeHCCtsZoDa0EMDM2rElpkXCRCItWasiSjgxGZ9bKXB2NtXSdVNCLRJWCPgbcncMRmp44V8faGo6enAi7DBFZgoI+BoYmZhmbmovUHDfZNq2tZWRyltGp2bBLEZFFKOhj4NWhlRFt0QcXiLvUqheJJAV9DGSmJ47a0MqMjrUKepEoU9DHwKH+U1SUGZvWRrPr5tWg1wVZkShS0MfAC8fHOK95VeSGVmY01FRSn6rQBVmRiIpmcshPOHh8jIta6sIu47Q2NdXSNaSgF4kiBX3EjU/P0XVykm0RD/r0WHoFvUgUKegj7sX+UwBc1BrxoF9bS9fQJAsLHnYpIpJDQR9xB/tGAbg4BkE/M7fAwCnNSy8SNQr6iDvYd4pUZVnkJjPL1dGYnq5YF2RFokdBH3EvBBdiy8os7FJOa5PG0otEloI+4uIw4gagrbEGM42lF4kiBX2EnRyfYWBsOvIjbgCqK8pprU+p60YkghT0EfbC8TEAtkX8QmxGR6PG0otEkYI+wuIW9O1ra9RHLxJBCvoIO9g3RkNNJevrqsMuJS+b1tbSNzrF9Nx82KWISBYFfYQd7BtjW0sdZtEecZPR0ViLO/QOa/1YkShR0EeUu6dH3LSuDruUvG0KFkbRBVmRaMkr6M3sOjM7aGaHzOzWRfZXm9kDwf7HzWxLsL3JzB4xs1Nm9teFLT3Z+kanGJuaY1trfdil5E0LkIhE07JBb2blwJ3A9cB24AYz255z2E3AkLtfANwBfCbYPgX8MfB7Bau4RBzsCy7ExmBoZcb6umqqKsoU9CIRk0+L/grgkLu/7O4zwP3ArpxjdgH3BI8fBK41M3P3cXf/HunAlzPwfBD0F7XEp+umrMxob9RC4SJRk0/QtwFdWc+7g22LHuPuc8AI0JRvEWZ2s5ntNbO9AwMD+b4s0fYdHWZzUy1raqvCLuWMXNC8+pVhoSISDfkE/WJDPnLnos3nmCW5+13uvtPddzY3N+f7ssRyd546OsRlmxrDLuWMbWuto3NwgqlZDbEUiYp8gr4b6Mh63g70LnWMmVUADcDJQhRYinqGJ+kfm+ayTWvCLuWMbWutY37BeWngVNiliEggn6B/ArjQzLaaWRWwG9iTc8we4Mbg8fuAh91dK1CcpaeODgPwpji26IOLx+q+EYmOiuUOcPc5M7sF+CZQDnze3Q+Y2aeAve6+B7gbuNfMDpFuye/OvN7MOoF6oMrM3gv8nLs/W/ivkhxPHx2iprI88ouNLGbLulVUlZe9cjFZRMK3bNADuPtDwEM5227LejwFvH+J1245h/pK0lNHh3lDewMV5fG7n62yvIzzmlfxgoJeJDLilyQJNzU7z7O9I7Hstsm4uLXulfsARCR8CvqI2d8zwuy8x/JCbMZFrXX0jkwxOjUbdikigoI+cp46OgTE80JsRubagrpvRKJBQR8xTx8dpmNtDc0xmZp4MZmlDw9q5I1IJCjoIyTON0pla1tTw+rqCvXTi0SEgj5CekemOD46HfugNzMualmtoBeJCAV9hDx1JN0/H/egB9jWWs/B42PovjmR8CnoI+Txw4PUVpVz8Yb43SiVa1vLaoYnZhkYmw67FJGSp6CPCHfn4ef6+ekL1lEZwxulcl0UjLzRHbIi4Yt/oiTEs8dG6R2Z4u2XtIRdSkFozhuR6FDQR8TDz/UD8LMXrw+5ksJoWl3NutXVatGLRICCPiL+7fl+3tixJtbj53O9bmM9z3QNh12GSMlT0EdA/9gUz3QN8/aEtOYzrjqviRf7T+mCrEjIFPQR8Mjz6W6baxPSP59x9fnp1SQffXkw5EpESpuCPgL+7bl+NjakuCQBwyqzvW5jPXWpCh596UTYpYiUNAV9yKZm5/neiye49pIWzBZbeje+KsrLuHJrEz94SS16kTAp6EP26MuDTM7Oc80lyeqfz7j6/CaODE7QMzwZdikiJUtBH7Kv/+gYtVXlvPm8prBLKYqrLwj66dWqFwmNgj5EJ8dn2PNML7/4pjZSleVhl1MUF62vY+2qKn6gfnqR0CjoQ3TfD48yM7fAh6/eEnYpRVNWZrz5vCYefWlQE5yJhERBH5K5+QW++NgR3nJBExe2JGu0Ta43n9/EsZEpOgcnwi5FpCQp6EPy7WeP0zsyxY1v3hJ2KUX3ynh69dOLhEJBH5K/+0En7Y01ibtJajFb162itT7F99VPLxIKBX0Injs2yuOHT/LBqzZTXpassfOLMTOuuWQ933nuOMMTM2GXI1JyFPQhuPORQ6Qqy/jA5R1hl7JiPnjVZqZmF3jgia6wSxEpOQr6FfbvB/v52o+OcfNbz2dNbVXY5ayYSzbUc8XWtdz72BHmFzT6RmQlKehX0MTMHH/0L/s5v3kVv/mz54ddzor78NVb6B6a5DvPHQ+7FJGSoqBfQXd8+wW6hyb5k196A9UVybxB6nR+bnsLGxpS3PNoZ9iliJQUBf0K2d8zwt3fO8wNV2ziiq1rwy4nFBXlZfzaVZv5/qFBXtQSgyIrRkG/ArpOTvAb9z5J0+pqbr3+4rDLCdXuyzuoqijjCz/oDLsUkZKhoC+ynuFJbvjsY4xNzfKFD19OQ01l2CWFqml1Ne/7qXbu/+FR/uOFgbDLESkJCvoiOjYyyQ13PcbI5Cz/8OtXsqOtIeySIuGP3nUJF7XU8bH7nqbrpKZFECk2BX0RuDv/tLeLd/7Fdxkan+Hem67kDe1rwi4rMmqrKvjbD/4U7s7N9z7J5Mx82CWJJJqCvsCe7R3lhs8+xn998Eec37yaL3/0ai7tUMjn2ty0ir+44U083zfKx+57mpPjumNWpFjyCnozu87MDprZITO7dZH91Wb2QLD/cTPbkrXvE8H2g2b284UrPTr6Rqb4wvcP8wt/9V3e+Zff5dneUf7kl17Pl37jzYmfmfJc/Oy29fzxu7bz7wf7edv/fIS/f7STufmFsMsSSRxbbo5wMysHXgDeAXQDTwA3uPuzWcd8FHiDu3/EzHYDv+juHzCz7cB9wBXARuDfgIvcfcm/q+/cudP37t17jl+rcNyd8Zl5RidnGZ2apX90mmMjk/QMT3Gwb5RnukboG50CYEdbPb98WTvvvbSNxlWlc9fruXrx+Bi3f/UA3z80SHtjDW85fx1XbF3LG9obWLe6moaaSspKYE4gkXNhZk+6+87F9lXk8forgEPu/nLwZvcDu4Bns47ZBdwePH4Q+GtLr3S9C7jf3aeBw2Z2KHi/R8/mi5zO832j3PKPTy+6L/uXmQf/58F2BxbccYeFBWcu+JmdW2B6boGZJVqYZrB5bS1XnreWN7av4S0XrGNbq1rvZ+PCljr+4aYr+eaBPh58spt/3X+MB/a+OidOeZlRn6qgqqKMyvIyKsqMMjNI/2/JRdX1q0Hi5m3bmvnDd20v+PvmE/RtQPZMVN3AlUsd4+5zZjYCNAXbH8t5bVvuB5jZzcDNAJs2bcq39p+Qqihn2+m6SewnH5oZBpTZq4/Ly4yK8nSIVFeUU11ZRlV5GaurK6hLVVCXqmTd6io2rqmhpT5FVYUucRSKmXHdjg1ct2MDCwvOC/1jHOwbY/DUDCfHZxiZnGV2Pv2Ld27+1V/QLPEXUl9qh0iEtdSnivK++QT9Yg2j3P+Kljomn9fi7ncBd0G66yaPml5jy7pV3Pmrl53NSyViysqMi1vrubi1PuxSRBIhnyZpN5A9n2470LvUMWZWATQAJ/N8rYiIFFE+Qf8EcKGZbTWzKmA3sCfnmD3AjcHj9wEPe7pjfA+wOxiVsxW4EPhhYUoXEZF8LNt1E/S53wJ8EygHPu/uB8zsU8Bed98D3A3cG1xsPUn6lwHBcV8ifeF2DvjN0424ERGRwlt2eOVKi9rwShGRODjd8EoNGxERSTgFvYhIwinoRUQSTkEvIpJwkbsYa2YDwJFzeIt1wIkClRMHpfZ9Qd+5VOg7n5nN7t682I7IBf25MrO9S115TqJS+76g71wq9J0LR103IiIJp6AXEUm4JAb9XWEXsMJK7fuCvnOp0HcukMT10YuIyE9KYoteRESyKOhFRBIuMUG/3ALmSWBmHWb2iJk9Z2YHzOy3g+1rzezbZvZi8M/GsGstJDMrN7OnzexrwfOtwSL0LwaL0idugV4zW2NmD5rZ88H5fnOSz7OZ/U7wZ3q/md1nZqkknmcz+7yZ9ZvZ/qxti55XS/vLINN+ZGZnvbJSIoI+WMD8TuB6YDtwQ7AwedLMAb/r7pcAVwG/GXzPW4HvuPuFwHeC50ny28BzWc8/A9wRfN8h4KZQqiquvwE63QoAAANZSURBVAC+4e4XA28k/f0TeZ7NrA34LWCnu+8gPR36bpJ5nv8OuC5n21Ln9XrSa3hcSHqp1b852w9NRNCTtYC5u88AmQXME8Xdj7n7U8HjMdL/8beR/q73BIfdA7w3nAoLz8zagXcBnwueG3AN6UXoIWHfF8DM6oG3kl7nAXefcfdhEnyeSa+NUROsUFcLHCOB59nd/4P0mh3Zljqvu4C/97THgDVmtuFsPjcpQb/YAuavWYQ8ScxsC/Am4HGgxd2PQfqXAbA+vMoK7n8Dvw8sBM+bgGF3nwueJ/FcnwcMAF8Iuqw+Z2arSOh5dvce4M+Bo6QDfgR4kuSf54ylzmvBci0pQZ/XIuRJYWargX8GPu7uo2HXUyxm9gtAv7s/mb15kUOTdq4rgMuAv3H3NwHjJKSbZjFBn/QuYCuwEVhFutsiV9LO83IK9mc9KUFfMouQm1kl6ZD/ort/Odh8PPNXuuCf/WHVV2BvAd5jZp2ku+OuId3CXxP8FR+Sea67gW53fzx4/iDp4E/qeX47cNjdB9x9FvgycDXJP88ZS53XguVaUoI+nwXMYy/on74beM7d/1fWruzF2W8EvrLStRWDu3/C3dvdfQvpc/qwu/8q8AjpReghQd83w937gC4z2xZsupb0usuJPM+ku2yuMrPa4M945vsm+jxnWeq87gE+FIy+uQoYyXTxnDF3T8QP8E7gBeAl4A/DrqdI3/GnSf/V7UfAvuDnnaT7rb8DvBj8c23YtRbhu78N+Frw+Dzgh8Ah4J+A6rDrK8L3vRTYG5zrfwEak3yegU8CzwP7gXuB6iSeZ+A+0tchZkm32G9a6ryS7rq5M8i0H5MelXRWn6spEEREEi4pXTciIrIEBb2ISMIp6EVEEk5BLyKScAp6EZGEU9CLiCScgl5EJOEU9CLLMLPLg/nAU2a2Kpg3fUfYdYnkSzdMieTBzD4NpIAa0vPQ/EnIJYnkTUEvkodgDqUngCnganefD7kkkbyp60YkP2uB1UAd6Za9SGyoRS+SBzPbQ3qq5K3ABne/JeSSRPJWsfwhIqXNzD4EzLn7PwbrE//AzK5x94fDrk0kH2rRi4gknProRUQSTkEvIpJwCnoRkYRT0IuIJJyCXkQk4RT0IiIJp6AXEUm4/w+FIHO42za+DwAAAABJRU5ErkJggg==\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "xs = uniform.index\n", "p_avg = 0.5\n", "k = 75\n", "n = 100\n", "likelihood1 = compute_likelihood(xs, p_avg, k, n)\n", "\n", "plt.plot(likelihood1)\n", "plt.xlabel('x')\n", "plt.title('Likelihood 1');" ] }, { "cell_type": "code", "execution_count": 5, "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXoAAAEWCAYAAABollyxAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjMsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+AADFEAAAgAElEQVR4nO3de3Rc51nv8e+j60i2JMuyLNmSfMnFTly3TYNzaQpdJWkhaWndQrtwgDblhBW6IIVyOZAeICft4Swoh3VygAYWKWkJoSQpoYe6JZCWJhzaNDcncVo7iRMnlq2LZcu6W/fLc/6YPcl0Itlje0b7Mr/PQmRm7z0zz3QnP71697vf19wdERFJrrKwCxARkeJS0IuIJJyCXkQk4RT0IiIJp6AXEUk4Bb2ISMIp6CW2zOzHzOxA1vNOM3v3WbzPbWb298HjDWZ20szKg+f/YWa/VLiql6zh42b23WJ/jpQmBb1E3lIB7u7fcfethfwsdz/i7ivdfb6Q73suzGyLmX3NzPrNbNDMHjKzgn5vSTYFvUj0rQJ2A1uBFuBJ4GuhViSxoqCX2DKzd5lZ9xL7LjKzQ2a2K3i+3sz+KWgVHzKzX1vidZvMzM2sImvzRjN71MzGzOybZrYm6/gPmNl+MxsOunkuztp3cbBtODjmA1n7msxst5mNmtmTwPlLfU93f9Ld73L3QXefBW4HtppZU77/W0lpU9BL4pjZpcA3gU+6+31mVgZ8HXgOaAOuAT5lZj+Z51v+HPCLwFqgCvjt4HO2APcCnwKagQeBr5tZlZlVBp/5zeB1nwS+nNXlcgcwBawD/kvwk693An3uPnAGr5ESpqCXpPkx0t0cN7j7N4JtlwHN7v5Zd59x91eBLwC78nzPL7n7S+4+CXwFuCTY/rPAv7j7t4KW9p8CNcBVwJXASuCPg898GPgGcH1wofdngFvdfdzd9wF351OImbWT/iXxm3nWLkLF6Q8RiZVPAP/P3R/J2rYRWG9mw1nbyoHv5PmefVmPJ0gHOMB64HBmh7svmFkX6b8a5oAud1/Ieu3hYF8z6f/2unL2nZKZNZP+C+Ev3f3ePGsXUYteEucTwAYzuz1rWxdwyN1XZf3Uuft7z/Gzekn/EgHAzAzoAHqCfR1Bt1HGhmBfP+lfBB05+5ZkZo2kQ363u//Pc6xbSoyCXuKi0sxSWT9L/TU6BlwLvNPM/jjY9iQwama/a2Y1ZlZuZtvN7LJzrOkrwPvM7JqgT/63gGnge8ATwDjwO2ZWaWbvAt4P3BcM3fwqcJuZ1ZrZNuCGpT7EzOqBh4BH3f2Wc6xZSpCCXuLiQWAy6+e2pQ5092HgPcB1ZvY/gmB9P+m+9UPACeBvgIZzKcjdDwC/APxF8J7vB94f9MnPAB8Argv2/SXwMXd/MXj5zaS7gPqAvwW+dIqP+hDp6wy/GNzMlfk55V8BIhmmhUdERJJNLXoRkYRT0IuIJJyCXkQk4RT0IiIJF7kbptasWeObNm0KuwwRkVh5+umnT7h782L7Ihf0mzZtYs+ePWGXISISK2a25N3V6roREUk4Bb2ISMIp6EVEEk5BLyKScAp6EZGEU9CLiCScgl5EJOEU9CVqYcH552d7GJmYDbsUESkyBX2Jun9PF5+6fy+37t4XdikiUmQK+hI0PDHDn/zbi6Qqy/ja3l6ePjwYdkkiUkQK+hL0p988wOjUHF/+pStZW1fNZ77+PAsLWoBGJKkU9CVmX88IX37iCB+9ciM/srGRW667iO93j/DVZ3vCLk1EikRBX0IWFpxbv7aPphVV/MZ7tgDwwUvaeGvHKv7k315kfHou5ApFpBgU9CXkBz0jPHNkmN94zxYaaioBKCsz/vv7t3F8bJp/eqY75ApFpBgU9CXk2SNDAFx90dof2n7phkaa66rZe2Q4jLJEpMgU9CXkue4R1tZV01qfesO+N7c1sK93JISqRKTYFPQlZG/XMJd0rMLM3rBv+/p6Dh4/ycSM+ulFkkZBXyKGJ2Y4dGKct3asWnT/9rYGFhxeODq2zJWJSLEp6EvEc93pbpm3nSLoIT38UkSSRUFfIvYeGcYM3tzesOj+dQ0pmlZUKehFEkhBXyL2dg1xQfNK6lKVi+43M97U1sC+3tFlrkxEik1BXwLcnee6R7hkiW6bjDe31fPysTGmZueXqTIRWQ4K+hLQNTjJ4PgMl2w4XdA3MLfgHOjTBVmRJFHQl4C93ekbod7afuqgf9P6dP/9D9RPL5IoCvoSsPfIMKnKMra21p3yuPbGGlbVVrJfN06JJIqCvgTs7Rpi+/oGKstPfbrNjO3rG9SiF0kYBX3Czc4vsK939LQXYjO2tzVwoG+MmbmFIlcmIstFQZ9wLx1Lh/ZSd8Tm2t5Wz+y889IxXZAVSQoFfcK90j8OwIUtK/M6/s1tuiArkjQK+oTrPJEO+o2rV+R1fEdjLdUVZbzaf7KYZYnIMsor6M3sWjM7YGYHzeyWRfZXm9n9wf4nzGxTzv4NZnbSzH67MGVLvjpPjLOuIUVNVXlex5eVGe2NNXQNTha5MhFZLqcNejMrB+4ArgO2Adeb2bacw24Ehtz9AuB24HM5+28H/vXcy5UzdWhgnE1N+bXmMzpW19I1NFGkikRkueXTor8cOOjur7r7DHAfsDPnmJ3A3cHjB4BrLJj03Mw+CLwK7C9MyXImOk+Ms2nNmQV9e2MN3UNq0YskRT5B3wZ0ZT3vDrYteoy7zwEjQJOZrQB+F/jMuZcqZ2p4YoahiVnOO8Og72isZWRyltGp2SJVJiLLKZ+gf+NyROB5HvMZ4HZ3P+WVPTO7ycz2mNme/v7+PEqSfBwKLsSeaYu+Y3UtAF2D6r4RSYJ8gr4b6Mh63g70LnWMmVUADcAgcAXwJ2bWCXwK+G9mdnPuB7j7ne6+w913NDc3n/GXkMV1DqSDfvOa2jN6XUdjJujVfSOSBBV5HPMUcKGZbQZ6gF3Az+Ucsxu4AXgM+DDwsLs78GOZA8zsNuCku3++AHVLHg6dmKDMXm+h56tjdQ0A3bogK5IIpw16d58LWuEPAeXAF919v5l9Ftjj7ruBu4B7zOwg6Zb8rmIWLfnpPDHO+lU1VFfkN7Qyo6GmkrrqCnXdiCREPi163P1B4MGcbbdmPZ4CPnKa97jtLOqTc9A5MM7mM+yfh/TkZu2ra+nSyBuRRNCdsQnl7hw6ceZj6DM6GmvUohdJCAV9Qg2OzzA2NXfGI24yOlbX0j00SfpSi4jEmYI+oTJDK890xE1GR2MNk7PznDg5U8iyRCQECvqEej3o85u1MtdrY+k18kYk9hT0CdU5ME55MEHZ2dBNUyLJoaBPqM4TE3Q01px2+cClZH5BaM4bkfhT0CfUobOYzCxbbVUFTSuq1KIXSQAFfQK5O51nMT1xrnZNVyySCAr6BOofm2ZiZv6sbpbK1qEFSEQSQUGfQGc7a2WujtW19A5PMr+gsfQicaagT6DMBdSOsxxxk9HRWMvcgnN0RK16kThT0CdQz3A6mNevOsegD2axVPeNSLwp6BOoe2iC5rpqUpVnNmtlrtfmpdcFWZFYU9AnUM/wJG3n2JqH9F8EZhpLLxJ3CvoE6hmaPOs7YrNVVZSxZmU1feqjF4k1BX3CLCw4vcNTtBUg6AFa61P0jU4X5L1EJBwK+oTpPznNzPwC7QXougFobUhxbGSqIO8lIuFQ0CdMpj+9vfHspifOlW7RK+hF4kxBnzCZBb0L1nXTkGJkcpbJmfmCvJ+ILD8FfcJkxtAXYtQNQEt9CkCtepEYU9AnTPfQJI21layozmvd99NqzQS9+ulFYktBnzA9Q5MF67YBaG2oBuCYWvQisaWgT5ie4UnaVxXmQiyo60YkCRT0CeLudA9NFLRFX5eqZEVVubpuRGJMQZ8gg+MzTM0uFOSu2GytDSl13YjEmII+QQo94iajtUFj6UXiTEGfIJmbpQrZdQPpfnrdHSsSXwr6BOkp8F2xGa31KY6PTbOglaZEYklBnyDdQxPUVVfQUFNZ0PdtbUgxt+CcGNfkZiJxpKBPkJ7hwo6hz2jRTVMisaagT5DuAs1Dn0t3x4rEm4I+QXqGCrOyVK7WhnTQa4ilSDwp6BNiZHKWsem5gl+IBVizspryMtMQS5GYUtAnRE+RhlYClJcZa+uq6RvRxViROFLQJ0SxbpbKaKnX3bEicZVX0JvZtWZ2wMwOmtkti+yvNrP7g/1PmNmmYPvlZrY3+HnOzD5U2PIlo6fAC47k0kpTIvF12qA3s3LgDuA6YBtwvZltyznsRmDI3S8Abgc+F2zfB+xw90uAa4G/NrPCTJQuP6R3ZIqqijKaVlQV5f21dqxIfOXTor8cOOjur7r7DHAfsDPnmJ3A3cHjB4BrzMzcfcLd54LtKUC3VhZJz3B6xI2ZFeX9W+pTjE3PMT49d/qDRSRS8gn6NqAr63l3sG3RY4JgHwGaAMzsCjPbD/wA+ERW8L/GzG4ysz1mtqe/v//Mv4XQO1ycoZUZmQVI1H0jEj/5BP1iTcTclvmSx7j7E+7+JuAy4NNmlnrDge53uvsOd9/R3NycR0mSq2dokvWr3vA/bcHo7liR+Mon6LuBjqzn7UDvUscEffANwGD2Ae7+AjAObD/bYmVx03PzHB+bZn0RW/TrGtLvraAXiZ98gv4p4EIz22xmVcAuYHfOMbuBG4LHHwYedncPXlMBYGYbga1AZ0Eql9ccC8a3FzPoW7WkoEhsnXYEjLvPmdnNwENAOfBFd99vZp8F9rj7buAu4B4zO0i6Jb8rePmPAreY2SywAPyKu58oxhcpZZkx9O1FDPqaqnLqUxUaSy8SQ3kNdXT3B4EHc7bdmvV4CvjIIq+7B7jnHGuU08gEfTFb9JDupz8+qrtjReJGd8YmQG8Q9JnJx4qlpT7FsTG16EXiRkGfAL3Dk6xZWU2qsryon7O2vlotepEYUtAnQLEWHMnVUp/i+NiUlhQUiRkFfQKkb5YqbrcNQEtdNbPzztDETNE/S0QKR0Efc+5Oz/Ak6xuK36JfW59ZgETdNyJxoqCPuaGJWaZmF4o+4gagpT49DYIuyIrEi4I+5nqXaWglwNq6dIu+Xy16kVhR0MfcazdLLcPF2LWZFr1umhKJFQV9zGWWEFyOFn11RTmNtZXquhGJGQV9zPUOT5KqLKOxtnJZPi+9pKC6bkTiREEfc70jk6wv4oIjudbWpziurhuRWFHQx1zP8FRRFxzJ1VJXrRa9SMwo6GOuZ6i4K0vlaqlP0X9ymnndHSsSGwr6GJuanefEyeIuOJKrpb6a+QVnYFytepG4UNDHWGa1p+UM+szdsZrcTCQ+FPQxlrlZajm7btbWaSy9SNwo6GOsO4SgzywSfnxMLXqRuFDQx1j30CRlVvwFR7I1q0UvEjsK+hjrHpqgtT5FVcXyncbK8jLWrKzSEEuRGFHQx1j30CTtjbXL/rlr63TTlEicKOhjrGdoclkmM8vVUl+t+W5EYkRBH1Oz8wscHQkr6DXfjUicKOhjqm9kigUnnK6b+hQnTk4zN7+w7J8tImdOQR9TXUMTwPLMQ5+rpb4adzhxUmvHisSBgj6muocyC44sf4u+pS6zdqz66UXiQEEfU2GMoc/QSlMi8aKgj6kwxtBn6O5YkXhR0MdU99AkbSH0zwM0raiizNBYepGYUNDHVE9IN0sBVJSXsWZlNX0KepFYUNDH0Nz8An2jU6GMuMlY15CiT2PpRWJBQR9DR0emmF/wUIO+tSFF38hkaJ8vIvlT0MdQmEMrM9Y11HB0RF03InGgoI+h7hBvlspobUgxNjXHyem50GoQkfwo6GOoe2gSs3SrOizrgvH7fWrVi0ReXkFvZtea2QEzO2hmtyyyv9rM7g/2P2Fmm4Lt7zGzp83sB8E/ry5s+aWpe2gytDH0Ga31CnqRuDhtUphZOXAHcB2wDbjezLblHHYjMOTuFwC3A58Ltp8A3u/ubwZuAO4pVOGlrHtoItRuG3j9jtyjuiArEnn5NAkvBw66+6vuPgPcB+zMOWYncHfw+AHgGjMzd3/W3XuD7fuBlJlVF6LwUhbWgiPZWtSiF4mNfIK+DejKet4dbFv0GHefA0aAppxjfgZ41t3fMPjazG4ysz1mtqe/vz/f2ktSFMbQA6Qqy1m9ooqjumlKJPLyCXpbZJufyTFm9ibS3Tm/vNgHuPud7r7D3Xc0NzfnUVLpisIY+ozW+hTH1KIXibx8gr4b6Mh63g70LnWMmVUADcBg8Lwd+L/Ax9z9lXMtuNRFYQx9xrqGlMbSi8RAPkH/FHChmW02sypgF7A755jdpC+2AnwYeNjd3cxWAf8CfNrdHy1U0aUss+BI26oItOgbUprvRiQGThv0QZ/7zcBDwAvAV9x9v5l91sw+EBx2F9BkZgeB3wQyQzBvBi4A/sDM9gY/awv+LUrI4YFxyssstJkrs61rSDE4PsPU7HzYpYjIKVTkc5C7Pwg8mLPt1qzHU8BHFnndHwJ/eI41SpbOgfTQysry8O91aw1u2Do2OsXGphUhVyMiSwk/LeSMdJ4YZ1NEQnXda2Pp1X0jEmUK+hhxdw4PTLCpKfwLsfD6TVMaSy8SbQr6GBkYn+Hk9Byb1kSjRZ+ZBkEtepFoU9DHSOeJcYDIdN2sqK6gLlWheelFIk5BHyOdA+mhlRsj0nUDmZWm1KIXiTIFfYx0nkgPrYzCzVIZrQ016qMXiTgFfYx0DozTtqom1OmJc62r192xIlEXncSQ0zo8MBGZC7EZrQ0p+k9OMzu/EHYpIrIEBX1MuHswhj463TaQ7qN3h+Njb5iUVEQiQkEfE4PjM4xNz0XuDtTXx9Jr5I1IVCnoY6JzID20cvOaqLXo09MgqJ9eJLoU9DHReSIztDKqLXoFvUhUKehj4vDAOGUGHREaWglQn6qgprJcLXqRCFPQx8ShgQnaGqM1tBLAzFi3KqVFwkUiLFqpIUs6PBCdWStzdTTW0jWooBeJKgV9DLg7hyI0PXGujtU1HBmcCLsMEVmCgj4GhiZmGZuai9QcN9k2rK5lZHKW0anZsEsRkUUo6GPg9aGVEW3RBxeIu9SqF4kkBX0MZKYnjtrQyoyO1Qp6kShT0MfAweMnqSgzNqyOZtfN60GvC7IiUaSgj4GXjo1xXvOKyA2tzGioqaQ+VaELsiIRFc3kkB9y4NgYW1rqwi7jlDY01dI1pKAXiSIFfcSNT8/RNTjJ1ogHfXosvYJeJIoU9BH38vGTAGxpjXjQr66la2iShQUPuxQRyaGgj7gDfaMAXBSDoJ+ZW6D/pOalF4kaBX3EHeg7SaqyLHKTmeXqaExPV6wLsiLRo6CPuJeCC7FlZRZ2Kae0QWPpRSJLQR9xcRhxA9DWWIOZxtKLRJGCPsIGx2foH5uO/IgbgOqKclrrU+q6EYkgBX2EvXRsDICtEb8Qm9HRqLH0IlGkoI+wuAV9++oa9dGLRJCCPsIO9I3RUFPJ2rrqsEvJy4bVtfSNTjE9Nx92KSKSRUEfYQf6xtjaUodZtEfcZHQ01uIOvcNaP1YkShT0EeXu6RE3rSvDLiVvG4KFUXRBViRa8gp6M7vWzA6Y2UEzu2WR/dVmdn+w/wkz2xRsbzKzR8zspJl9vrClJ1vf6BRjU3Nsba0Pu5S8aQESkWg6bdCbWTlwB3AdsA243sy25Rx2IzDk7hcAtwOfC7ZPAX8A/HbBKi4RB/qCC7ExGFqZsbaumqqKMgW9SMTk06K/HDjo7q+6+wxwH7Az55idwN3B4weAa8zM3H3c3b9LOvDlDLwYBP2Wlvh03ZSVGe2NWihcJGryCfo2oCvreXewbdFj3H0OGAGa8i3CzG4ysz1mtqe/vz/flyXa3iPDbGyqZVVtVdilnJELmle+NixURKIhn6BfbMhH7ly0+RyzJHe/0913uPuO5ubmfF+WWO7OM0eGuHRDY9ilnLGtrXV0DkwwNashliJRkU/QdwMdWc/bgd6ljjGzCqABGCxEgaWoZ3iS42PTXLphVdilnLGtrXXMLziv9J8MuxQRCeQT9E8BF5rZZjOrAnYBu3OO2Q3cEDz+MPCwu2sFirP0zJFhAN4WxxZ9cPFY3Tci0VFxugPcfc7MbgYeAsqBL7r7fjP7LLDH3XcDdwH3mNlB0i35XZnXm1knUA9UmdkHgZ9w9+cL/1WS49kjQ9RUlkd+sZHFbFqzgqrystcuJotI+E4b9ADu/iDwYM62W7MeTwEfWeK1m86hvpL0zJFh3tLeQEV5/O5nqywv47zmFbykoBeJjPglScJNzc7zfO9ILLttMi5qrXvtPgARCZ+CPmL29YwwO++xvBCbsaW1jt6RKUanZsMuRURQ0EfOM0eGgHheiM3IXFtQ941INCjoI+bZI8N0rK6hOSZTEy8ms/ThAY28EYkEBX2ExPlGqWxtq2pYWV2hfnqRiFDQR0jvyBTHRqdjH/RmxpaWlQp6kYhQ0EfIM4fT/fNxD3qAra31HDg2hu6bEwmfgj5Cnjg0QG1VOReti9+NUrm2tqxkeGKW/rHpsEsRKXkK+ohwdx5+4Tg/esEaKmN4o1SuLcHIG90hKxK++CdKQjx/dJTekSnefXFL2KUUhOa8EYkOBX1EPPzCcQB+/KK1IVdSGE0rq1mzslotepEIUNBHxL+/eJy3dqyK9fj5XG9aX89zXcNhlyFS8hT0EXB8bIrnuoZ5d0Ja8xlXntfEy8dP6oKsSMgU9BHwyIvpbptrEtI/n3HV+enVJB97dSDkSkRKm4I+Av79heOsb0hxcQKGVWZ70/p66lIVPPbKibBLESlpCvqQTc3O892XT3DNxS2YLbb0bnxVlJdxxeYmvveKWvQiYVLQh+yxVweYnJ3n6ouT1T+fcdX5TRwemKBneDLsUkRKloI+ZP/y/aPUVpXz9vOawi6lKK66IOinV6teJDQK+hANjs+w+7lePvS2NlKV5WGXUxRb1taxekUV31M/vUhoFPQhuvfJI8zMLfDxqzaFXUrRlJUZbz+vicdeGdAEZyIhUdCHZG5+gS8/fph3XNDEhS3JGm2T6+3nN3F0ZIrOgYmwSxEpSQr6kHzr+WP0jkxxw9s3hV1K0b02nl799CKhUNCH5G+/10l7Y03ibpJazOY1K2itT/Go+ulFQqGgD8ELR0d54tAgH71yI+VlyRo7vxgz4+qL1/LtF44xPDETdjkiJUdBH4I7HjlIqrKMn72sI+xSls1Hr9zI1OwC9z/VFXYpIiVHQb/M/uPAcb7x/aPc9M7zWVVbFXY5y+bidfVcvnk19zx+mPkFjb4RWU4K+mU0MTPH7//zPs5vXsGv/vj5YZez7D5+1Sa6hyb59gvHwi5FpKQo6JfR7d96ie6hSf7op99CdUUyb5A6lZ/Y1sK6hhR3P9YZdikiJUVBv0z29Yxw13cPcf3lG7h88+qwywlFRXkZv3DlRh49OMDLWmJQZNko6JdB1+AEv3zP0zStrOaW6y4Ku5xQ7bqsg6qKMr70vc6wSxEpGQr6IusZnuT6LzzO2NQsX/r4ZTTUVIZdUqiaVlbz4R9p574nj/CfL/WHXY5ISVDQF9HRkUmuv/NxRiZn+ftfuoLtbQ1hlxQJv/++i9nSUscn732WrkFNiyBSbAr6InB3/nFPF+/9s+8wND7DPTdewVvaV4VdVmTUVlXw1x/9Edydm+55msmZ+bBLEkk0BX2BPd87yvVfeJz/+sD3Ob95JV/9lau4pEMhn2tj0wr+7Pq38WLfKJ+891kGx3XHrEix5BX0ZnatmR0ws4Nmdssi+6vN7P5g/xNmtilr36eD7QfM7CcLV3p09I1M8aVHD/FTf/Ed3vvn3+H53lH+6KffzFd++e2Jn5nyXPz41rX8wfu28R8HjvOu//UIf/dYJ3PzC2GXJZI4dro5ws2sHHgJeA/QDTwFXO/uz2cd8yvAW9z9E2a2C/iQu/+smW0D7gUuB9YD/w5scfcl/1bfsWOH79mz5xy/VuG4O+Mz84xOzjI6Ncvx0WmOjkzSMzzFgb5RnusaoW90CoDtbfX8zKXtfPCSNhpXlM5dr+fq5WNj3Pb1/Tx6cID2xhrecf4aLt+8mre0N7BmZTUNNZWUlcCcQCLnwsyedvcdi+2ryOP1lwMH3f3V4M3uA3YCz2cdsxO4LXj8APB5S690vRO4z92ngUNmdjB4v8fO5oucyot9o9z8D88uui/7l5kH/8+D7Q4suOMOCwvOXPAzO7fA9NwCM0u0MM1g4+parjhvNW9tX8U7LljD1la13s/GhS11/P2NV/DQ/j4eeLqbf913lPv3vD4nTnmZUZ+qoKqijMryMirKjDIzSP/fkouq61eDxM27tjbze+/bVvD3zSfo24Dsmai6gSuWOsbd58xsBGgKtj+e89q23A8ws5uAmwA2bNiQb+0/JFVRztZTdZPYDz80Mwwos9cfl5cZFeXpEKmuKKe6soyq8jJWVldQl6qgLlXJmpVVrF9VQ0t9iqoKXeIoFDPj2u3ruHb7OhYWnJeOj3Ggb4yBkzMMjs8wMjnL7Hz6F+/c/Ou/oFniD1JfaodIhLXUp4ryvvkE/WINo9z/ipY6Jp/X4u53AndCuusmj5reYNOaFdzx85eezUslYsrKjIta67motT7sUkQSIZ8maTeQPZ9uO9C71DFmVgE0AIN5vlZERIoon6B/CrjQzDabWRWwC9idc8xu4Ibg8YeBhz3dMb4b2BWMytkMXAg8WZjSRUQkH6ftugn63G8GHgLKgS+6+34z+yywx913A3cB9wQXWwdJ/zIgOO4rpC/czgG/eqoRNyIiUninHV653KI2vFJEJA5ONbxSw0ZERBJOQS8iknAKehGRhFPQi4gkXOQuxppZP3D4HN5iDXCiQOXEQal9X9B3LhX6zmdmo7s3L7YjckF/rsxsz1JXnpOo1L4v6DuXCn3nwlHXjYhIwinoRUQSLolBf2fYBSyzUvu+oO9cKvSdCyRxffQiIvLDktiiFxGRLAp6EZGES0zQn24B8yQwsw4ze8TMXjCz/Wb268H21Wb2LTN7OfhnY9i1FpKZlZvZs2b2jeD55mAR+peDRekTt0Cvma0yswfM7MXgfL89yefZzH4j+Hd6n5nda2apJJ5nM/uimR03s31Z2xY9r5b250Gmfd/MznplpbECLWgAAAOXSURBVEQEfbCA+R3AdcA24PpgYfKkmQN+y90vBq4EfjX4nrcA33b3C4FvB8+T5NeBF7Kefw64Pfi+Q8CNoVRVXH8G/Ju7XwS8lfT3T+R5NrM24NeAHe6+nfR06LtI5nn+W+DanG1LndfrSK/hcSHppVb/6mw/NBFBT9YC5u4+A2QWME8Udz/q7s8Ej8dI/8ffRvq73h0cdjfwwXAqLDwzawfeB/xN8NyAq0kvQg8J+74AZlYPvJP0Og+4+4y7D5Pg80x6bYyaYIW6WuAoCTzP7v6fpNfsyLbUed0J/J2nPQ6sMrN1Z/O5SQn6xRYwf8Mi5EliZpuAtwFPAC3ufhTSvwyAteFVVnD/B/gdYCF43gQMu/tc8DyJ5/o8oB/4UtBl9TdmtoKEnmd37wH+FDhCOuBHgKdJ/nnOWOq8FizXkhL0eS1CnhRmthL4J+BT7j4adj3FYmY/BRx396ezNy9yaNLOdQVwKfBX7v42YJyEdNMsJuiT3glsBtYDK0h3W+RK2nk+nYL9u56UoC+ZRcjNrJJ0yH/Z3b8abD6W+ZMu+OfxsOorsHcAHzCzTtLdcVeTbuGvCv7Eh2Se626g292fCJ4/QDr4k3qe3w0ccvd+d58FvgpcRfLPc8ZS57VguZaUoM9nAfPYC/qn7wJecPf/nbUre3H2G4CvLXdtxeDun3b3dnffRPqcPuzuPw88QnoRekjQ981w9z6gy8y2BpuuIb3uciLPM+kumyvNrDb4dzzzfRN9nrMsdV53Ax8LRt9cCYxkunjOmLsn4gd4L/AS8Arwe2HXU6Tv+KOk/3T7PrA3+Hkv6X7rbwMvB/9cHXatRfju7wK+ETw+D3gSOAj8I1Addn1F+L6XAHuCc/3PQGOSzzPwGeBFYB9wD1CdxPMM3Ev6OsQs6Rb7jUudV9JdN3cEmfYD0qOSzupzNQWCiEjCJaXrRkRElqCgFxFJOAW9iEjCKehFRBJOQS8iknAKehGRhFPQi4gknIJe5DTM7LJgPvCUma0I5k3fHnZdIvnSDVMieTCzPwRSQA3peWj+KOSSRPKmoBfJQzCH0lPAFHCVu8+HXJJI3tR1I5Kf1cBKoI50y14kNtSiF8mDme0mPVXyZmCdu98cckkieas4/SEipc3MPgbMufs/BOsTf8/Mrnb3h8OuTSQfatGLiCSc+uhFRBJOQS8iknAKehGRhFPQi4gknIJeRCThFPQiIgmnoBcRSbj/D0vwkgyaRQ5CAAAAAElFTkSuQmCC\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "xs = uniform.index\n", "p_avg = 0.5\n", "k = 25\n", "n = 100\n", "likelihood2 = compute_likelihood(xs, p_avg, k, n)\n", "\n", "plt.plot(likelihood2)\n", "plt.xlabel('x')\n", "plt.title('Likelihood 2');" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Now we can do the Bayesian update in the usual way, multiplying the priors and likelihoods:" ] }, { "cell_type": "code", "execution_count": 6, "metadata": {}, "outputs": [], "source": [ "posterior = uniform * likelihood1 * likelihood2" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Computing the total probability of the data:" ] }, { "cell_type": "code", "execution_count": 7, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "0.00015958597615453273" ] }, "execution_count": 7, "metadata": {}, "output_type": "execute_result" } ], "source": [ "total = posterior.sum()\n", "total" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "And normalizing the posterior." ] }, { "cell_type": "code", "execution_count": 8, "metadata": {}, "outputs": [], "source": [ "posterior /= total" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Here's what it looks like." ] }, { "cell_type": "code", "execution_count": 9, "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYgAAAEWCAYAAAB8LwAVAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjMsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+AADFEAAAgAElEQVR4nO3de3xddZnv8c8396aXNGkDvULBgshFUQvqHMUL6oCO4JwBLePMgIO38eB4dGbO6DleEHUc5nhXPMoII4IOoDOOnbGKKBcFAVsuCgUKBQotLW3a9JKmTdIkz/ljrd3u7u4kO21WdvbO9/165dW11/qttZ61d5pn/y7rtxQRmJmZFaopdwBmZjYxOUGYmVlRThBmZlaUE4SZmRXlBGFmZkU5QZiZWVFOEJYJSbskHZvRsV8jaX3e61WSXjNGx36HpJ/nvQ5Ji8fi2OnxMntfxkr++6nEv0jaJum3ZQ6tJJJ+KunCcsdRDeT7ICqbpLXAkcAA0A0sBz4QEbsO8XiLgKeA+ojoH5sox1b6x+u6iFgwin0WcQjXJSmA4yJizSjDRNJtJHF+e7T7ThSSXgX8K/D8iOgudzw2vlyDqA5viYhpwEuA04CPlSsQSXXl3L9Szz2BHQ2sPZTkMN7vZ1rbOeS/af78D+YEUUUi4lngp8DJAJLmSVomqVPSGknvzpWVdLqklZJ2Stok6Yvppl+l/25Pm0NekZb/S0mPpE0NN0k6Ou9YIel/SHoceDxv3eJ0uUXSdyV1SHpa0sdy/5ElXSTpTklfktQJXFp4XZKmSPpOeu6HSZJg/va1kl4/2usqdu503R0FIbxJ0pOStkj6v3mxXyrpurw4FqXXXSfps8CrgK+n5/v6Ibwvd0j6fHrdT0k6e+hP/4D344AmuCLv0aWSbkzP3ZU2KS0pLCvpYuDbwCvSa/hUuv3d6e9TZ/r7NS9v36F+F94v6fH0fJ+W9DxJd6Wf042SGoa4ltxn9DVJOyQ9KunMvO23SfqspDuB3cCx6bp3pdtr0vf1aUmb02tuKfi8Lpb0DHBLKe/vpBIR/qngH2At8Pp0eSGwCvh0+vp24BtAE3Aq0AGcmW67C/jzdHka8PJ0eREQQF3eOd4KrAFeANSR1FB+k7c9gJuBNmBK3rrF6fJ3gR8D09PjPwZcnG67COgHPpAee0qRa/xH4Nfp8RcCDwHrh3gPRnNdB507XXdHwbXdmp77qDT2d6XbLiVpQqLYOYDbcmULjlfq+7IXeDdQC/wVsIG0WXiE34nX5L8/Rd6jS4Ee4E3psT8H3D1E2cL343XAFpLaaiPwNeBXJfwuLANmACcBvcAvgWOBFuBh4MIhriX3GX0IqAfeDuwA2vLe42fS49alZfa978BfkvzuHpv+Pvw7cG3B5/VdYCpFfvcm+49rENXhPyRtB+4gSQr/IGkh8Erg7yOiJyIeIPk2+OfpPnuBxZJmR8SuiLh7mOO/F/hcRDwSSfv9PwCn5tci0u2dEbEnf0dJtST/qT8aEV0RsRb4Ql4cABsi4msR0V+4f+ptwGfT468DvjpMrKO5rlLODXB5eu5ngC8DF4xwzBGV+L48HRH/HBEDwDXAXJL+prFwR0QsT499LfCiEvd7B3B1RNwXEb3AR0lqGIvyyhT7Xbg8InZGxCqSBP/ziHgyInaQ1HpfPMw5NwNfjoi9EXEDsBp4c97270TEqvQz3Fsk3i+m59qVxrtUBzYnXRoR3cN8/pOWE0R1eGtEzIyIoyPi/ekv+jygMyK68so9DcxPly8GjgcelbRC0h8Nc/yjga9I2p4mok5AeccCWDfEvrOBhvTcxeIYbt+ceQVlnh6qIKO7rlLOXVjm6TSew1XK+/JcbiEidqeL08bg3Accm6RppkmltcHPIy/m9I/uVkb+PDflLe8p8nq463o20q/8qcLPYLjP8IB40+U6Dky0pfwOTEpOENVrA9AmaXreuqOAZwEi4vGIuAA4Argc+KGkqSRV7kLrgPemSSj3MyUifpNXZqjhcFtIvtXn1zb2xTHCvjkbSZqW8vcvapTXVcq5KXLuDelyN9Cct23OKI5dyvtyqA6IK62ttI/BcSG59vz+p6nALEb3eY7WfEnKe53/GYx0vgPiTfft58AE5aGcQ3CCqFJpU8xvgM9JapL0QpJv198DkPRnktojYhDYnu42QNJPMUjSZpvzTeCjkk5K922RdH6JcQwANwKflTQ9bZb6MHDd8Hse4Mb0/K2SFpD0GRQ1yusq1d+l514IfBC4IV3/AHCGpKPSjs+PFuy3aajzHe77knY03zbE5sdIagRvllRP0mfUWMpxS/B94J2STpXUSNLceE/aRJaVI4C/llSf/t69gGQ4dyn+FfiQpGMkTSOJ94aYoEO4JxoniOp2AUlH3AbgR8AnI+LmdNtZwCpJu4CvAEvTvordwGeBO9MmpZdHxI9Ivo1fL2knSRtySSNqUh8g+Vb7JEk/yfeBq0ex/6dImgaeAn5O0mY+lJKvaxTn/zFwL0lC+AlwFUD6Xt4A/D7d/l8F+30FOC8dhVSs3+Rw3peFwJ3FNqTt+u8n6XN6Nj3H+mJlRysifgl8HPg3kprd84ClY3HsYdwDHEdS6/oscF5EbC1x36tJfl9+RfL708MwXzDsQL5RzqwCSXqAZERaqX8oK5Kki0hGJL2y3LFMRr4xxKwCRcSp5Y7Bqp+bmMzMrCg3MZmZWVGuQZiZWVFV0wcxe/bsWLRoUbnDMDOrKPfee++WiCh6n0zVJIhFixaxcuXKcodhZlZRJA05M4GbmMzMrCgnCDMzK8oJwszMinKCMDOzopwgzMysKCcIMzMrygnCzMyKcoKwUVv+4Eae3e6nM5pVOycIG5WNO/bw/u/dx5dvfqzcoZhZxpwgbFR++chmAG57rIPBQU/0aFbNnCBsVH75SPIo346uXlZt2FnmaMwsS04QVrLdff3c+cRWzj11HhLcunpzuUMysww5QVjJ7nh8C339g7x9yUJeuGAmtzzqBGFWzZwgrGS/eGQT05vqOO2YNl77/HZ+t347W3f1ljssM8uIE4SVZHAwuOXRDl59fDv1tTW87oQjiIDbH+sod2hmlhEnCCvJ79ZvZ8uuXl7/giMBOHleC7OnNXLraicIs2rlBGEl+cUjm6itEa95fvLgqZp0+fbVm+kfGCxzdGaWBScIK8kvH9nMkqNbmdncsG/da59/BDt7+rl/3fYyRmZmWXGCsBHt7NnLo891ccbxBz629lXHz6a2Rvza/RBmVckJwka0IZ136ehZzQesn9FUz9yWJp7p3F2OsMwsY04QNqKN23sAmNsy5aBtc2Y08dzOnvEOyczGQaYJQtJZklZLWiPpI0W2nyHpPkn9ks7LW3+qpLskrZL0e0lvzzJOG96GHUkNYt7MpoO2HdnSxKadvhfCrBplliAk1QJXAGcDJwIXSDqxoNgzwEXA9wvW7wb+IiJOAs4CvixpZlax2vA2bN9DbY04YvrBCWLOjCae29FDhCfuM6s2dRke+3RgTUQ8CSDpeuBc4OFcgYhYm247YJxkRDyWt7xB0magHfBwmTLYuL2HI6c3Ulujg7bNmdHEnr0D7Ozpp2VKfRmiM7OsZNnENB9Yl/d6fbpuVCSdDjQATxTZ9h5JKyWt7OjwSJqsbNixh7kzD+5/gKSJCWCT+yHMqk6WCeLgr5swqnYISXOBa4F3RsRBd2NFxJURsSQilrS3tx98ABsTG3f0MLfl4OYlSGoQAM/tcIIwqzZZJoj1wMK81wuADaXuLGkG8BPgYxFx9xjHZiWKCDbu6GHeEDWIfQnCNQizqpNlglgBHCfpGEkNwFJgWSk7puV/BHw3In6QYYw2gq3dffT1DzJviBrEETMaAdjkGoRZ1cksQUREP3AJcBPwCHBjRKySdJmkcwAknSZpPXA+8C1Jq9Ld3wacAVwk6YH059SsYrWh7bsHYogaRFN9La3N9a5BmFWhLEcxERHLgeUF6z6Rt7yCpOmpcL/rgOuyjM1Ks+8eiCI3yeUcOaPJndRmVch3UtuwctNszC1yk1zOnJYmNrqJyazqOEHYsDbu6KGhroZZUxuGLDPHNQizquQEYcPasH0Pc1uakIqNWk4cOaOJLbuSzmwzqx5OEDas4e6ByJmTbt/c5VqEWTVxgrBhbdy+Z8h7IHJy90K4mcmsujhB2JAGBoNNXb3DjmCCpIkJ4LkdntXVrJo4QdiQNnf1MDAYw45ggv1NTL4Xwqy6OEHYkDakN8mNVINoba6noa7GTUxmVcYJwoa0ccfI90AASOLIGY2esM+syjhB2JD23SQ3Qg0CYO6MKW5iMqsyThA2pA3be5jaUMuMppFnZEkePeoEYVZNnCBsSBt3JENch7tJLmdO2sTkR4+aVQ8nCBvSxh09Q87iWujIGU309g+yY8/ejKMys/HiBGFD2rC9Z8jnQBTyUFez6uMEYUX19g+wZVdvSR3U4EePmlUjJwgrasuuPgCOTJ8YN5IjPd2GWdVxgrCitu5Kps2YPW10CcLTbZhVDycIK2prWoOYNW3o50Dkyz0zYpNndDWrGk4QVtSWUdYgAFqnNrCtuy+rkMxsnDlBWFFb0z/0bcM8Sa5QW3MDnU4QZlXDCcKK2rqrl6b6Gpobakvep3VqPdt2O0GYVQsnCCtq664+Zk1tLOku6py2qQ1s2+0b5cyqRaYJQtJZklZLWiPpI0W2nyHpPkn9ks4r2HahpMfTnwuzjNMOtrW7j9kldlDnzGxO+iA83YZZdcgsQUiqBa4AzgZOBC6QdGJBsWeAi4DvF+zbBnwSeBlwOvBJSa1ZxWoH29rdy6xRdFBD0gfRPxh09fZnFJWZjacsaxCnA2si4smI6AOuB87NLxARayPi98Bgwb5/CNwcEZ0RsQ24GTgrw1itQNLENLoaRGta3iOZzKpDlgliPrAu7/X6dN2Y7SvpPZJWSlrZ0dFxyIHagSKCrbv6aBtlE1Nrcz2A+yHMqkSWCaJY72apjdMl7RsRV0bEkohY0t7ePqrgbGhdvf30DQwye+romphcgzCrLlkmiPXAwrzXC4AN47CvHabR3kWd09aclPe9EGbVIcsEsQI4TtIxkhqApcCyEve9CXijpNa0c/qN6TobB7l5mEbbSb2vBuF7IcyqQmYJIiL6gUtI/rA/AtwYEaskXSbpHABJp0laD5wPfEvSqnTfTuDTJElmBXBZus7GQW4m19F2Us9oqqO2Rk4QZlVi5IcNH4aIWA4sL1j3ibzlFSTNR8X2vRq4Osv4rLhcE9Fo5mECkERrcz2d3e6kNqsGvpPaDpJrYhrNPEw5rc0NbHcNwqwqOEHYQbZ29zG9qY6GutH/erR6wj6zquEEYQfZsqt31M1LOZ6wz6x6OEHYQQ7lLuqctqkN7oMwqxJOEHaQZB6mQ0sQuT4IT9hnVvmcIOwgW3f1jfoeiJxWT9hnVjWcIOwAA4PBtt19zD7EJiZPt2FWPZwg7ADbd/cxGKO/izqnbaon7DOrFk4QdoBDeRZ1vpnNrkGYVQsnCDvAln3zMB3iKCZP2GdWNZwg7AC5mVwP/T4IT9hnVi2cIOwA+2ZyPcQmJk/YZ1Y9nCDsAFu7+6jR/r6E0fKEfWbVwwnCDrBlVx9tUxuorSn2UL/StDY3uJParAo4QdgBOrt7mTXKR40Wap3a4CYmsyrgBGEH2JrWIA5Ha7Mn7DOrBk4QdoCt3X2HPMQ1xxP2mVUHJwg7wOFM9Z3jCfvMqoMThO3T2z9AV0//IQ9xzWmb6gn7zKqBE4Ttk7v7ue0wm5g83YZZdXCCsH1yd1Ef7iim3IR9nm7DrLJlmiAknSVptaQ1kj5SZHujpBvS7fdIWpSur5d0jaQHJT0i6aNZxmmJ3ER9sw+zBtGa1iC2e0ZXs4pWUoKQ9G+S3iyp5IQiqRa4AjgbOBG4QNKJBcUuBrZFxGLgS8Dl6frzgcaIOAV4KfDeXPKw7HR2J9NsHP4wV0/YZ1YNSv2D//+APwUel/SPkk4oYZ/TgTUR8WRE9AHXA+cWlDkXuCZd/iFwpiQBAUyVVAdMAfqAnSXGaodorJqYPGGfWXUoKUFExC8i4h3AS4C1wM2SfiPpnZLqh9htPrAu7/X6dF3RMhHRD+wAZpEki25gI/AM8PmI6Cw8gaT3SFopaWVHR0cpl2LD2NrdR12NmDGl7rCO4wn7zKrDaJqMZgEXAe8C7ge+QpIwbh5qlyLrCgfGD1XmdGAAmAccA/yNpGMPKhhxZUQsiYgl7e3tpVyGDWPrrl7apjaQVOIO3f4J+5wgzCpZSV8VJf07cAJwLfCWiNiYbrpB0sohdlsPLMx7vQDYMESZ9WlzUgvQSdKc9bOI2AtslnQnsAR4spR47dB0dvcd8qNGCyUT9rmT2qySlVqD+HZEnBgRn8slB0mNABGxZIh9VgDHSTpGUgOwFFhWUGYZcGG6fB5wSyS33z4DvE6JqcDLgUdLvio7JFu7+w77Jrmc1mZP2GdW6UpNEJ8psu6u4XZI+xQuAW4CHgFujIhVki6TdE5a7CpglqQ1wIeB3FDYK4BpwEMkieZfIuL3JcZqh2gsJurLmdlc72GuZhVu2CYmSXNIOpKnSHox+/sMZgDNIx08IpYDywvWfSJvuYdkSGvhfruKrbdsdY7BRH05bVMbuH/d9jE5lpmVx0h9EH9I0jG9APhi3vou4H9nFJOVQc/eAXb1Hv48TDkz8ybsO9xObzMrj2ETRERcA1wj6U8i4t/GKSYrg9yIo7HrpK5n70DQ3TfAtMbDGzZrZuUxUhPTn0XEdcAiSR8u3B4RXyyym1WgfRP1jVUn9dT9E/Y5QZhVppH+505N/52WdSBWXlt2JdNsHO48TDm56Ta27e5jYduI3VVmNgGN1MT0rfTfT41POFYu+2sQY9fEBLDNI5nMKtZITUxfHW57RPz12IZj5ZKbh2msm5i2+14Is4o1UhPTveMShZXd1u4+6mvFjKax6S/wjK5mla+UUUw2CYzVPEw5LVPqkdzEZFbJRmpi+nJE/E9J/8nBE+0REecU2c0qUGd332FP852vtka0TKl3E5NZBRupPeHa9N/PZx2IldfWMbyLOqe1ucFNTGYVbKQmpnvTf29PJ9w7gaQmsTp9CJBVia3dvRw9a2yHo3o+JrPKVuojR98MPAF8Ffg6sEbS2VkGZuOrc9fYNjEBtHlGV7OKVuqQlS8Ar42INQCSngf8BPhpVoHZ+OnZO0B338CYNzHNbG7g0ee6xvSYZjZ+Sp3ue3MuOaSeBDZnEI+VwdbcPExjdA9Ejp8qZ1bZRhrF9N/TxVWSlgM3kvRBnE/ynAarAp1jfJNcTuvUBvbsHaBn7wBN9bVjemwzy95ITUxvyVveBLw6Xe4AWjOJyMbdlu5kHqYsRjEBbN+9lzktThBmlWakUUzvHK9ArHxyNYix7qTOzcfU2d3HnJamMT22mWWvpE5qSU3AxcBJwL7/6RHxlxnFZeNoa1qDaMugkxo8H5NZpSq1k/paYA7JE+ZuJ3nCnIenVImt3X001NYwfYyf25Dr0/B0G2aVqdQEsTgiPg50p/MzvRk4JbuwbDxt3dU3pvMw5eyf8ts1CLNKVGqCyH0F3C7pZKAFWJRJRDbuOjOYZgP2NzFt81BXs4pUaoK4UlIr8HFgGfAwcPlIO0k6S9JqSWskfaTI9kZJN6Tb75G0KG/bCyXdJWmVpAfTfhDLwNbuvjEf4grQUFfDtMY6NzGZVaiSGp0j4tvp4u3AsaXsI6kWuAJ4A7AeWCFpWUQ8nFfsYmBbRCyWtJQk6bxdUh1wHfDnEfE7SbPYX4uxMbZ1Vy/HjPE8TDnJfEyuQZhVolLnYpol6WuS7pN0r6Qvp3+0h3M6sCYinkwn9rseOLegzLlA7pkTPwTOVNIQ/kbg9xHxO4CI2BoRA6VelI1O0sQ0tkNcc1qbG+h0gjCrSKU2MV1PMrXGnwDnAVuAG0bYZz6wLu/1+nRd0TIR0Q/sAGYBxwMh6aY0Kf2vYieQ9B5JKyWt7OjoKPFSLN+evgF29w1k0sQEyd3UbmIyq0ylJoi2iPh0RDyV/nwGmDnCPsWGxBQ+dGioMnXAK4F3pP/+saQzDyoYcWVELImIJe3t7SNfhR0kdw/E7Aw6qSEZyeQmJrPKVGqCuFXSUkk16c/bSGZzHc56YGHe6wXAhqHKpP0OLUBnuv72iNgSEbuB5cBLSozVRmFLRndR57Q2N3gUk1mFGjZBSOqStBN4L/B9oC/9uR740AjHXgEcJ+mY9GFDS0lGQOVbBlyYLp8H3BIRAdwEvFBSc5o4Xk0ycsrG2OadPQAcMSO7BLGzp5/+gcFMjm9m2RlpLqbph3rgiOiXdAnJH/ta4OqIWCXpMmBlRCwDrgKulbSGpOawNN13m6QvkiSZAJZHxEg1FjsEm7uSJqYjpmczirh1anKz3PY9e5mdUUe4mWWj5LkVJJ0DnJG+vC0i/mukfSJiOUnzUP66T+Qt95BMHV5s3+tIhrpahjZ39SJl1weRPx+TE4RZZSl1mOs/Ah8kaeZ5GPhgus4qXEdXD7OmNlBXW2p31Ojsn9HVI5nMKk2pNYg3AadGxCCApGuA+4GD7o62yrJ5Zy/tGTUvwf5nQng+JrPKM5qvjfnDWlvGOhArj81dvRwxPbumn9apnvLbrFKVWoP4HHC/pFtJ7l04A/hoZlHZuNnc1cMJcw55LMKI9s/o6iYms0ozYoJIp764A3g5cBpJgvj7iHgu49gsYwODwZZdfZkNcQWYUl9LY12N74Uwq0AjJoiICEn/EREv5eD7GKyCdXb3MTAYmQ1xBZCU3CznJiazilNqH8Tdkk7LNBIbd5u70pvkMuyDgKQfwqOYzCpPqX0QrwXeJ2kt0E3SzBQR8cKsArPs7btJLsMmJkjusdiyqzfTc5jZ2Cs1QZydaRRWFh07s72LOqd9WiNPdnRneg4zG3vDJoj0KW7vAxYDDwJXpdNyWxXINTG1Z9zENHt6I1t29RIRY/7cazPLzkh9ENcAS0iSw9nAFzKPyMbN5q5eZjTV0VRfm+l52qc10ts/yK5ef7cwqyQjNTGdGBGnAEi6Cvht9iHZeNm8s5cjZmT/qO/Z05Ob5Tq6epneVJ/5+cxsbIxUg9g39MRNS9Vnc1dP5iOYANqnJUko9+wJM6sMI9UgXpQ+DwKSkUtT0te5UUwzMo3OMrW5q5clR7dmfp78GoSZVY6RngeRbeO0lU1EJPMwjUcTUzrNt4e6mlWWbOZ4tglvZ08/ff2D49LE1NrcQG2NnCDMKowTxCTVMU5DXAFqa0Tb1AY3MZlVGCeISWrzON0klzN7WqNrEGYVxglikhqvaTZy2qc3ugZhVmGcICap8ZqoLyeZj8nDXM0qiRPEJLV5Zy9T6muZ1ljqdFyHp31aIx3pdBtmVhkyTRCSzpK0WtIaSQc9v1pSo6Qb0u33SFpUsP0oSbsk/W2WcU5GyRDXxnGbG6l9eiN9/YPs7PH9lmaVIrMEIakWuIJkDqcTgQsknVhQ7GJgW0QsBr4EXF6w/UvAT7OKcTIbr7uoc3wvhFnlybIGcTqwJiKejIg+4Hrg3IIy55JMCAjwQ+DM9BGnSHor8CSwKsMYJ63NXb3jNoIJ8hKEO6rNKkaWCWI+sC7v9fp0XdEy6VxPO4BZkqYCfw98argTSHqPpJWSVnZ0dIxZ4JNBx87ecbkHIid3rg7XIMwqRpYJoljjdmEP5VBlPgV8KSJ2DXeCiLgyIpZExJL29vZDDHPy2dM3QFdv/7gNcYVkFBO4BmFWSbIcwrIeWJj3egGwYYgy6yXVAS1AJ/Ay4DxJ/wTMBAYl9UTE1zOMd9LYP8R1/JqYctNtuAZhVjmyTBArgOMkHQM8CywF/rSgzDLgQuAu4DzglkjGQb4qV0DSpcAuJ4ex89yOJEEcOY41iJoaMWtqA1u6fC+EWaXILEFERL+kS4CbgFrg6ohYJekyYGVELAOuAq6VtIak5rA0q3hsv2c6dwOwsLV5XM/r6TbMKkumd0lFxHJgecG6T+Qt9wDnj3CMSzMJbhJbt20PNYJ5M6eM63nbpze6icmsgvhO6klofedu5rZMoaFufD/+2dMa3UltVkGcICahZzp3s7BtfGsPkDxZbsuuPk+3YVYhnCAmoXXbdo97/wMk8zH1DQyyc4+n2zCrBE4Qk0zP3gE27exlYVsZEoRvljOrKE4Qk8z6bXsAOKoMCSI33YafC2FWGZwgJpl129IhrmXog8jVIDzU1awyOEFMMuvKdA8EuAZhVmmcICaZdZ27aayrGdeJ+nJmTqmnrkauQZhVCCeISSYZ4to8bg8KyldTI2ZNa3CCMKsQThCTzLrOPWXpoM5pn97IZjcxmVUEJ4hJJCJY17mbha3j30GdM69lCs+mI6nMbGJzgphEduzZS1dvf1nugcg5qq2ZZzp3+25qswrgBDGJrOtMvrmXNUHMaqa3f9AjmcwqgBPEJFKuab7z5fo/nk5jMbOJywliEinnTXI5uQTxzFYnCLOJzgliEnmmczetzfVMb6ovWwzzW6cg7a/NmNnE5QQxiaxL74Eop8a6WubOaNp3R7eZTVxOEJPI+m17yp4gIOmodg3CbOJzgpgkBgaD9WV6DkSho9qa3UltVgGcICaJTTt72DsQZb2LOueotmY6unrZ0zdQ7lDMbBhOEJPE2q3dQHlHMOXkmrlyo6rMbGJygpgkHt6wE4AT5swocyRw9KypgIe6mk10mSYISWdJWi1pjaSPFNneKOmGdPs9khal698g6V5JD6b/vi7LOCeDVRt2cuSMxrJM811o370Q7ocwm9AySxCSaoErgLOBE4ELJJ1YUOxiYFtELAa+BFyert8CvCUiTgEuBK7NKs7J4qFnd3DK/JZyhwFAa3M90xrrnCDMJrgsaxCnA2si4smI6AOuB84tKHMucE26/EPgTEmKiPsjYkO6fhXQJKn8X30r1O6+fp7o2MVJ8yZGgpDEwjYPdTWb6LJMEPOBdXmv16fripaJiH5gBzCroMyfAPdHxEGzu0l6j6SVklZ2dHSMWeDV5pGNOxkMOHmC1CAAjnaCMJvwskwQxR5ZVjjH87BlJJ1E0uz03mIniIgrI2JJRCxpb28/5ECr3UPPJh3UJ88vfwd1zpEq9sQAAA+bSURBVFGzmlnXuZvBQU/7bTZRZZkg1gML814vADYMVUZSHdACdKavFwA/Av4iIp7IMM6q99CzO5g1tYE5M5rKHco+C9uSab/9dDmziSvLBLECOE7SMZIagKXAsoIyy0g6oQHOA26JiJA0E/gJ8NGIuDPDGCeFhzbs5OT5LWV5DvVQPJLJbOLLLEGkfQqXADcBjwA3RsQqSZdJOictdhUwS9Ia4MNAbijsJcBi4OOSHkh/jsgq1mrWs3eAxzd1TajmJUj6IMAJwmwiq8vy4BGxHFhesO4Tecs9wPlF9vsM8JksY5ssVj/XRf9gcPIEGcGUM2/mFGo87bfZhOY7qavcQxt2ABNrBBNAQ10Nc1umeNpvswnMCaLKPfTsTlqm1LOgtfxzMBVaNLuZJzp2lTsMMxuCE0SVW7VhByfPnzGhOqhzXrRgJg9v2Mnuvv5yh2JmRThBVLG+/kEe3dg14fofck5b1Eb/YPDAuu3lDsXMinCCqGKPb+6ib2CQkyZY/0POS45uRYKVa7eVOxQzK8IJoorduWYLAC85amaZIymuZUo9zz9yOivWdpY7FDMrwgmiii1/8DlOmd/CggnwmNGhnLaojfue3kb/wGC5QzGzAk4QVerZ7Xt4YN12zj5lTrlDGdaSRa109w3w6HNd5Q7FzAo4QVSpnz30HABnnzy3zJEM7/Rj2gDczGQ2ATlBVKmfPriRF8ydwTGzp5Y7lGHNbZnC/JlT3FFtNgE5QVShTTt7WPn0Ns4+eWI3L+WctqiVFWs7ifDU32YTiRNEFbppVdK89KYJ3v+Qs2RRG5u7elnXuafcoZhZHieIKrT8wY0cd8Q0Fh8xvdyhlCTXD/Fb90OYTShOEFWmo6uX3z7VydmnTOzO6XyL26fRMqWelU4QZhOKE0SVufautQwGvLmCEkRNjTj9mDZueXQzvf0D5Q7HzFJOEFVk7ZZuvvmrJ3nrqfN4/pzKaF7KuegPFrG5q5cbV6wrdyhmlnKCqBIRwaX/uYqG2hr+95teUO5wRu0PnjeL0xa18o3bnnAtwmyCcIKoEjc/vInbVnfwoTcczxEzmsodzqhJ4oNnHs/GHT3cuHJ9ucMxM5wgqsKevgE+9Z8Pc8Kc6Vz4iqPLHc4h+2+LZ7Hk6Fa+cesa1yLMJgAniArX2d3HO7/zW57dvodPnXMSdbWV+5FK4oOvP46NO3r4gWsRZmVXV+4A7NA9tqmLi69ZwaadvXzp7S/iZcfOKndIh+2Vi2fz0qNbufxnjzK3pYkzX3BkuUMym7Qy/bop6SxJqyWtkfSRItsbJd2Qbr9H0qK8bR9N16+W9IdZxllpnt2+h6/+8nH++Io76dk7yA3veTl//OIF5Q5rTEjiK0tP5ehZzVx8zUo+f9NqBgY9BYdZOWRWg5BUC1wBvAFYD6yQtCwiHs4rdjGwLSIWS1oKXA68XdKJwFLgJGAe8AtJx0fEpGqY7h8YpKunn81dvTzRsYsnNu/i7qe28psnthIBrzpuNv903guZ2zKl3KGOqQWtzfzwfX/AJ3+8iq/fuobbH+vgdSccwcuOaeOUBS1Ma6ybkM/YNqs2WTYxnQ6siYgnASRdD5wL5CeIc4FL0+UfAl9X8j//XOD6iOgFnpK0Jj3eXWMd5PbdfZz/zTE/LMN9542I/dsDBiMYDBgYDPoHB9k7EPTsHWB338H58NjZU/mfZx7Pf3/JfBa2TdwHAR2upvpaLj/vhSxZ1Mo1d63la7c8zlfSN622RkxvqmNqQx21NaK2RkiQnzLyE4hTiVW7E+bO4GsXvHjMj5tlgpgP5N/1tB542VBlIqJf0g5gVrr+7oJ95xeeQNJ7gPcAHHXUUYcUZE2NOO7IaYe070g03J+mvD9otTVJyZoa0VBbQ12taKyrZUZTPTOm1NE2tYHntU/j2PapNDdMrm6j85cs5PwlC9nZs5eVazt5bNMuunr2snNPP919/QwOBgNpkt3ngEU3T1n1W9iaTStCln9tiv11LPzfOlSZUvYlIq4ErgRYsmTJIf0lmNFUzzfe8dJD2dXG0Yymel53wpG87gR3WpuNlyw7qdcDC/NeLwA2DFVGUh3QAnSWuK+ZmWUoywSxAjhO0jGSGkg6nZcVlFkGXJgunwfcEslTY5YBS9NRTscAxwG/zTBWMzMrkFkTU9qncAlwE1ALXB0RqyRdBqyMiGXAVcC1aSd0J0kSIS13I0mHdj/wPybbCCYzs3JTtTzmccmSJbFy5cpyh2FmVlEk3RsRS4ptq9x5GczMLFNOEGZmVpQThJmZFeUEYWZmRVVNJ7WkDuDpwzjEbGDLGIVTKSbbNU+26wVf82RxONd8dES0F9tQNQnicElaOVRPfrWabNc82a4XfM2TRVbX7CYmMzMrygnCzMyKcoLY78pyB1AGk+2aJ9v1gq95ssjkmt0HYWZmRbkGYWZmRTlBmJlZUZMqQUg6S9JqSWskfaTI9kZJN6Tb75G0aPyjHFslXPOHJT0s6feSfinp6HLEOZZGuua8cudJCkkVPySylGuW9Lb0s14l6fvjHeNYK+F3+yhJt0q6P/39flM54hwrkq6WtFnSQ0Nsl6Svpu/H7yW95LBPGhGT4odkyvEngGOBBuB3wIkFZd4PfDNdXgrcUO64x+GaXws0p8t/NRmuOS03HfgVyaNtl5Q77nH4nI8D7gda09dHlDvucbjmK4G/SpdPBNaWO+7DvOYzgJcADw2x/U3AT0meyPly4J7DPedkqkGcDqyJiCcjog+4Hji3oMy5wDXp8g+BMyVV8jPvR7zmiLg1InanL+8meXpfJSvlcwb4NPBPQM94BpeRUq753cAVEbENICI2j3OMY62Uaw5gRrrcQoU/lTIifkXy3JyhnAt8NxJ3AzMlzT2cc06mBDEfWJf3en26rmiZiOgHdgCzxiW6bJRyzfkuJvkGUslGvGZJLwYWRsR/jWdgGSrlcz4eOF7SnZLulnTWuEWXjVKu+VLgzyStB5YDHxif0MpmtP/fR5TZE+UmoGI1gcIxvqWUqSQlX4+kPwOWAK/ONKLsDXvNkmqALwEXjVdA46CUz7mOpJnpNSS1xF9LOjkitmccW1ZKueYLgO9ExBckvYLk6ZUnR8Rg9uGVxZj//ZpMNYj1wMK81ws4uMq5r4ykOpJq6XBVuomulGtG0uuB/wOcExG94xRbVka65unAycBtktaStNUuq/CO6lJ/t38cEXsj4ilgNUnCqFSlXPPFwI0AEXEX0EQyqV21Kun/+2hMpgSxAjhO0jGSGkg6oZcVlFkGXJgunwfcEmnvT4Ua8ZrT5pZvkSSHSm+XhhGuOSJ2RMTsiFgUEYtI+l3OiYhKfl5tKb/b/0EyIAFJs0manJ4c1yjHVinX/AxwJoCkF5AkiI5xjXJ8LQP+Ih3N9HJgR0RsPJwDTpompojol3QJcBPJCIirI2KVpMuAlRGxDLiKpBq6hqTmsLR8ER++Eq/5/wLTgB+k/fHPRMQ5ZQv6MJV4zVWlxGu+CXijpIeBAeDvImJr+aI+PCVe898A/yzpQyRNLRdV8hc+Sf9K0kQ4O+1X+SRQDxAR3yTpZ3kTsAbYDbzzsM9Zwe+XmZllaDI1MZmZ2Sg4QZiZWVFOEGZmVpQThJmZFeUEYWZmRTlBWNlJGpD0gKSHJP1AUvMo9981yvLfkXRekfVLJH01Xb5I0tfT5fdJ+ou89fNGc75h4nhVOrPqA5Km5K1fNNSMnWNwzrXpfRCllv+ypDNGKPMLSa2HH51NNE4QNhHsiYhTI+JkoA94X/7G9MafzH9XI2JlRPx1kfXfjIjvpi8vAsYkQQDvAD6fXvueMTrmmJHUBrw8nSRuONeSzIRsVcYJwiaaXwOL02/Rj0j6BnAfsFDSBZIeTGsal+fvJOkLku5T8kyL9nTduyWtkPQ7Sf9WUDN5vaRfS3pM0h+l5V8j6aAJ/CRdKulv01rHEuB76bf+N0v6UV65N0j69yL7n6nkmQQPKpnTv1HSu4C3AZ+Q9L0i70OtpH9Oaxg/z9UwJD1P0s8k3ZvGf0K6/i1KnmFyf/qN/sh0/ax0//slfYt0vh5JUyX9JH1vHpL09iIxnAf8LC3fouTZC89PX/+rpHen5ZaRzHtkVcYJwiYMJfNfnQ08mK56Psn0xS8G9gKXA68DTgVOk/TWtNxU4L6IeAlwO8kdpgD/HhGnRcSLgEdI5ubJWUQyMeGbgW9Kahopvoj4IbASeEdEnEpy5+oLcgmJ5M7Vfym4pibgO8DbI+IUktkL/ioivk3yh/XvIuIdRU53HMn03CcB24E/SddfCXwgIl4K/C3wjXT9HSTf9l9MMvX1/0rXfxK4I12/DDgqXX8WsCEiXpTW3H5WJIb/BtybXvsO4BLgO5KWkjxX4p/TbduARkmVPPOxFeEEYRPBFEkPkPzxfYZkyhOAp9N57QFOA26LiI50KvbvkTxABWAQuCFdvg54Zbp8cvot+0GS5pyT8s55Y0QMRsTjJHMSnTDaoNNpG64lmVJ6JvAKDp4u/fnAUxHxWPr6mry4h/NURDyQLt8LLJI0DfgDkmlRHiCZQys33/8C4Kb0Wv+O/dd6Bsl7QkT8BNiWrn+QpBZ1uaRXpQmg0Fzy5i6KiJvT/a4A3lVQdjNj1/RmE8SkmYvJJrQ96TfyfdJ5obrzV43ieLn5Y74DvDUififpIpJ5bArLDPW6VP8C/CfJg4d+kCavfIf6wKn8WXUHgCkkX+i2F75Xqa8BX4yIZZJeQ/IshJyDri0iHpP0UpK5ez4n6ecRcVlBsT0kE9wB+6ZKf0G6vo1k9tCcpnS9VRHXIKxS3AO8WtJsSbUkbd63p9tqSNrLAf6UpLkFkqm9N0qqJ6lB5DtfUo2k55E8tnJ1iXF0pccFICI2kEyp/DGShFToUZJv/4vT13+eF/eoRMRO4ClJ58O+zvsXpZtbgGfT5QvzdvsV6bVLOhtoTZfnAbsj4jrg8ySPsiz0CLA47/WH0nUXAFen7ytKsvkcYO2hXJdNXE4QVhHSaYs/CtxK8vzh+yLix+nmbuAkSfeS9FHkvgl/nCSx3EzyhzrfapI/1D8F3hcRpT569DskfRb5Q1O/B6yLiIeLxN1D0jfxg7T5ZxD4ZonnKuYdwMWSfgesYv9jNi9Nz/FrYEte+U8BZ0i6D3gjSRMewCnAb9Omqv8DfKbIuX5CWuuSdDxJs9LfRMSvSRLPx9JyLwXuLlJ7sgrn2VzNDpOS+yXuj4irRixcYSTdAfzRcE+ek/QVYFlE/HL8IrPx4BqE2WFIay0vJO0IrkJ/w/6RT0N5yMmhOrkGYWZmRbkGYWZmRTlBmJlZUU4QZmZWlBOEmZkV5QRhZmZF/X/9xXeXNV5BFQAAAABJRU5ErkJggg==\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "posterior.plot()\n", "\n", "plt.xlabel('Probability of heads (x)')\n", "plt.ylabel('Probability')\n", "plt.title('Posterior distribution, uniform prior');" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The posterior mean is pretty close to the value I used to construct the data, 0.25." ] }, { "cell_type": "code", "execution_count": 10, "metadata": {}, "outputs": [], "source": [ "def pmf_mean(pmf):\n", " \"\"\"Compute the mean of a PMF.\n", " \n", " pmf: Series representing a PMF\n", " \n", " return: float\n", " \"\"\"\n", " return np.sum(pmf.index * pmf)" ] }, { "cell_type": "code", "execution_count": 11, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "0.2475247524738673" ] }, "execution_count": 11, "metadata": {}, "output_type": "execute_result" } ], "source": [ "pmf_mean(posterior)" ] }, { "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.6" } }, "nbformat": 4, "nbformat_minor": 1 }