{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "# Coin flipping\n", "\n", "***\n", "\n", "The ideas behind flipping a coin are important in computing.\n", "\n", "Computers rely on bits - a bit is a variable that can take on one of two values, one or zero.\n", "\n", "Flipping a coin results in one of two outcomes, heads or tails.\n", "\n", "***" ] }, { "cell_type": "code", "execution_count": 1, "metadata": {}, "outputs": [], "source": [ "# We'll use numpy and scipy.stats to analyse flipping a coin.\n", "import numpy as np\n", "import scipy.stats as ss\n", "\n", "# We'll use this for visualisation.\n", "import matplotlib.pyplot as plt\n", "import seaborn as sns" ] }, { "cell_type": "code", "execution_count": 2, "metadata": {}, "outputs": [], "source": [ "# This just sets the default plot size to be bigger.\n", "plt.rcParams['figure.figsize'] = (12, 8)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "***\n", "\n", "A fair coin will give a head fifty percent of the time and a tail fifty percent of the time.\n", "\n", "We say the probability of a head is 0.5 and the probability of a tail is 0.5.\n", "\n", "We can use the following function to simulate this - giving a 1 for a head and 0 for a tail.\n", "\n", "***" ] }, { "cell_type": "code", "execution_count": 3, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "array([1])" ] }, "execution_count": 3, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# (Number of times to flip a coin, probability of a head, number of times to do this)\n", "np.random.binomial(1, 0.5, 1) # returns the head(1) or tail(0) value" ] }, { "cell_type": "code", "execution_count": 4, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "array([486])" ] }, "execution_count": 4, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# Flip a fair coin 1000 times - how many heads?\n", "np.random.binomial(1000, 0.5, 1) # returns the number of heads" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "***\n", "\n", "How likely are we to see a certain number of heads when flipping a coin however many times?\n", "\n", "***" ] }, { "cell_type": "code", "execution_count": 5, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "0.025225018178380496" ] }, "execution_count": 5, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# (No. of heads, no. of flips, probability of a head)\n", "ss.binom.pmf(500, 1000, 0.5) # ss - scipy stats - " ] }, { "cell_type": "code", "execution_count": 6, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "" ] }, "execution_count": 6, "metadata": {}, "output_type": "execute_result" }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAs8AAAHSCAYAAAAT0iZvAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjAsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+17YcXAAAgAElEQVR4nOzdd3hV14Hu/3epF9R7RaJ303HcG244tuMSlxQ7kzyeFE/zTHI9xbmZTGYmmfubOL9MPE4cJ47jxC2OnSFu2I6NC9iA6B0kISQhhCqSkFA7Z90/OPgSImAjJK1Tvp/n0cPRPutsvecPdF42a69lrLUCAAAAcGZRrgMAAAAAoYLyDAAAAHhEeQYAAAA8ojwDAAAAHlGeAQAAAI8ozwAAAIBHMV4GGWOulfT/S4qW9Li19rsnPR8v6ZeSFkhqlXSHtbbGGLNY0mPHh0n6lrX2JS/nHEp2drYtKyvzEhkAAAAYlvXr17dYa3OGes6caZ1nY0y0pD2Slkqql7RO0l3W2h0njPmqpDnW2i8bY+6U9Clr7R3GmCRJ/dbaQWNMgaTNkgol2TOdcygLFy60FRUVnt40AAAAMBzGmPXW2oVDPedl2sZiSZXW2mprbb+kZyXddNKYmyQ9GXj8gqQrjTHGWttjrR0MHE/QsdLs9ZwAAABAUPFSnosk1Z3wfX3g2JBjAmW5Q1KWJBljlhhjtkvaKunLgee9nBMAAAAIKl7Ksxni2MlzPU45xlq7xlo7U9IiSX9vjEnweM5jJzbmPmNMhTGmorm52UNcAAAAYHR4Kc/1kkpO+L5YUsOpxhhjYiSlSWo7cYC1dqekbkmzPJ7z+Oses9YutNYuzMkZct42AAAAMCa8lOd1kiYbY8qNMXGS7pS0/KQxyyXdE3h8m6S3rbU28JoYSTLGjJc0VVKNx3MCAAAAQeWMS9UFVsq4X9IKHVtW7ufW2u3GmG9LqrDWLpf0M0lPGWMqdeyK852Bl18k6UFjzIAkv6SvWmtbJGmoc47wewMAAABG1BmXqgsmLFUHAACA0XauS9UBAAAAEOUZAAAA8IzyDAAAAHhEeQYAAAA8ojwDAAAAHlGeAQAAAI8ozwAAAIBHlGcAAADAI8ozAAAA4BHlGQAAAPCI8gwAAAB4RHkGAAAAPIpxHQAAMHqeXlM7aue+e0npqJ0bAIIVV54BAAAAjyjPAAAAgEeUZwAAAMAjyjMAAADgEeUZAAAA8IjyDAAAAHhEeQYAAAA8ojwDAAAAHlGeAQAAAI8ozwAAAIBHlGcAAADAI8ozAAAA4BHlGQAAAPCI8gwAAAB4RHkGAAAAPKI8AwAAAB5RngEAAACPKM8AAACAR5RnAAAAwCPKMwAAAOAR5RkAAADwiPIMAAAAeER5BgAAADyiPAMAAAAeUZ4BAAAAjyjPAAAAgEeUZwAAAMAjyjMAAADgEeUZAAAA8IjyDAAAAHhEeQYAAAA8ojwDAAAAHlGeAQAAAI8ozwAAAIBHlGcAAADAI8ozAAAA4BHlGQAAAPCI8gwAAAB4RHkGAAAAPKI8AwAAAB5RngEAAACPKM8AAACAR5RnAAAAwCPKMwAAAOAR5RkAAADwiPIMAAAAeER5BgAAADyiPAMAAAAeUZ4BAAAAjyjPAAAAgEeeyrMx5lpjzG5jTKUx5sEhno83xjwXeH6NMaYscHypMWa9MWZr4M8rTnjNysA5NwW+ckfqTQEAAACjIeZMA4wx0ZIekbRUUr2kdcaY5dbaHScM+6KkdmvtJGPMnZK+J+kOSS2SPmmtbTDGzJK0QlLRCa/7jLW2YoTeCwAAADCqvFx5Xiyp0lpbba3tl/SspJtOGnOTpCcDj1+QdKUxxlhrN1prGwLHt0tKMMbEj0RwAAAAYKx5Kc9FkupO+L5ef3z1+I/GWGsHJXVIyjppzK2SNlpr+0449kRgysZDxhhzVskBAACAMealPA9Vau3ZjDHGzNSxqRx/fsLzn7HWzpZ0ceDrc0P+cGPuM8ZUGGMqmpubPcQFAAAARoeX8lwvqeSE74slNZxqjDEmRlKapLbA98WSXpL0eWtt1fEXWGsPBP7skvS0jk0P+RPW2sestQuttQtzcnK8vCcAAABgVJzxhkFJ6yRNNsaUSzog6U5Jd580ZrmkeyR9KOk2SW9ba60xJl3SK5L+3lq76vjgQMFOt9a2GGNiJd0g6a1zfjcAMMqeXlM7Kue9e0npqJwXADCyznjlOTCH+X4dWyljp6TnrbXbjTHfNsbcGBj2M0lZxphKSQ9IOr6c3f2SJkl66KQl6eIlrTDGbJG0ScdK+U9H8o0BAAAAI83LlWdZa1+V9OpJx755wuNeSbcP8brvSPrOKU67wHtMAAAAwD12GAQAAAA8ojwDAAAAHlGeAQAAAI8ozwAAAIBHlGcAAADAI8ozAAAA4JGnpeoAADgZG8YAiERceQYAAAA8ojwDAAAAHlGeAQAAAI8ozwAAAIBHlGcAAADAI8ozAAAA4BHlGQAAAPCIdZ4BIMxYa9XW3a8Dh49q24EOHe7p1+GjA0qKi9b4rGSVZCQpLoZrJwAwHJRnAAgTHT0D+u2Gev16zX5VNXf/0XNx0VEa8PllJUUZqTA9UWVZyVpcnqnscfFuAgNACKI8A0AIs9ZqU91h/XpNrX6/uUF9g37NK03XPy2brpLMJG2t71B6UqwSY6PVO+BXbVu3alp7tL+1Rx9Vt+qj6lZdMS1XF0/OUXSUcf12ACDoUZ4BIETVtvbom8u3aeXuZiXHReu2BcW6e0mpZhamfTym9Uj/x48T46I1NT9VU/NTJUmdvQN6eXOD3thxSFvqO/SpeUUqyUwa8/cBAKGE8gwAIaZv0KfH3q3Wj96pVGx0lP7x+um6a0mpxsWf3a/01IRY3b1kvHYe7NT/bDqgH79bpfMnZum6mfmKiWZONAAMhfIMACFkdVWL/ul321Td3K1lcwr00LIZyk9LOKdzTi9IVXl2st7Y0agPq1rV0TOguxaXMo0DAIZAeQaAENA/6Ne/v7ZTT6yq0fisJD35Z4t16ZScETt/Qmy0bjyvSFnJ8Xpl60G9uKFety4oVpShQAPAiSjPABDk6tp6dP/TG7S5vkP3XlCmB6+bpoTY6FH5WRdOylbfoE9v7WxSfGyUPjmnUIYCDQAfozwDQBB7fVujvv7CZknSjz+7QNfOyh/1n3n51Fz1Dfj1fmWLEmKidfXM0f+ZABAqKM8AEIQGfH7926vHpmmcV5ymH909f8xWwjDG6NpZ+eod9GvlnmbFx0aP6BQRAAhllGcACDKHe/r11V9v0OqqVt17QZn+4frpY74joDFGN80tVN+gTyu2N6ooPVGTcseNaQYACEasRQQAQaSyqUs3P7JKFTXt+s/bz9O3bpzpbCvtKGN06/xiZY+L04sb6tU74HOSAwCCCeUZAILEyt1N+tQjq3Wkb1DP3He+bl1Q7DqSYqOjdNuCEnUcHdArWw+6jgMAzlGeAcAxa61+/sE+/dkv1qkkM0n/c/9FWjA+w3Wsj5VmJumSKTlav79duxo7XccBAKcozwDgkN9avbz1oL798g4tnZGnF77yCRWlJ7qO9SeunJarvNR4vbThgHr6B13HAQBnKM8A4Ej/oF9Pr6nVh1Wt+tJF5Xr0MwuUFBec93HHREfp9gUl6u4f1O83N7iOAwDOUJ4BwIEjfYP62QfV2nmwUzfMKdA/3TBDUUG+HXZheqIun5arzfUd2nagw3UcAHCC8gwAY6z1SJ9+8m6VDnb06u4lpbpgYrbrSJ5dNiVXhekJ+v3mBvUP+l3HAYAxR3kGgDHUcPiofvxulY4O+PSli8o1szDNdaSzEh1l9Mk5herqG9TqqhbXcQBgzFGeAWCMVDcf0U/fr1ZsdJS+fMlElWYlu440LOOzkjUtP0Xv7mlWTx83DwKILJRnABgDOxo69IvVNUpLjNWfXzpR2SnxriOdk6tn5qs/sH03AEQSyjMAjLKKmjb9ek2tCtISdN/FE5SWGOs60jnLT03QvNIMfVjdqsM9/a7jAMCYoTwDwCj6YG+zXtx4QJNyx+mLF01QUnxwLkU3HFdNz5WR9NbOJtdRAGDMUJ4BYJS8vatJr25r1KyiNH3uE+MVFxNev3LTk+J0/oQsbaxtV2Nnr+s4ADAmwus3OQAEAWut3tjeqLd2HtK8knTdsbBEMVHh+ev2sik5iouJ0hvbG11HAYAxEZ6/zQHAEWutXt16UCv3NGtRWYZuXVCs6CDf/ORcJMXH6NIpOdrV2KWalm7XcQBg1FGeAWCE+K3V8s0NWlXVqk9MyNLNc4sUZcK3OB93wcRspcTH6A+7DrmOAgCjjvIMACPAWqvfb27Qmn1tunhytm6YUyATAcVZkuJionTBpGxVNXer4fBR13EAYFRRngHgHB2fqnG8OF87Mz9iivNxi8syFRcTpQ8q2XUQQHijPAPAObDWasX2Q1pV1aoLJmZFZHGWpMS4aC0cn6Et9YfVcXTAdRwAGDWUZwA4B3/Y1aT39jZrcXmmls2OnKkaQ7lwYraslVZXcfUZQPiiPAPAMK3c3aS3dzVpwfgM3XheYUQXZ0nKSI7TrKI0rd3Xpt4Bn+s4ADAqKM8AMAxr9rXqjR2HNLckXZ+aFxmranhx8eRs9Q36VVHT5joKAIwKyjMAnKVtBzq0fFODpuWn6Nb5xRTnExRnJKk8O1mrqlrl81vXcQBgxFGeAeAsVLcc0fMVdSrJTNKdi0rDegOU4bp4UrY6jg5o64EO11EAYMRRngHAo50HO/XUh/uVkRynz58/XnEx/AodypT8FOWMi9cHe5tlLVefAYQXfvMDgAd1bT265+drlRAbrS9cUKak+BjXkYJWlDG6aHK2Gjp6Vc2W3QDCDOUZAM6go2dA9z6xVn2Dft17QZnSk+JcRwp6c0vSlRQXrTXVra6jAMCIojwDwGkM+Pz62tMbVNvWo8c+t0B5qQmuI4WE2OgozS/N0I6DnerqZdMUAOGD8gwAp2Ct1f9evl0fVLbo3z41W0smZLmOFFIWlWXKb6UN+9tdRwGAEUN5BoBT+PmqGj29plZfuWyibl9Y4jpOyMlJiVd5drLW7W+XnxsHAYQJyjMADOEPOw/pO6/s0DUz8/T1q6e6jhOyFpdlqq27X9XN3DgIIDxQngHgJDsPduovn9momYWpeviOuYpiLedhm1mYqqS4aK3dx42DAMID5RkATnC4p1/3PVWhcQkxevzzi5QUx5J05yKGGwcBhBnKMwAE+PxWf/nsJjV29OrRzy5Qfhora4wEbhwEEE4ozwAQ8PCbe/Tenmb9842zNL80w3WcsMGNgwDCCeUZACSt2N6oH71TqTsWluiuxaysMdKO3zhY1XzEdRQAOCeUZwARr7LpiP72+c06rzhN/3zTTBnDDYIj7fiNg+v2tbmOAgDnhPIMIKId6RvUnz9VofiYKD362QVKiI12HSksceMggHDhqTwbY641xuw2xlQaYx4c4vl4Y8xzgefXGGPKAseXGmPWG2O2Bv684oTXLAgcrzTG/NBwqQfAGLPW6sHfblFNa4/+6+55KkxPdB0prC0sy5DfSpvrO1xHAYBhO2N5NsZES3pE0nWSZki6yxgz46RhX5TUbq2dJOlhSd8LHG+R9Elr7WxJ90h66oTXPCrpPkmTA1/XnsP7AICz9vTaWr285aAeWDpFF0zMdh0n7OWmJKg4I1Eba1l1A0Do8nLlebGkSmtttbW2X9Kzkm46acxNkp4MPH5B0pXGGGOt3WitbQgc3y4pIXCVukBSqrX2Q2utlfRLSTef87sBAI+2N3Ton3+/Q5dMydFXLp3oOk7EmFeSroMdvTrYcdR1FAAYFi/luUhS3Qnf1weODTnGWjsoqUNS1kljbpW00VrbFxhff4ZzSpKMMfcZYyqMMRXNzc0e4gLA6R3pG9T9T29URlKsHv70eewgOIbmFKcrykibag+7jgIAw+KlPA/1qXLyQp2nHWOMmaljUzn+/CzOeeygtY9Zaxdaaxfm5OR4iAsAp2at1T+8uFX7W7v1wzvnKWtcvOtIESU5PkZT81O1qf4waz4DCEleynO9pBMXPS2W1HCqMcaYGElpktoC3xdLeknS5621VSeMLz7DOQFgxD27rk7LNzfogaVTtGTCyf9BhrEwryRdXb2DqmpizWcAocdLeV4nabIxptwYEyfpTknLTxqzXMduCJSk2yS9ba21xph0Sa9I+ntr7arjg621ByV1GWPOD6yy8XlJ/3OO7wUATmt3Y5e+tXy7Lp6cra9eNsl1nIg1LT9FibHR2ljH1A0AoeeM5Tkwh/l+SSsk7ZT0vLV2uzHm28aYGwPDfiYpyxhTKekBSceXs7tf0iRJDxljNgW+cgPPfUXS45IqJVVJem2k3hQAnKx3wKe/fGajUhJi9P1Pz2Wes0Mx0VGaXZym7Q0d6hvwuY4DAGclxssga+2rkl496dg3T3jcK+n2IV73HUnfOcU5KyTNOpuwADBc331tl3Yf6tITX1iknBTmObs2vyRda/e1aVtDpxaMz3AdBwA8Y4dBAGHvnd1N+sXqGt17QZkun5p75hdg1JVkJikrOY41nwGEHMozgLDWcqRPX//NFk3NS9GD101zHQcBxhjNK01XdUu32nv6XccBAM8ozwDClrVW33hhizp7B/TDu+YpITbadSScYG7Jsekam7lxEEAIoTwDCFtPfbRfb+9q0j9cN01T81Ncx8FJMpPjVJaVpA21h2VZ8xlAiKA8AwhLlU1H9K+v7NRlU3N0zwVlruPgFOaWZKjlSJ8aOnpdRwEATyjPAMLOoM+vv/3NZiXGRes/bpujY8vJIxjNKkxVlJG21DN1A0BooDwDCDs/ea9am+sO6zs3z1JuSoLrODiNpPgYTc5N0Zb6DrbrBhASKM8AwsqOhk794K09umFOgW6YU+g6DjyYU5ymjqMDqmvrcR0FAM6I8gwgbPQN+vTA85uUnhSnf7mJPZhCxYyCVMVEGW1m6gaAEEB5BhA2fviHvdrV2KXv3jJbGclxruPAo/jYaE3LT9HWA53y+Zm6ASC4edqeGwCC3Ybadj26skqfXlisK6fn6ek1ta4j4SzMKU7XtoZOVbcccR0FAE6LK88AQl7vgE9/95vNKkhL1EM3zHAdB8MwNT9F8TFR2lLf4ToKAJwW5RlAyPuvt/equrlb/37LbKUkxLqOg2GIjY7SjIJUbW/oUN+gz3UcADglpm0ACGnbGzr043erdduCYl0yJcd1nGEbrWkma/e1aXF55qice6TNKU7XxrrDem9Pi5bOyHMdBwCGxJVnACFrwOfXN17YoszkOD20jOkaoW5S7jglxUVr+eYG11EA4JQozwBC1k/fr9b2hk79y00zlZbEdI1QFx1lNKswTW/tOKSe/kHXcQBgSEzbAODUcKcrNHf16b/e3quZhalq6x5gdY0wMac4TWtr2vTWzibdeB6b3AAIPlx5BhBy/NbqxY31iok2FKwwU5adrLzUeC3fxNQNAMGJ8gwg5Kzd16b9rT1aNruQ1TXCTJQxWja7UO/taVZX74DrOADwJyjPAEJKV++AVmxv1KSccZpfmu46DkbBsjn56vf59dbOQ66jAMCfoDwDCCmvbD0on9/qxrmFMsa4joNRMK8kQwVpCXply0HXUQDgT1CeAYSMvU1d2lLfoUum5Ch7XLzrOBglUVFG188u0Ht7WtTJ1A0AQYbyDCAkDPj8Wr6pQVnJcbo0hDdDgTfL5hQcm7qxg6kbAIIL5RlASHhvT7Nau/t149xCxUbzqyvczStJVyFTNwAEIT6BAAS9liN9WrmnWXOK0zQ5N8V1HIwBYwJTN/Y2q+MoUzcABA/KM4CgZq3V8k0NignMg0XkWDanQAM+qzeZugEgiFCeAQS1rQc6VNl8RFfPyFMqazpHlLkl6SpKT9SrW5m6ASB4UJ4BBK3+Qb9e29aogrQELZmQ5ToOxpgxRsvmFOj9vc3q6GHqBoDgQHkGELSOz3e9YU6holjTOSJdP/vY1I03djS6jgIAkijPAIJUe0+/3tvTrNlFaSrPTnYdB46cV5ymovREvcLUDQBBgvIMICi9tvWgjJGum5XvOgocMsbohjkF+mBvC1M3AAQFyjOAoFPVfETbGjp1yZQcpSfFuY4Dx5bNKdCg32oFUzcABAHKM4Cg4vNbvbLloNKTYnXJZHYShDS7KE3FGYlsmAIgKFCeAQSVdTVtauzs1XWzCthJEJL+36obqypbdLin33UcABGOTyYAQeNov09v7jikCdnJmlWY6joOgsgNsws16Ld6YzsbpgBwi/IMIGis3N2k3gGfls0pkGFpOpxgVlGqSjOT9DKrbgBwjPIMICi0d/drdXWr5pdmqCAt0XUcBJkTp260dzN1A4A7lGcAQeHNnYdkJF01I891FASpZbML5PNbrdjOqhsA3KE8A3DuQPtRbao7rAsnZSstMdZ1HASpmYWpGp+VxIYpAJyiPANwylqr17YdVFJctC6dwtJ0ODVjjJbNLtDqqla1MXUDgCOUZwBO7TnUpeqWbl0xLVcJsdGu4yDILZvD1A0AblGeATjj81u9tq1RWclxWlye6ToOQsCMglSVZyezYQoAZyjPAJz57fp6NXX16eqZ+YqJ4tcRzswYo+tn52t1VYtaj/S5jgMgAvFpBcCJ3gGf/vPN3SrJSGRDFJyVZbML5bfSCjZMAeAA5RmAE7/6aL8Odfbpmln5bIiCszK9IEUTspP1ytYG11EARCDKM4Axd6RvUP+9skoXT87WhOxxruMgxBzfMOXDqla1MHUDwBijPAMYc79YtU9t3f3626unuo6CEHX97AL5rfT6NlbdADC2KM8AxlRHz4B+8l61rpqep7kl6a7jIERNy0/RhJxkvcqGKQDGGOUZwJj66fvV6uod1ANLp7iOghBmjNENswv0UXWrmruYugFg7FCeAYyZliN9+vmqfbphToFmsMIGztGyOcdW3XidDVMAjCHKM4Ax8+OVVeod8Omvr+KqM87dlLxxmpQ7Tq9sYdUNAGOH8gxgTDR29OqXH+3XLfOLNSmXFTZw7o5tmFKgtfva1NTV6zoOgAhBeQYwJn70zl5Za/VXV052HQVh5IY5x1bdWMGqGwDGCOUZwKg72HFUz6+r120LSlSSmeQ6DsLIlLwUTc4dp5e3sOoGgLFBeQYw6n7ybrX81uqrl010HQVhaNmcAq2taVNTJ1M3AIw+yjOAUdXU2atn1tbqlvlFXHXGqFg2u0CWVTcAjBHKM4BR9dh71Rr0W33t8kmuoyBMTc5L0ZQ8pm4AGBuUZwCjpuVIn369plY3nVeo8VnJruMgjC2bXah1NW06xNQNAKOM8gxg1Dz+/j71Dvr0tSu46ozRtWxOvqyVXmO7bgCjjPIMYFS0d/frqQ9rdMOcQk3MYV1njK5JuSmalp+iV7cy7xnA6KI8AxgVP1+1T939Pv0FV50xRq6fXaB1+9vU2MHUDQCjh/IMYMR1HB3QL1bV6LpZ+ZqSl+I6DiLE9YFVN17bxtQNAKOH8gxgxD25ukZdfYO6n6vOGEOTcsdpWn6KXmHVDQCjyFN5NsZca4zZbYypNMY8OMTz8caY5wLPrzHGlAWOZxlj3jHGHDHG/Oik16wMnHNT4Ct3JN4QALeO9vv0i9U1umJarmYWprmOgwhzw5wCVexvZ+oGgFFzxvJsjImW9Iik6yTNkHSXMWbGScO+KKndWjtJ0sOSvhc43ivpIUl/d4rTf8ZaOzfw1TScNwAguDxfUae27n59hd0E4cD1swskSa+y6gaAUeLlyvNiSZXW2mprbb+kZyXddNKYmyQ9GXj8gqQrjTHGWtttrf1Ax0o0gDA36PPrp+9Xa8H4DC0qy3QdBxFoQs44TS9I1SuUZwCjxEt5LpJUd8L39YFjQ46x1g5K6pCU5eHcTwSmbDxkjDFDDTDG3GeMqTDGVDQ3N3s4JQBXXtl6UPXtR/XlS7nqDHdumFOg9fvb1XD4qOsoAMKQl/I8VKm1wxhzss9Ya2dLujjw9bmhBllrH7PWLrTWLszJyTljWABuWGv16MoqTc4dpyuncQsD3Dk+deO1baz5DGDkeSnP9ZJKTvi+WFLDqcYYY2IkpUlqO91JrbUHAn92SXpax6aHAAhRK/c0a1djl/780omKihryP5KAMVGenayZhal6ZcvJH1UAcO68lOd1kiYbY8qNMXGS7pS0/KQxyyXdE3h8m6S3rbWnvPJsjIkxxmQHHsdKukHStrMNDyB4/HhllQrSEnTjeYWuowC6fnaBNtQe1gGmbgAYYWcsz4E5zPdLWiFpp6TnrbXbjTHfNsbcGBj2M0lZxphKSQ9I+ng5O2NMjaTvS7rXGFMfWKkjXtIKY8wWSZskHZD005F7WwDG0obadq3Z16YvXTxBcTEsHw/3lh2fusGNgwBGWIyXQdbaVyW9etKxb57wuFfS7ad4bdkpTrvAW0QAwe7HK6uUlhirOxeVnHkwMAbKspM1qyhVL285qC9dPMF1HABhhEtEAM5JZdMRvbnzkO75xHglx3v69zgwJpbNLtSmusOqb+9xHQVAGKE8AzgnP/tgn+Kio3TPBWWuowB/5P9N3WDVDQAjh/IMYNjauvv14oZ63TK/SFnj4l3HAf5IaVaSZhel6WVW3QAwgijPAIbt1x/tV9+gX392YbnrKMCQbphToM31Hdrf2u06CoAwQXkGMCx9gz49+eF+XTolR5PzUlzHAYZ049xCGSP9biNXnwGMDMozgGH5/eaDajnSpy9dzFVnBK+CtESdX56l3206oNNsPwAAnlGeAZw1a60ef79aU/NSdNGkbNdxgNP61Lwi7Wvp1ub6DtdRAIQByjOAs7a6qlW7Grv0xYvKZQxbcSO4XTs7X3ExUfrdxgOuowAIA5RnAGft8ferlT0uTjfOZStuBL/UhFhdNT1Xv9/coAGf33UcACGO8gzgrFQ2HdE7u5v12fPHKyE22nUcwJOb5xaptbtfH1S2uI4CIMRRngGclZ+v2qe4mCh99vzxrqMAnl02NVfpSbFM3QBwzthLFwgzT6+pHZXz3r2kVO3d/frt+np9am6RstkUBSEkLiZKy2YX6MUNB9TdN8hW8gCGjSvPADx7rqJOfYN+feGiMg95CpAAACAASURBVNdRgLN287wiHR3w6Y0dbNcNYPgozwA88fmtnvpwv5aUZ2pafqrrOMBZW1CaoeKMRL3EhikAzgHlGYAnf9h5SAcOH9W9F5S5jgIMS1SU0c1zi/TB3mY1dfW6jgMgRFGeAXjy5Ic1KkhL0NIZea6jAMN287xC+e2xHTIBYDgozwDOqKmzV6sqW/XZ88crJppfGwhdk3JTNLsoTS9trHcdBUCI4nZjAGf0YXWroqOMYqOjRm01D2Cs3Dq/SN/6/Q7taOjUjELm7wM4O1xCAnBavQM+baw9rDlFaRrH8l4IAzfNLVJcdJR+s77OdRQAIYjyDOC0NtS2q9/n1ycmZrmOAoyIjOQ4LZ2Rp99tPKD+QbbrBnB2KM8ATslvrT6qblVJRqKKM5JcxwFGzO0Li9XeM6A/7DzkOgqAEEN5BnBKlU1H1HKkn6vOCDsXT85RfmqCfrOeGwcBnB3KM4BT+qi6VcnxMZpVmOY6CjCioqOMbl1QpJW7m3SokzWfAXjH3T8AhtTe06/djV26dGoOy9NhTI3Wii53Lyn9o+9vW1CiR96p0osbDugrl00clZ8JIPzwiQhgSOtq2iRJi8syHScBRkd5drIWlWXoNxV1sta6jgMgRFCeAfwJn99qfU27puanKD0pznUcYNTcvrBE1S3d2lDb7joKgBDBtA0Af2LnwU519Q1qcTlXnRE+hpoO0jfoU1x0lP791V26ZX7xsM578nQQAOGNK88A/sTafW1KT4zVlLwU11GAURUfE63ZRWnacqCDNZ8BeEJ5BvBHWo70qbL5iBaWZSrKGNdxgFE3f3yG+gf92nagw3UUACGA8gzgj6zb16YoIy0sy3AdBRgTZVlJyh4Xp7WBm2QB4HQozwA+Nujza31tu6YXpCo1IdZ1HGBMGGO0uCxTtW09Othx1HUcAEGO8gzgY9saOtXT7+NGQUSc+eMzFBNltHYfV58BnB7lGcDH1u5rVWZynCbmjHMdBRhTSXExmlOcpo11h9U34HMdB0AQozwDkCQd6uxVTWuPFnOjICLU4vIs9Q/6tan+sOsoAIIY5RmAJGltTZuio4zmj+dGQUSmkoxEFaQlaE11GzsOAjglyjMADfj82lR7WDMKUjUunr2TEJmMMVpSnqXGzl7VtfW4jgMgSFGeAWhHQ6eODvi0qIwbBRHZzitJU3xMlNZw4yCAU6A8A9C6/W3KSIrVhJxk11EAp+JjojW3JF1bD3Sop2/QdRwAQYjyDES41iN9qm7u1oLx3CgISNKSCVka9Futr213HQVAEKI8AxFufW27jKQF3CgISJLyUxM0PitJa/e1yc+NgwBOQnkGIpjPb7Vhf7um5KUoLZEdBYHjlpRnqrW7X1VNR1xHARBkKM9ABNt7qEudvYNaWMZVZ+BEswrTNC4+RquqWlxHARBkKM9ABKvY367k+BhNy091HQUIKjHRUTp/Qqb2HDqips5e13EABBHKMxChunoHtKuxU/NL0xUdxY2CwMkWl2cpJspoVVWr6ygAggjlGYhQG2oPy2+lheNZ2xkYyrj4GM0tSdfG2nZ1s2wdgADKMxCBrLWqqGlTWVaSclLiXccBgtaFk7I16LdaW8OmKQCOoTwDEaimtUet3f1cdQbOIC81QZNzx+mjqlYN+vyu4wAIApRnIAKt39+m+JgozSpKcx0FCHoXTspWV9+gthzocB0FQBCgPAMRpm/Qp20HOjW7KE1xMfwKAM5kcu445abEa1VliyybpgARj09OIMJsO9Cpfp+fHQUBj4wxunBStg529Kq6pdt1HACOUZ6BCLOhtl1ZyXEqzUxyHQUIGXNL0pUUF61VlWyaAkQ6yjMQQdq6+7WvpVvzx2fIGNZ2BryKjY7S+ROytKuxS01dbJoCRDLKMxBBNta2y0iaV5LuOgoQcs6fkKXYaKP39jS7jgLAIcozECH81mpDbbsm5oxTelKc6zhAyBkXH6NFZZnaVHdY7d39ruMAcITyDESImtZutfcMaP54rjoDw3Xx5BwZGb1fydVnIFJRnoEIsWH/YcXHRGlGAWs7A8OVlhireaXpqqhpV1fvgOs4ABygPAMR4Njazh2s7QyMgEun5Mjnt6y8AUQoPkWBCLA9sLbz/FLWdgbOVda4eM0uTtNH+9rU0z/oOg6AMUZ5BiLA+tp2ZSbHaXwWazsDI+GyKbnqH/Trw+pW11EAjDHKMxDmPl7buZS1nYGRkp+WoOn5KVpd2aruPq4+A5GE8gyEuY11x9Z2nl/KKhvASLp0aq6ODvj09Jpa11EAjCHKMxDG/NZqw/52TchJZm1nYISVZiZpQk6yHnu/Wr0DPtdxAIwRyjMQxva39hxb25kbBYFRceW0PDV39elXH+13HQXAGPFUno0x1xpjdhtjKo0xDw7xfLwx5rnA82uMMWWB41nGmHeMMUeMMT866TULjDFbA6/5oWEyJjDiNuxvV3xMlGYWsrYzMBrKs5N10aRsPbqyirnPQIQ4Y3k2xkRLekTSdZJmSLrLGDPjpGFflNRurZ0k6WFJ3wsc75X0kKS/G+LUj0q6T9LkwNe1w3kDAIbWN+jTVtZ2BkbdA1dPUWt3v578sMZ1FABjwMsn6mJJldbaamttv6RnJd100pibJD0ZePyCpCuNMcZa222t/UDHSvTHjDEFklKttR9aa62kX0q6+VzeCIA/tr2BtZ2BsTC/NEOXT83RY+9Vs+sgEAG8lOciSXUnfF8fODbkGGvtoKQOSVlnOGf9Gc4pSTLG3GeMqTDGVDQ3N3uIC0A6NmWDtZ2BsfHA0qk63DOgn39Q4zoKgFHmpTwPNRfZDmPMsMZbax+z1i601i7Myck5zSkBHNfe3a9q1nYGxszs4jQtnZGnxz+oVkcPV5+BcOalPNdLKjnh+2JJDacaY4yJkZQmqe0M5yw+wzkBDNOGunZJ0jzWdgbGzANLp6ird1A/fb/adRQAo8hLeV4nabIxptwYEyfpTknLTxqzXNI9gce3SXo7MJd5SNbag5K6jDHnB1bZ+Lyk/znr9AD+xIlrO2ewtjMwZqYXpGrZ7AI9sWqf2rr7XccBMErOWJ4Dc5jvl7RC0k5Jz1trtxtjvm2MuTEw7GeSsowxlZIekPTxcnbGmBpJ35d0rzGm/oSVOr4i6XFJlZKqJL02Mm8JiGzH13ZewI2CwJj766smq2fAp5+8W+U6CoBREuNlkLX2VUmvnnTsmyc87pV0+yleW3aK4xWSZnkNCsCbDbXtimNtZ8CJyXkpunlukX6xukb3XlimgrRE15EAjDAWfwXCSE//IGs7A449sHSKrJV+8OZe11EAjAI+XYEw8vq2RvUPsrYz4FJJZpI+94nx+s36Ou091OU6DoARRnkGwsgL6+uVmRynMtZ2Bpz62uWTlBwXo++9vtt1FAAjjPIMhIn69h6trmrV/NJ01nYGHMtMjtOXL5uot3Ye0rqa063cCiDUUJ6BMPHihgOSpHlM2QCCwhcuLFNuSry++9ounWb1VgAhhvIMhAFrrV5YX68LJmaxtjMQJJLiYvQ3S6do/f52vbnjkOs4AEYI5RkIA+tq2lXb1qPbFhSfeTCAMXP7gmJNyEnW917fpUGf33UcACOA8gyEgRfW1yk5LlrXzsp3HQXACWKio/SNa6apqrlbL6yvdx0HwAigPAMhrqd/UK9sOajrZxcoKc7TvkcAxtA1M/M0vzRdD7+1R0f7fa7jADhHlGcgxK3Y3qjufh9TNoAgZYzRg9dN16HOPj2xep/rOADOEeUZCHEvrK9XaWaSFpVluo4C4BQWl2fqqum5enRlldq7+13HAXAOKM9ACDu+tvOt84sVFcXazkAw+/o109TdN6hH3ql0HQXAOaA8AyHspQ0HZK10y/wi11EAnMHU/BTdOr9Yv/xwv+rbe1zHATBMlGcgRFlr9cKGep0/IVMlmWzHDYSCv1k6RcZI339zj+soAIaJ8gyEqDX72rS/tUefXljiOgoAjwrTE3XvhWV6aeMB7TzY6ToOgGGgPAMh6vmKOqXEx+i6WQWuowA4C1+9dJJS4mP0vdd3uY4CYBgoz0AI6uwd0KtbD+qTcwuVGBftOg6As5CWFKuvXT5JK3c3a3VVi+s4AM4S5RkIQS9vPqjeAb/uYMoGEJLuuaBMBWkJ+t7ru2WtdR0HwFmgPAMh6LmKOk3NS9Gc4jTXUQAMQ0JstP5m6RRtrjus17c1uo4D4CxQnoEQs7uxS5vrDuvTi0pkDGs7A6Hq1vnFmpI3Tv9nxW4N+Pyu4wDwiPIMhJjnK+oUG21089xC11EAnIPoKKNvXDNN1S3der6iznUcAB5RnoEQ0j/o10sbD+iq6XnKGhfvOg6Ac3Tl9FwtHJ+hH7y1Vz39g67jAPCA8gyEkD/sPKS27n59ehE3CgLhwBijB6+bpuauPj2xqsZ1HAAeUJ6BEPJ8RZ3yUxN0yeQc11EAjJCFZZlaOiNPP15ZpbbuftdxAJwB5RkIEY0dvXp3T7NuW1Cs6ChuFATCyTeumaru/kE98k6l6ygAzoDyDISI326ol99Kty8sdh0FwAibnJei2xYU66kP96u+vcd1HACnQXkGQoDfb/XM2lqdPyFT47OSXccBMAr++qopMkb6/pt7XEcBcBqUZyAEvLe3WfXtR/WZJeNdRwEwSgrTE3XvhWV6aeMB7TzY6ToOgFOgPAMh4Ok1tcpKjtM1M/NdRwEwir566SSlxMfoP17f5ToKgFOgPANB7mDHUf1hV5M+vahEcTH8lQXCWVpSrL56+SS9s7tZH1W3uo4DYAh8EgNB7rl1dfL5re5aVOo6CoAxcO8FZcpPTdB3X9sla63rOABOQnkGgtigz6/n1tXpkik5Ks1Kch0HwBhIiI3WA0unaFPdYa3Y3ug6DoCTUJ6BIPbO7mYd7OjVZ5Zw1RmIJLfML9Kk3HH6jxW7Nejzu44D4AQxrgMAOLVfr9mvvNR4XTkt13UUAKfw9JraUTnvN66ZqvueWq/frK/XXYv5BzQQLLjyDASpurYevbunWXcsKlVMNH9VgUizdEaeFozP0MNv7tHRfp/rOAAC+EQGgtSz62plJN25qMR1FAAOGGP0v66dpqauPj31UY3rOAACKM9AEBrw+fXcunpdMS1PhemJruMAcGRxeaYumZKj/15Zpa7eAddxAIjyDASlN7YfUsuRPm4UBKCvXz1Vh3sG9Pj7+1xHASDKMxCUfrF6n0oyE3XJlBzXUQA4Nrs4TdfNytfj71errbvfdRwg4lGegSCztb5D62rade8F5YqOMq7jAAgCDyydoqMDPv343SrXUYCIR3kGgswTq/YpOS5aty8sdh0FQJCYnJeiT80r1pOra9TY0es6DhDRKM9AEGnq6tXvtzTo9oUlSk2IdR0HQBD566smy2+t/uvtva6jABGN8gwEkV99VKtBv9W9F5S5jgIgyJRkJumuxaV6bl2d9rd2u44DRCx2GASCRO+AT0+v2a8rpuaqLDvZdRwAjg21c2FheqKMkf7q2U369MLhrQF/N6v4AOeEK89AkPj95ga1HOnXn11U7joKgCCVmhCrT0zI0ua6w2rsZO4z4ALlGQgC1lo9sapGU/NSdMHELNdxAASxSybnKC4mSm/tOOQ6ChCRKM9AEFizr007DnbqCxeWyRiWpwNwaknxMbp4crZ2HOxUXVuP6zhAxKE8A0HgiVX7lJEUq5vnFbmOAiAEXDgxW0lx0XpzJ1efgbFGeQYcq23t0Zs7DunuJaVKiI12HQdACIiPjdZlU3NV2XREVc1HXMcBIgrlGXDssferFBMVpc9/osx1FAAhZEl5ptISY/XmjkOy1rqOA0QMyjPgUFNnr56vqNdtC4uVl5rgOg6AEBIbHaXLp+aqtq1Huxu7XMcBIgblGXDo8Q/2adDn15cvmeg6CoAQtGB8hjKT4/TGjkPyc/UZGBOUZ8CRwz39+tVH+3XjeYUqzUpyHQdACIqOMrpqep4aO3u19UCH6zhARKA8A478YnWNevp9+splk1xHARDC5hSnKT81QW/tOCSfn6vPwGijPAMOHOkb1BOrarR0Rp6m5qe4jgMghEUZo6Uz8tTa3a8Nte2u4wBhj/IMOPDMmlp1HB3QVy9jrjOAczctP0UlGYl6e1eTBnx+13GAsEZ5BsZY74BPj71frQsnZWleaYbrOADCgDFGV8/MV8fRAa3d1+Y6DhDWKM/AGHthfb2au/r0NeY6AxhBE3PGaUJOslbublLfoM91HCBsUZ6BMTTg8+sn71Vpbkm6PjExy3UcAGHm6hn56u73aXVVq+soQNiKcR0AiCTPratTXdtR/fONM/XM2jrXcQCEmdLMJE3PT9H7e5u1pDxTSXF8zAMjjSvPwBg52u/TD/+wVwvHZ+jyqbmu4wAIU1fNyFPfgF/v721xHQUIS5RnYIw8+WGNmrr69I1rp8kY4zoOgDBVkJaoOcVpWl3Vos7eAddxgLDjqTwbY641xuw2xlQaYx4c4vl4Y8xzgefXGGPKTnju7wPHdxtjrjnheI0xZqsxZpMxpmIk3gwQrDqODujRlVW6bGqOFpdnuo4DIMxdNT1PPr/V27uaXEcBws4Zy7MxJlrSI5KukzRD0l3GmBknDfuipHZr7SRJD0v6XuC1MyTdKWmmpGsl/XfgfMddbq2da61deM7vBAhij71XpY6jA/r6NVNdRwEQAbLGxWtxeaYqatrU3NXnOg4QVrxceV4sqdJaW22t7Zf0rKSbThpzk6QnA49fkHSlOfb/0jdJetZa22et3SepMnA+IGI0dfXq5x/U6JPnFWpmYZrrOAAixBXT8hQTHaUV2xtdRwHCipfyXCTpxGUB6gPHhhxjrR2U1CEp6wyvtZLeMMasN8bcd/bRgdDwo7cr1e/z64GlU1xHARBBxsXH6JLJ2dpxsFP7W7tdxwHChpfyPNSdTdbjmNO99kJr7Xwdmw7yNWPMJUP+cGPuM8ZUGGMqmpubPcQFgkdta4+eWVurOxaVqDw72XUcABHmokk5SomP0evbGmXtyR/dAIbDS3mul1RywvfFkhpONcYYEyMpTVLb6V5rrT3+Z5Okl3SK6RzW2sestQuttQtzcnI8xAWCx8Nv7VGUMfrLKya7jgIgAsXFROnK6Xna39ajnQc7XccBwoKX8rxO0mRjTLkxJk7HbgBcftKY5ZLuCTy+TdLb9tg/cZdLujOwGke5pMmS1hpjko0xKZJkjEmWdLWkbef+doDgsaG2XS9tPKAvXFiu/LQE13EARKgF4zOUMy5er28/JJ+fq8/AuTpjeQ7MYb5f0gpJOyU9b63dboz5tjHmxsCwn0nKMsZUSnpA0oOB126X9LykHZJel/Q1a61PUp6kD4wxmyWtlfSKtfb1kX1rgDt+v9W3lm9Xbkq87r9ikus4ACJYdJTRNTPz1XKkTxX721zHAUKep307rbWvSnr1pGPfPOFxr6TbT/Haf5X0rycdq5Z03tmGBULF8xV12lLfoR/cMVfj4tkeF4Bb0wtSND4zSW/vbFJ336CS+b0EDBs7DAIjrOPogP5jxW4tHJ+hm+YWuo4DADLG6PrZBerqG9R/r6x0HQcIaZRnYIQ9/OYetff061s3zmQbbgBBoyQzSfNK0vXT9/eptrXHdRwgZFGegRG0u7FLT320X3cvLtWsIjZEARBcrp6Zr5goo397dafrKEDIojwDI8TaYzcJpiTE6O+uZhtuAMEnLTFWX71sol7f3qjVlS2u4wAhiTsGgDN4ek2tp3FbD3Tow+pW3XheoV7bxna4AILTly6eoGfX1enbL+/Qy39xkWKiuY4GnA3+xgAj4Gi/Ty9vblBBWoIWl2e6jgMAp5QQG61/vH66djV26Zl1da7jACGH8gyMgFe2HlR3/6BunV+sKG4SBBDkrp2Vr/MnZOr7b+zW4Z5+13GAkEJ5Bs7RnkNd2lDbrksm56gwPdF1HAA4I2OMvnnDTHUcHdAP3trrOg4QUijPwDnoHfDppY0HlJMSr8un5bqOAwCezShM1d1LSvXLD2u07UCH6zhAyKA8A+dgxfZGdR4d0K3zixXLTTcAQszXr5mmzOR4/f2LWzXo87uOA4QEPu2BYapuPqI1+9p04aRslWYmuY4DAGctLTFW//uTM7T1QId++eF+13GAkEB5Boahf9CvFzceUGZynK6anuc6DgAM2w1zCnTZ1Bz95xu71XD4qOs4QNCjPAPD8Nq2g2rr7tct84sUF8NfIwChyxijf7lplnzW6n8v3+46DhD0+NQHztKOhg6t2demiyZla0L2ONdxAOCclWQm6W+umqI3dxzS62zyBJwW5Rk4C4d7+vXbDQdUlJ6oq2cyXQNA+Pizi8o1LT9F31q+XV29A67jAEGL8gx45LdWz1fUy+e3umNRiWKi+OsDIHzERkfp32+ZrUNdvfre67tcxwGCFp/+gEcrdzepprVbN84tVPa4eNdxAGDEzSvN0BcuKNevPqrVu3uaXccBghLlGfBgf2u3/rCzSXNL0jWvJN11HAAYNd+4dqom5Y7T13+zma27gSFQnoEz6Okf1HPr6pSRHKcbzyuUMcZ1JAAYNQmx0frBHXPV1t2vf/zdNllrXUcCggrlGTgNn9/quXV16uod1B0LS5QQG+06EgCMullFafrrqybrlS0HtXxzg+s4QFChPAOn8X9W7NbepiO68bxClbCLIIAI8uVLJ2p+aboe+t02Nk8BTkB5Bk7h5S0N+vG7VVpclqlF5Zmu4wDAmIqJjtL3Pz1Xg36rr7+wWX4/0zcAifIMDGnnwU59/TdbtGB8hm44r8B1HABwoiw7Wf+0bIZWVbbq8Q+qXccBgkKM6wBAsDnc06/7nqpQamKMHv3MfL21s8l1JAAYMU+vqT2r8dZazSxM1Xdf26Xmrn6VZycPOe7uJaUjEQ8Ielx5Bk4w6PPrL57ZqEMdfXr0swuUm5rgOhIAOGWM0a3zi5WRFKdn19aqk90HEeEoz0CAtVYP/c82vb+3Rf9y80zNL81wHQkAgkJCbLQ+c/549Q769MzaWvmY/4wIRnkGAh55p1LPrK3T/ZdP0h2L+O9HADhRfmqCPjWvWPtbe7Rie6PrOIAzlGdA0m/X1+v/e2OPbplXpL+9eorrOAAQlOaWpOv8CVn6oLJFWw90uI4DOEF5RsT7YG+L/tdvt+jCSVn67q1z2EEQAE7j+tn5Ks1M0m831Kuxs9d1HGDMUZ4R0XYe7NSXf7Vek3LH6dHPLlBcDH8lAOB0YqKidNfiUiXEROnJ1TXqOMoNhIgsNAVErJqWbt3z87UaFx+jJ76wSKkJsa4jAUBISEuM1ec/UaajAz798sMa9Q74XEcCxgzlGRGpvr1Hd//0Iw36rX75xcUqSEt0HQkAQkpheqLuXlyqQ529emZtrQZ8fteRgDFBeUbEaezo1d0/XaMjfYN66ouLNSUvxXUkAAhJU/JSdPPcIu1tOqJ/fGmrrGUJO4Q/dhhERGnu6tPdj3+ktu5+/epLSzSzMM11JAAIaQvLMnX46ICer6hXUXqS/uqqya4jAaOK8oyI0d7dr8/9bI0OHu7VL7+4WHNL0l1HAoCwcOW0XGUkxenht/YoOT5aX7p4gutIwKihPCMiNHf16XM/W6Pqlm49ce8iLSrLdB0JAMKGMeb/tnfv8VFWdx7HP7+5JZmQO5CQCwTkJghFRNRS70rxyu6qq1Wr1nZtu7bdbbe166v7al233dVet9qLL6u1tdYq6rqltqVKdb0j4gWVSxC5RgkhCQnkNslkzv7xPGgIEwmayYTM9/16zet55syZ4TC/OTO/PM9zzuGmC2bR0R3n239cR8CMqz82Md3NEkkJJc8y4r3d3MHld7xAXUsnd145jwWTR6e7SSIiI044GODHlxxNIvEKNz6yloDBVQuUQMvIowGDMqJtbmjjop8/R0NrjHs+M58Tp4xJd5NEREascDDArZcezcdnlnLDH9Zy9/Nb0t0kkUGnI88yYtz7wrb97u9o6eCuZ7eQcI6rF0ykpq6VmrrWNLVORCQzhIMBbv3EXK6992W++fs1AFxxQnV6GyUyiHTkWUakrY1t3PH0ZgIG15w0ifJCzeMsIjJUIqEAP710LmccWco3f7+GHz5ao2nsZMRQ8iwjzurtzdz5zGaikSCfPekIxuZlp7tJIiIZJxIKcNvlc/n7eZXc8vhGvv7Qa1pIRUYEXbYhI4Zzjidq6lm+rp7qklwuP2480Sx9xEVE0iUUDHDzBbMpy8/mlsc3smtvjJ9eNpdoRN/NcvjSkWcZEWLxHh58qZbl6+o5uqqQqxdUK3EWERkGzIyvLJzGd/72KJ7csItP3L6ChtZYupsl8oEpeZbDXkNrjE/euZJXtjdzxpGlXHhMJaGgPtoiIsPJZcdN4LbLj2F93V4W/+RZXq9tSXeTRD4QZRhyWHtxSxPn3PI0q7c3c/G8Kk6bPhYzS3ezREQkiYUzy3jgcyfgnOOC255jyart6W6SyCHTeW05LDnn+MXTm7h5WQ1VRTnc9Y/zeXV7c7qbJSKSsfpOF/p+rlowkftf3MZ1D77Ggy/Vcu7scYQCyY/nXXrc+MFqosig0JFnOey0dHRzzW9e4j//tJ6FM0pZ+sWPMaM8P93NEhGRARqVFeKqj07kpCmjWbm5iV88tYnm9q50N0tkQHTkWQ4rKzY18tUHVlPX0sk3z53BpxZU6zINEZHDUDBgLDpqHBVFUR56uZZbHn+T8z9SwUcqC/S9LsOakmc5LLR3xfnushp+9dwWqkui3P/ZEzhmQlG6myUiIh/SrIoCyguyeeClWpas2s66HXtYPKdc09nJsKVPpgx7Kzc38bUHV7O1sZ2rPlrNdYum6UtVRGQEKRmVxTUnTeKpDbtYvm4nWxvbuGBuJVNK89LdNJEDKAORYaulvZsfLd/Ar5/fQlVRlPuuOZ7jJ5Wku1kiIpICATNOmTaWqaV5LFm1nbue28LRVYWcOaOUMXlZIeQykAAADkBJREFU6W6eyLuUPMuw05Nw/G7lNn7waA0tHd1ccfwErls0nVwteiIiMuKVF+Zw7amTeaKmnqffbOC07/8f/7JwKpcfP0Fz+MuwoGxEhpXn32rk3/+whvV1ezluYjHfOm+mZtIQEckw4WCAhTPKmFtVxItbm7jhD2u5f1Ut/7F4JvOqi9PdPMlwSp5lWHhl225ufXwjj6+vp6Iwh59dNpezjirTiGsRkQw2Oi+Lu6+ez7I36rjxkbVceNvznDZ9LF8+YyqzKgvS3TzJUEqeZcj1nkh/c0MbT9TUs7G+lZxwkIUzSlkweTTN7d38bqVWnhIRyXRmxlmzxnHytDHc9ewWbn9qE+f95BnOnFHKP58xhZnlwyuJPpTFYg6FFosZPpQ8y5DrSThq6vbw7FuNbG5oIzcrxKKZZRw3sZiscDDdzRMRkWEoGglx7amTueKECdz17BZ+8fQmzrllJ2ccWconT5jAiZNHEwjobKWknpJnGTLbm9pZsmo7v3puC3s74xTkhDln1jiOrS4mEtIgEBERObi87DBfOn0KV55QzZ3PbOKeF7axfN1OJpREuXT+eC6aV0VxbiTdzZQRTMmzpFRTWxfL1+7kkdd38PSbuzBgytg85s8pZmppHkEdJRARkQ+gIBrmKwunce1pk1n2Rh33rNjKf/15PT94bAOnThvDmTPKOG36WCXSMuiUPMug29HSwWNrd7LsjTpe2NxET8JRUZjDF0+bwsXHVvFkza50N1FEREaIrFCQxXMqWDyngvV1e7hv5XaWvVHHX9bsJGAwr7qYhTNKOba6mOnj8sgK6fJA+XCUPMuH1tTWxYpNjTz3VgPPvdXIpl1tABwxJpfPn3wEi44qY2Z5vmbOEBGRlJpels8N58/kW+fN4PW3W1i+diePrt3Jt/+4DoBw0DhyXD6zKgqYVVFAVXGU8sIcxhVkk52CMTcJ54j3OOKJBPGEI5FwBAJGwIygGYGANy1fQL+PhxUlz9KvZCOGu3sS1LV0Uru7ndrdHdQ2d7BrbwyASCjAxJJczjqqjKmleZTmZwPwWm0Lr9W2DGnbRUQkc5kZsysLmV1ZyFcWTuOd5g5e3d7M6tpmXtvewtJX3+G3fX7jRo+KMHpUFtFIkGgkRHY4SDQSJBQwepyjJ+FI+NvuHkcs3kOsO0EsnqArnvDuxxO0dHQTTzh6ehw9zg2ovZFggOxwgKxwkOxQgNysEHnZYfKyQ+Rlh8jPDlNTt5fxxVFyIjpynm5KniWpeE+ChtYY9Xs62bk3xs49ndTviVG/t5OE/10wKitEZVEOR1cVMmnMKCoKc3QNs4iIDDvlhTmUF+Zw9qxxACQSjtrdHbzd3ME7+24tnezaG6Ozu4f2rjiNbV10dMXpca7XkWJvGwkFyAoFiIQC5GWHyAoFyQp7ZVsb2wkFjFAwQChohAMBggEjFPSOOCecdwQ64fAT8QSd3V7i3dndQ2d3gub2brbv7qA9Fmdf+v2bFVsBKM3PYkJxLhNKov4t991tQU44Te9wZlHynME6u3vY1tTO1sZ2tja2sbWxnS2NbWxrauft3R3EE+/9xVwUDVOan820sjwqCnOoLMqhICesSzFEROSwEwgY40uijC+JDvprD+Y8zz0JR2sszp6Obo4sz2dbYxtbGtvZ1tjOkxt2Ue+f+d2nKBpm0phRTBqdy6QxozhijLcdXxzVrFaDaEDJs5ktAn4MBIE7nHM39Xk8C7gbOAZoBC52zm3xH7se+DTQA3zJOfeXgbymfDjOOZrautjR0sk7zR3etqWDupZOdjR3sq2pnbo9nfs9Jz87RPXoXGZVFHDe7HJ2tHQyNi+LsflZGmAhIiIyxIIBoyAnTEFOmPM/Un7A4+1d8f0Ogm1uaGfTrlaeqNnFAy/V7vc644ujflKdy/iSXCqLcqgszKGiKIdoRMdSD8VB3y0zCwI/Bc4EaoEXzWypc25tr2qfBnY75yab2SXAzcDFZjYDuASYCZQDy81sqv+cg72m9OGco62rh8bWGA2tXTS2xmhs673ton5vp5cgt3QSiyf2e344aJQVZDOuIIcFk0fvd8qnuiRKYXT/6XxStUqSiIjIQOm3qH/RSIjpZflML8s/4LE9nd1s2tXGpl2t3rbB2z6zseGA/KAkN+Il00VRKoq8s8vlBTkUj4pQHI1QlBshPzuks82+gfypMR/Y6JzbBGBm9wGLgd6J7mLgBn//QeAn5r3Di4H7nHMxYLOZbfRfjwG85rDinHd9UsJ5Awbcu/ve1iXeeyzhvPo9ztH17kCCBF09vfb33Xp66OhK0BaLszcWp7UzTmusm9ZYnL2dcVrfLYvT1NZ1wAd+n7ysECWjIpSMyuKoigIWzixjnJ8olxd625LciFZfEhERyQD52WHmVBUyp6pwv/JEwrGrNeYN+t83+N/fX7djD4+t20lXklwjFDCKcvcl02GKcyMURSPkZYfJCQfJiQTIiYS8ff9+tr+fHfYGXgYC5m3NCAa8W+/9YOC9WUiCfvlwTNgHkjxXANt73a8FjuuvjnMubmYtQIlfvqLPcyv8/YO95rBw3YOrWbKq9uAVB8morJB3y/a2edkhyvKzGZUVojg34iXIuVmU+KOCi3MjFOdGUjLFjoiIiIwsgYBRmp9NaX42x0woOuDxRMLR0BZjR3MnTe1d7G7rosm/7W73t23dbNjZSlNbF62xeNJke7BcPK+Kmy+cnbLX/yAGkjwnS/n7zr3SX53+ypNdtZ50Phczuwa4xr/bamY1/bQzVUYDDUP8b8r+FIP0UwzSTzFIP8Ug/TI2BpeluwHvGdIYfNe/pcGE/h4YSPJcC1T1ul8JvNNPnVozCwEFQNNBnnuw1wTAOXc7cPsA2pkSZrbKOTcvXf++KAbDgWKQfopB+ikG6acYpJ9ikPwIcF8vAlPMbKKZRfAGAC7tU2cpcKW/fyHwuHPO+eWXmFmWmU0EpgArB/iaIiIiIiLDykGPPPvXMH8B+AvetHK/dM6tMbMbgVXOuaXAncBv/AGBTXjJMH69JXgDAePAtc65HoBkrzn4/z0RERERkcFjboBLR2YqM7vGv3RE0kQxSD/FIP0Ug/RTDNJPMUg/xUDJs4iIiIjIgGmtRhERERGRAVLyjLeKopm9YmaP9Cm/1cxae93PMrP7zWyjmb1gZtVD3daRqm8MzPMdM9tgZuvM7Eu9ym/xY/Camc1Nb8tHjiQxON3MXjazV83sGTOb7JerH6SAmW0xs9f993uVX1ZsZo+Z2Zv+tsgvVz9IgX5i8D0zW++/zw+bWWGv+tf7Magxs4+nr+UjR7IY9Hrsq2bmzGy0f1/9IAX6i4GZfdH/rK8xs+/2Ks+4fqDk2fNPwLreBWY2DyjsU+/dZciBH+EtQy6Do28MrsKbznC6c+5I4D6//Cy8WVum4M3//fMhbONI1zcGPwcuc87NAe4F/s0vVz9InVOdc3N6TQP1r8BfnXNTgL/690H9IJX6xuAx4Cjn3GxgA3A9gJnNwBscPxNYBPzMzLRa1eDoGwPMrAo4E+i9Vrf6QersFwMzOxVvJejZzrmZwPf98ozsBxmfPJtZJXAOcEevsiDwPeC6PtUXA7/29x8ETjcbhutGHmaSxQD4PHCjcy4B4Jyr98sXA3c7zwqg0MzGDWmDR6B+YuCAfH+/gPfmYlc/GDq93+tfA3/Tq1z9YAg45x51zsX9uyvw1iUALwb3OedizrnNwEZgfjramCF+hPeb3HuglvrB0Pk8cJNzLgYH/CZnXD/I+OQZ+G+8Dtl7bckvAEudczv61N1vGXJg3zLk8uEki8ERwMVmtsrM/mxmU/zyZMvFVyAfVrIYfAb4k5nVAp8EbvLL1Q9SwwGPmtlL5q2sClC673vI3471y9UPUiNZDHq7Gvizv68YpMYBMTCz84G3nXOr+9RVDFIjWT+YCpzoX6r3pJkd65dnZAwGssLgiGVm5wL1zrmXzOwUv6wcuAg4JdlTkpRpupIPIVkMfFlAp3Nunpn9HfBL4EQUg0H3PjH4MnC2c+4FM/sa8EO8hFoxSI0Fzrl3zGws8JiZrX+fuopBahwQA+fcUwBm9g289Qp+69dVDFIjWT/4BrAwSV3FIDWSxSAEFAHHA8cCS8xsEhkag4xOnoEFwPlmdjaQjXeKeg0QAzb6Z6KjZrbRv76zv2XI5YM7IAZmdg/ee/2QX+dh4C5/fyDLxcuhSRaDP+Jdb/6CX+d+YJm/r36QAs65d/xtvZk9jHfqc6eZjXPO7fBPR+87Vap+kAL9xOApM7sSOBc43b03v6tikAJJYnAyMBFY7f8mVwIvm9l8FIOU6Kcf1AL/43/+V5pZAhhNhsYgoy/bcM5d75yrdM5V413w/rhzrsg5V+acq/bL2/3EGfpfhlw+oH5icDnwv8BpfrWT8QbqgBeDK/xR1scDLUkur5FDkCwGeNexFZjZVL/ambw3mFD9YJCZWa6Z5e3bxzvK9gb7v9dXAr/399UPBll/MTCzRcDXgfOdc+29nrIUuMS82Wcm4g1aWznU7R5J+onBi865sb1+k2uBuc65OtQPBt37fBe9+5vs/y5EgAYytB9k+pHnQ5V0GXJJiZuA35rZl4FWvMsFAP4EnI03KKEd+FR6mjeyOefiZvYPwEP+EYbdeNd7gvpBKpQCD/tH1kLAvc65ZWb2It7p0U/jzTJwkV9f/WDw9ReDjXiXkT3mP7bCOfc559waM1sCrMW7nONa51xPmto+UiSNwfvUVz8YfP31gwjwSzN7A+gCrvQPmmRkP9AKgyIiIiIiA5TRl22IiIiIiBwKJc8iIiIiIgOk5FlEREREZICUPIuIiIiIDJCSZxERERGRAVLyLCIiIiIyQEqeRUREREQGSMmziIiIiMgA/T+l5DuSz5+pywAAAABJRU5ErkJggg==\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "sns.distplot(np.random.binomial(1000, 0.5, 1000))" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "***\n", "\n", "What about an unfair coin?\n", "\n", "***" ] }, { "cell_type": "code", "execution_count": 7, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "array([2])" ] }, "execution_count": 7, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# Flip an unfair coin 10 times - how many heads?\n", "np.random.binomial(10, 0.2, 1)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "***\n", "\n", "Suppose we flip an unfair coin ($p = 0.3$) ten times, what is the probability that the flips are as follows?\n", "\n", "$$ HHTTHHHTTT $$\n", "\n", "***" ] }, { "cell_type": "code", "execution_count": 8, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "0.0004084101" ] }, "execution_count": 8, "metadata": {}, "output_type": "execute_result" } ], "source": [ "(0.3)*(0.3)*(1.0-0.3)*(1.0-0.3)*(0.3)*(0.3)*(0.3)*(1.0-0.3)*(1.0-0.3)*(1.0-0.3)" ] }, { "cell_type": "code", "execution_count": 9, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "0.0004084101" ] }, "execution_count": 9, "metadata": {}, "output_type": "execute_result" } ], "source": [ "0.3*0.3*0.7*0.7*.3*.3*.3*.7*.7*.7" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "***\n", "The probability of $r$ heads when flipping an unfair coin $n$ times is\n", "\n", "$$ P(r \\mid n , p) = {n \\choose r} p^r (1-p)^{(n-r)} $$\n", "\n", "***\n", "$$ {n \\choose r} = \\frac{n!}{r!(n-r)!} $$\n", "\n", "$$ {10 \\choose 3} = \\frac{10!}{3!(7)!} $$\n", "$$ = 120$$\n", "\n", "$$ \\therefore P = {120} (0.3^3)(0.7^7) $$\n", "$$ = 0.266$$\n", "\n", "***" ] }, { "cell_type": "code", "execution_count": 10, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "0.00040841009999999976" ] }, "execution_count": 10, "metadata": {}, "output_type": "execute_result" } ], "source": [ "(0.3)**5*(0.7)**5" ] }, { "cell_type": "code", "execution_count": 28, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "[0.028247524900000005,\n", " 0.12106082100000018,\n", " 0.2334744405,\n", " 0.26682793200000016,\n", " 0.20012094900000013,\n", " 0.10291934520000007,\n", " 0.03675690899999999,\n", " 0.009001692000000002,\n", " 0.0014467004999999982,\n", " 0.00013778100000000015,\n", " 5.904899999999995e-06]" ] }, "execution_count": 28, "metadata": {}, "output_type": "execute_result" } ], "source": [ "noflips = 10\n", "\n", "p = 0.3\n", "d = [ss.binom.pmf(i, noflips, p) for i in range(noflips+1)]\n", "d" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "***\n", "\n", "$ {n \\choose r} $ is the number of ways to select $r$ items from $n$, ignoring the order you select them in.\n", "\n", "***" ] }, { "cell_type": "code", "execution_count": 24, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "120.0" ] }, "execution_count": 24, "metadata": {}, "output_type": "execute_result" } ], "source": [ "import math\n", "\n", "n = 10\n", "r = 3\n", "\n", "choose = lambda x, y: math.factorial(x) / (math.factorial(y) * math.factorial(x-y))\n", "\n", "choose(n, r)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "***\n", "\n", "Note the following for ${n \\choose 0}$ and ${n \\choose n}$.\n", "\n", "***" ] }, { "cell_type": "code", "execution_count": 13, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "1.0" ] }, "execution_count": 13, "metadata": {}, "output_type": "execute_result" } ], "source": [ "choose(10, 0)" ] }, { "cell_type": "code", "execution_count": 14, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "1.0" ] }, "execution_count": 14, "metadata": {}, "output_type": "execute_result" } ], "source": [ "choose(n, n)" ] }, { "cell_type": "code", "execution_count": 15, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "120.0" ] }, "execution_count": 15, "metadata": {}, "output_type": "execute_result" } ], "source": [ "choose(10,3)" ] }, { "cell_type": "code", "execution_count": 26, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "70.0" ] }, "execution_count": 26, "metadata": {}, "output_type": "execute_result" } ], "source": [ "choose(8,4)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "***\n", "\n", "Even though the chances are, with $p = 0.3$ and $10$ flips, that there are three heads, the most probable outcome is all tails.\n", "\n", "**Clarification:** The one instance of tails beats a single instance of three heads, not all 120 possible options of three heads. \n", "\n", "***" ] }, { "cell_type": "code", "execution_count": 18, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "0.02824752489999998" ] }, "execution_count": 18, "metadata": {}, "output_type": "execute_result" } ], "source": [ "(1-0.3)**10" ] }, { "cell_type": "code", "execution_count": 1, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "0.02824752489999998" ] }, "execution_count": 1, "metadata": {}, "output_type": "execute_result" } ], "source": [ "0.7**10" ] }, { "cell_type": "code", "execution_count": 2, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "5.9048999999999975e-06" ] }, "execution_count": 2, "metadata": {}, "output_type": "execute_result" } ], "source": [ "0.3**10" ] }, { "cell_type": "code", "execution_count": 3, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "0.012106082099999993" ] }, "execution_count": 3, "metadata": {}, "output_type": "execute_result" } ], "source": [ "0.7**9*0.3**1" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "***\n", "\n", "What has all of this got to do with computers and bits?\n", "\n", "Would you consider the following a data set?\n", "\n", "***" ] }, { "cell_type": "code", "execution_count": 1, "metadata": {}, "outputs": [], "source": [ "import itertools\n", "#[\"\".join(seq) for seq in itertools.product(\"01\", repeat=10)]" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## End" ] } ], "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.3" } }, "nbformat": 4, "nbformat_minor": 4 }