{ "metadata": { "name": "", "signature": "sha256:ba59f7ab69bbecea2f5e21a8e5e7168936c30d63b344d80d46379ada525969a9" }, "nbformat": 3, "nbformat_minor": 0, "worksheets": [ { "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "**Scneario**: We toss a coin and get 8 heads and 2 tails. Our objective is to determine probability of getting head in the \n", " next toss. Use bayesian approach to estimate parameter i.e. probability of getting head ( $\\theta$ )\n", " " ] }, { "cell_type": "heading", "level": 2, "metadata": {}, "source": [ "Theoretical Approach:" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Assuming $\\theta$ has a beta distribution with $\\alpha = \\beta = 2.0$, expected value for $\\theta$ is given as\n", "\n", "$E[\\theta] = \\frac{n_H + \\alpha}{\\alpha + \\beta + N}$\n", "\n", "$E[\\theta] = \\frac{8 + 2}{2 + 2 + 10} = \\frac{10}{14} = 0.71$\n", "\n", "\n", "\n", "\n" ] }, { "cell_type": "heading", "level": 2, "metadata": {}, "source": [ "Using Metropolis Algorithm" ] }, { "cell_type": "code", "collapsed": false, "input": [ "%matplotlib inline\n", "import random\n", "import scipy.stats\n", "import scipy.special \n", "import matplotlib.pyplot as plt\n", " \n", "def integral(theta):\n", " \"\"\" This is likelihood X prior after removing any constant terms.\"\"\" \n", " return theta**9*(1-theta)**3\n", " \n", " \n", "n = 100 # Number of points to sample\n", "samples = [random.random()] # Random starting point\n", " \n", "for i in range(n):\n", " #Grab last sample point \n", " theta = samples[-1]\n", " \n", " #Create new sample by randomly selecting a point from a normal distribution \n", " #of mean = 0 and sd = 0.1. if the new sample is outside of the \n", " #domain then ignore it and use existing sample\n", " newTheta = theta + random.normalvariate(0, 0.1)\n", " if newTheta < 0 or newTheta > 1: \n", " newTheta = theta\n", " \n", " #If the probability of new sample as compared to last sample is less than uniform distribution\n", " #then ignore it. \n", " acceptanceRatio = integral(newTheta)/integral(theta)\n", " if acceptanceRatio > random.random(): # accept only if going uphill\n", " samples.append(newTheta)\n", " else:\n", " samples.append(theta)\n", " \n", "print \"Estimate: \", sum(samples)/n\n", " \n", "# Plot how sample theta varies with each iteration\n", "ylab = [i for i in xrange(len(samples))]\n", "pylab.plot(samples, ylab)\n", "pylab.title('Random Walk Visualization')\n", "pylab.xlabel('Theta Value')\n", "pylab.ylabel('Time')\n", "pylab.show() " ], "language": "python", "metadata": {}, "outputs": [ { "output_type": "stream", "stream": "stdout", "text": [ "Estimate: 0.649221876089\n" ] }, { "metadata": {}, "output_type": "display_data", "png": "iVBORw0KGgoAAAANSUhEUgAAAYgAAAEZCAYAAACNebLAAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJzt3XmYXFWd//H3N0tn7XQ2spINyEJCIBB2EBpkG0VUFBQX\nRNFnHsZtFB1hHDWz6I/lp+gouI0/wA1ZFEZQQRSaRZYYZAmEBIJkI0mHJJ21O+nt+/vj3Erdrq7u\nrl6qblX15/U857lL3ar6VnVyvnXPufccc3dEREQyDUg6ABERKU5KECIikpUShIiIZKUEISIiWSlB\niIhIVkoQIiKSlRKEJMLMlpjZz5KOozvM7BYz+89ovdrM1vfR637QzB7oi9fq5D0OfN9mNt3MdpuZ\n9fF7vMXMVvbla0qylCDkADNbY2b1UeWx2cx+Zmaj8vR2BbkBx8xWmdnFse1TzKw1y75dZtbV/wen\nm3Gb2YlmtsfMRmR57Fkz+yd3/4W7n9ud1+2BA3G7+zp3r/Re3gQVfY+HxF73MXef15vXlOKiBCFx\nDpzv7pXAUcBC4N+SDanXHgFOi22fBqzMsu8Jd2/N4fW69avb3Z8CNgDvbfMiZkcAhwO3def1eqFP\nzxYK8LpSBJQgJCt3rwX+CCxI7TOzq8xsdfRr+yUze1fsscvM7HEzu97MtpvZ383svNjjs8zskei5\nfwTGx9/PzC6IXrPOzB42s3mxx9aY2RfM7IXo7OYnZjbRzP5gZjvN7EEzG93BR3mUtsngVODajH1v\niY7DzO40s01mtiOKd34u35eZfSaKf0qWh28FLs3YdynwO3evi767x6LXMTO7wcxqo8/2QioGM6sx\ns8tj73ngedH2d8xsXfS8ZWZ2agexzox+/Q8ws5Oi7zRV9pnZ69Fxx5vZk9HfZKOZfdfMBkePPRq9\n3PPR8y7KbHYzs8OjmOvM7EUze0fssVvM7EYzuy/6N/FU/GxEioMShGQyADM7GDgPeDr22GrgVHcf\nBfw78HMzmxh7/HjCr/NxwHXAT2KP/RL4a/TYfwIfIWr2MLM50eOfISSO3wP3mtmg6LkOXAi8FZgL\nnA/8AbgKmED4d/yZDj7PY8ACMxsdNSEdC9wOjI7tO5koQQC/Aw4DDgL+Bvyi868LzOyrhAr/NHff\nmOWQnwOnRd8p0XteQkgcmc4hJKzZ7l4FXARsj30PnTULLSWc+Y0hfJ93mllFZ7G7+5NRc1Nl9Lyn\noucCNAOfJfzNTiJ8//8UPS+VYI+Mnn9n/HWjRHIvcD/hu/w08Ivob53yPmBJ9L6rga93FqsUnhKE\nxBlwj5ntAtYBrwH/lXrQ3e9y983R+h3Aq8AJseevdfefRG3bPwUmm9kEM5tOqJi/4u5N7v4YofJI\neR9wn7v/2d1bgP8LDCNU3Cnfdfc3owr4MeBJd3/e3fcDdwNHZ/tA7r42+iynESrPV919H/CX2L4K\nokTo7re4+153byIkwaPMrLKj78vMvgWcBZzh7ts6iGE9UAN8ONr1VmAIIRllagIqgcPNbIC7r0p9\n512J+jLq3L3V3b8VvcfcXJ4b+S6wy92/HL3e39x9afR6a4EfAafn+FonAiPc/Rp3b3b3h4H7CIkx\n5Tfuviz6m/8CWNSNWKUAlCAkzoF3RmcI1cCZhIodADO7NOpYrTOzOuAIwq/LlAMVmbvXR6sjgSlA\nnbs3xI5dG1ufQqjEU891YD0wNXZMbWy9IWN7X/Q+HUk1Mx1oSgIej+172t2bzGygmV0TNaPtBF6P\njh3f7hWD0cDHgWvcfXcn7w/hbCGVID4M3BZVjG24+0PA94AbgVoz+2EnCaqNqBluRdQ8VgdUdRJ7\n5nP/kfB9fCC2b07UBLQp+j6+Ttu/d2emEP6GcWuj/RD+rWX+TTv7G0oClCAkK3d/lPCL8loAM5tB\n+AX5SWCsu48BXiS3TspNwBgzGx7bNyO2/kZ828wMmBbt70h3OkfjCSLVZv8Y7ZPGB4ALgLdGzTuz\nsrxXvImnjtDcdbOZxc92srkbONjMzgDeTfbmpfAG7t9192OB+cAc4IvRQ3uB+NVQk1IrZvaW6LiL\n3H109PfZSQ7fU/Tc/yD8ONgTe+j7wArgsOj7+DK51xkbgWnR3zJlBp3/TaXIKEFIZ74NHG9mJxAq\nJge2AgPM7KOEM4guRc0Ty4B/N7PBUefp+bFD7gTebmZnRm3XVxLOCp7oo8/xKHAMISH8Jdq3HDgE\nOIN0ghgJ7Ae2W7gs9RsZr2NkVLhRIv0g8BszO66jANx9L3AXcDOwxt3/lu04MzvWzE6Ivod6wveQ\nOtN4DrjQzIaZ2WHA5aQTViWhz2CrmVVE/SJdXqJsZtOAO4APu/vqjIdHAruBegsXDVyR8XgtcGgH\nL/10FP+/RH/zasLf/Fept+4qNkmeEoR0yN23En7pfsndVwDfBJ4kNCUdQWimOXA47TtQ49sfIPRX\nbAe+SuwXtLuvAj5EOGN5E3g78A53b+4svC7eO/45XgW2AJvcfVe0zwmVWCXpRPRTQjPIG4Szoye7\neB+PXutPwMcIHeudtaPfCkyP3ifzs6RedxThTG07sIaQkK+PHrsBaCRUzDcTOr9T7o/KK9HzGog1\n23UUO6E/ZALw69iVTMujx75A+LvtimL6VcZrLAFujZoc3xt/D3dvBN4B/APhb/o9QhJ6pYN4yLIt\nCbN8TRhkZv+P8B99i7svjPaNJVxBMoPwj/hid98RPXY14T9ZC/AZd/9jXgITEZGc5PMM4mbCZZJx\nVwEPuvsc4M/RNtF13u8jtLmeB9xkXd/VKiIieZS3Sji6lLEuY/cFpJsWbgVSN1q9k3BVR5O7ryFc\nE318vmITEZGuFfpX+sToDl0I7aipm6ymEIYjSNlA20scRUSkwBJrxok6CTvrAFGHlYhIggZ1fUif\nqjWzSe6+2cwmE64sgXDVyLTYcQeT5XppM1PSEBHpAXfv9qXFhT6D+C1hDB6i5T2x/e+Prt+eBcwm\njCvTjrsXffna176WeAyKU3EqzvaludnZvNn561+dX//aueEG53Ofc97zHue445yJE1MNG23L6NHO\nYYc5J57onH++s3JlaX2XPZW3Mwgzu40wbsv4aITHrwLXAHdEI1KuAS4GcPcVZnYH4a7NZuCfvDef\nSkRKUnMz1NeHUlcHy5ent+vrYe/e3m3X14f3GDEChg8PJbU+YgQsXtz2scxjUutTso3ZW4byliDc\n/ZIOHjqrg+O/Qfs7V0WkyO3bB7W1sGVLWNbWhsq9JxV4vPLevx8eeqjjijpVxo6Fgw/u+JjM7YoK\n6Nu59MpXofsg+oXq6uqkQ8iJ4uxb5RKnO+zZk67s4xV/tn379sGECTBxYrqMGRMq5t5U3jU11RT7\nV1oqf/Oeytud1PlgZmp5EumF+npYurTrit+sbYU/cWL7JJAqVVX6RV7szAzvQSe1EoRIP3LjjXDN\nNXDCCdkr+1QSGKmBt8tKTxOEmphE+pHmZrjwQvjOd5KOREqBEoRIP1JVBT/4Adx1F4wa1fNSWRn6\nC6S8qYlJpB9xh927Q9m1q3dl0KD2SaMnyWbYMPVh5Jv6IESkYNzD1Uu9TTK7dkFTU/ukcc01cOqp\nSX/K8qE+CBEpGLPwy3/YsNCp3RuNjW3PaJYsCTfIKUEkTwlCRBJVUQHjxoUCMGlS58dL4WhSHhER\nyUoJQkREslKCEBGRrJQgREQkKyUIERHJSglCRESy0mWuItJtqRvlGhralmz7spXOjnvllTBxjyRP\nCUKkxLmHm836spLu6tj9+8P9C6mb5TorQ4e23R4xAsaP7/y4hQuT/lYFlCBECq65GXbsCLOubd+e\nvezdm3uFvm8fDBzY/Yo6VaqqcjsuXoYMgQFqoC57GotJpIf27QuVeWcVfbbH9uwJlfLYsdnLmDFh\nPoZcKupUhT5wYNLfhhQzDdYn0kdeeAGefLLzSn779nAmMG5cqNA7quyzVf5VVfr1LYWlwfpE+siy\nZfDtb8PKlWF79mx497vhrLPCjGupyn74cA1TLeVNZxAiHaithUcfhUceCcu1a+Gkk+D00+G00+DE\nE9W0I6VBTUwiebZtGzz2WEgYt98O110HH/pQ0lGJdE1NTCJ5Nm4cvOtdoezaFS71FClnShAiPTBk\nCFx7bTijWLAgXaZPV7+ElA81MYn0QEMDPPssvPQSrFgRli+9FM4s5s8PJZ44pk1T4pDkqA9CpAjU\n1aUTRjxx7NmTThyHHw7z5oUyaxYM0nm85JkShEgR2749JIwVK2DVqnAJ7cqVsHEjHHJISBZz56YT\nx9y54X4Jkb6gBCFSghoaYPXqdMJ4+eWwXLUqJIhUwpgzByZODGMYjR8fOszHjw/3Yoh0RQlCpIy0\ntsKGDenE8cor8Oab4VLbrVvTZcCA9kkjs2TuHzo06U8nhaYEIdLPuIdB/TKTxtat2fel9g8enD2Z\nZO4bNy6U1F3jUrqUIESkS+6hw7yrZJLat317WDdLDzGSShpdLceO1dlKsVCCEJG8cA99Jdu2pQcq\nTK13ti91ttJVEsm2r6Ii6U9dXpQgRKSopJrAukok2fadeGK4CVH6hhKEiJSF1avh3HPhtdeSjqR8\n9DRBaFR6ESkqAwfCli3w1a/C//5vuJpLvwuTkcgZhJldDXwIaAWWAx8FRgC3AzOANcDF7r4j43k6\ngxApc62t8LvfwdKl8MwzYX4OM1i8GI49NiwXL4apUzV8Sa5KponJzGYCDwGHu/t+M7sd+D2wANjq\n7teZ2ZeAMe5+VcZzlSBE+hn3cBbxzDNtC4REceqpcPXVShadKaUEMRZ4EjgR2A3cDfw38F3gdHev\nNbNJQI27z8t4rhKEiOAOb7wRzi4uuigMkjhsWNJRFa+S6YNw9+3AN4F1wEZgh7s/CEx099rosFpg\nYqFjE5HSYBZu5jviCM3vnU8FH0fSzA4F/hmYCewE7jSzNvNyububWdZThSVLlhxYr66uprq6Ol+h\nikhC3MPIuOvWhalesy23bw/9EGedpfsmMtXU1FBTU9Pr10miiel9wNnu/vFo+8OE5qYzgTPcfbOZ\nTQYeVhOTSHlqaQlNRJ0lgIEDwwRMM2ZkX06apDnBc1VKfRBHAb8AjgP2AbcASwlXL21z92vN7Cpg\ntDqpRcrTpZfC/ffDYYd1nAA03HnfKZk5qd39eTP7KbCMcJnr34AfAZXAHWZ2OdFlroWOTUQKY9cu\n+NGPwvzeUrx0J7WI5F1rK7z6avoS1V/+En78Yzj//KQj6x9KpompN5QgRIpfPBksWxaWzz4bBuFL\n3eS2eDGceWYYzE/yTwlCRBJ1/fVw330hGYwb1zYZHHNMuCxVkqEEISKJOvhg+OY3w2Wn48YlHY3E\nlcyNciJSvk45RcmhnChBiIhIVkoQIiKSlRKEiIhkpQQhIiJZ6SomEekThx4axkaqqoLhw9NlxIjs\n67luDxmiuR56S5e5ikiitm6F2lqor4e9e8MyVbqznflYY2P3k0t3E1G5JyElCBEpS83N0NDQs+SS\n63Zzc9fJZOZMuOwyWLAg6W+k+0pmsD4Rke4YNAgqK0PJl+bmtokjWzJ57jk455ww2uwnPgEXXxyS\nRznTGYSISI6am+F3vwsDDT7xBLz//VBdHZLGjBkwYUJxznCnJiYRkQJavx5uvTWMPbV2bSi7d8O0\naW3ntkiV6dPDY0nMfqcEISKSsPr69rPjxdc3bQqDFmZLHgsWwKxZ+YlLCUJEpIi5h2lWH3oI/vhH\nePBB2LIlPFZREZqqHnggP++tBCEiUgQaG2H1ali5Ml1WrQrLIUNg3rxQ5s5Nr8+cmd/5tZUgREQK\noLERNmyANWtCs1HmctOm0GyUqvxTyWDu3ORGulWCEBHpA/v2hT6DjhLAli0weXL41T9jRvtlUh3R\nnVGCEBHpoXvvha9/PSSAurow+VFHCWDq1HBvRinRjXIiIj20eHG4Ce6222DYMHjf+8I9DgsXJh1Z\nsnQGISIScQ/3NfzqV6FUVsKnPgVXXJF0ZL2jMwgRkV4yg0WLQj/Diy/CU0/Bzp1JR5UcJQgREcLV\nSd//Ptx4Yxhj6dOfhrvuCgP19VdKECIiwH33wQ9/CDffDCefXN7Df+dKCUJEBGhqCp3Sp5ySdCTF\nowjHHRQRkWKgBCEiIlmpiUlE+oX6+nB1UkflxRdh/vykoywuug9CREpSczNs29Z5pR8vTU1hQp+J\nE8MyWznqKJg0KelP1vc01IaIlKWdO2HJkjAIXm1tusKvq4OxYzuu7DNLZWX/vTJJN8qJSFlKjXv0\nyCPhjOGSS+Aznwk3tOVziGxRJ7WIFLkRI+CGG8IQ2/ffD4MHh3GTzjgjzA2tRoX8UROTiJSc/fvD\n7GuXXw6PPRbmXJCO9bSJSWcQIlJyhgyBCy4IfQstLUlHU74SSRBmNtrM7jKzl81shZmdYGZjzexB\nM3vFzP5oZqOTiE1ERIKkziC+A/ze3Q8HjgRWAlcBD7r7HODP0baIiCSk4H0QZlYFPOvuh2TsXwmc\n7u61ZjYJqHH3eRnHqA9CRA5YtChM/3nWWWEMpWOOKb7pPotBydwHYWaLgB8CK4CjgGeAfwY2uPuY\n6BgDtqe2Y89VghCRAzZtgocfhr/8JZTVq0OSOOWUUE46CcaNSzrK5JVSgjgWeBI42d3/ambfBnYD\nn4onBDPb7u5jM56rBCEiHdq1K0zy88QTIWE8/XSYQ/rUU+H006G6Osw33d+U0o1yGwhnC3+Ntu8C\nrgY2m9kkd99sZpOBLdmevGTJkgPr1dXVVFdX5zdaESkZo0aFeyTOOSdsNzfD8uXw+ONwzz3wuc/B\n6NEhUaTK1KkJBpwnNTU11NTU9Pp1ErkPwsweBT7u7q+Y2RIgNWfTNne/1syuAka7+1UZz9MZhIj0\nSEsL7N4dziruvTeUdevgsMNCorjkEjjzzKSjzI+SaWICMLOjgP8BKoDXgI8CA4E7gOnAGuBid9+R\n8TwlCJF+orExVOh79oRlfL0n+/btC3dlV1bCyJFhGV8/7TS44oqkP3V+lFSC6CklCJHi19AAq1b1\nvkJvbU1X4vGKPHOZ677hw2FAP701uJT6IESkjN19N/zrv4ZKfteu0A/QkbFjYcqUUObNC/0Bqe2D\nDgp9CpWV6eXgwYX7HKIzCBHJI/cwbtKuXemE0dEyl2MGD26bMOLLbPs6Ww7qRz+P1cQkImXNPfQj\ndDepdLQcMgQ++Um4/vqkP1n+qYlJRMqaGQwbFsrEib17LXf4+c/hvvv6JrZy1U+7bESkPzPTkBy5\nUIIQEZGs1MQkIv3Grl2wfn2Yne6RR5KOpvgpQYhIWaivD5V/ZtmwIb3e0gLTpoXxmKZNg0svTTrq\n4tblVUzR0NtfB6a6+3lmNh84yd1/UogAM2LRVUwi/ZQ7PPkkvPZa9kTQ0BAq/lTln62MHh36H/qb\nvF3mamb3AzcDX3b3I81sMGE+hyN6FmrPKUGI9F9794Z5H5YvD+uQngvi7LPDiK0zZvTfu6U7k88E\nsczdjzWzZ9396Gjfc+6+qIex9pgShIi4wxtvwMqVoaxalV7ftg3mzAl3Zc+dG5bz5oV9I0YkHXly\n8nkfxB4zOzDlhpmdCOzs7huJiPQFs3RT0llntX1s92545ZV00rj77rBcvRoWLoR3vCOUI4/sn01N\n3ZXLGcRi4LvAAuAl4CDgve7+fP7DaxeLziBEpNsaG+Gxx9LDfDc3p5NFdXW4q7qc5XWojajfYQ5g\nwCp3b+p+iL2nBCEi3dXaGpqeamtD2bwZamrgtttCX8bIkXDDDfDxjycdaf7ksw9iEPB2YCbpJil3\n92919816SwlCRCBcrrp1a7rSzyybN6fXt24NA/RNnNi+TJoUlscdBxMmJP2p8ieffRD3Ag3AcqC1\nu28gItJTLS2wZAmsWdM2AWzbBlVV7Sv6iRPh8MPbJoEJEzSsRk/lcgbxgrsfWaB4OqUzCJH+Zfv2\n0Bl9001tE8FBB2luiO7I5xnEH83sXHd/oAdxiYi009oa5onILI2Nbbe3boWhQ+Gyy5KOuH/KJUE8\nAdxtZgOAVOe0u/uo/IUlIn2lpaXzSjiXirqvj2tpCc0+Q4Z0XS65JOlvsP/KpYlpDXAB8KK7J9oH\noSYmKRUNDbBzZ7KVcKq0tmaveHOtoLtzbK7HDRqk+xAKKZ9NTOuAl5JODiKl5KKL4M9/DjOg9VRV\nVRg7KLNMnpxeHzOm7WOjRqkylr6TyxnErcAs4A9AY7Rbl7mK5KC1NSSJhoa2pb6+/b7u7s+2r7Ex\ntNmnZl5LleHD2+/raH+u+4YNCwlIyaf45fMM4vWoVETFANXSIjkYMCBUrsOHF+b9WluzJ5hcks+b\nb3Y/UTU3t01II0bAnXeGYS2k9HWZINx9SQHiEJE+MGBAqKQLNTBdS0s4Q0oljg99KAykpwRRHjpM\nEGb2PXf/lJndm+Vhd/cL8hiXiBS5hoaQDDZsSJf165OOSvpSh30QZrbb3SvNrDrLw+7uBZ+wT30Q\nIoWxe3f7yj+z7NkDU6eGEp+o5yMfCZ3lUjz6fCym+PwPxUIJQqTvNDXBL34Br7/evvJvbk5PzRlP\nAPEyfrw6qEtFPjqpDzKzzxM6pTMlchWTiPSd/fvhoYfCXAnr1sHGjWH/ggVhvoQZM0KZPj29HDky\n2ZilsDo7g9gE/KCjJ7r7v+crqI7oDEIkfxobQ7PS2rUhYcSXqfXhw9snjVQimTwZxo4NVzPpzKK4\nqIlJRPLKPVwKmy15rF2bHmUVQqLobqmsVGLJFyUIESkKDQ1hFNbuln37epZYqqrC5b3SsXwkiHHu\nvq3XkfUhJQiR8tXYCHV13U8su3eHIUcefxzmzk36UxSnPu+kLrbkICLlraIiPclPdzQ3w6mnhqHB\nlSD6Vi5DbYiIFI0dO9r2g6xbFy7VVf9F3+tysL5ioiYmkfLW0gKbNrXvCI8vW1vbX0k1cyZceGEY\nPFDa6/M+iGKkBCFSetzDXdfbt4ernLZtS69v2tQ2AWzcGDqes11Km1qOHq2zhe4quQRhZgOBZcAG\nd3+HmY0FbgdmAGuAi919R8ZzlCBEErRvX/aKvqP11HLwYBg3LpSxY9PrEye2vSFv2jSdBeRDKSaI\nzwOLgUp3v8DMrgO2uvt1ZvYlYIy7X5XxHCUIkT7kDsuWwebNuVX0TU3ZK/ps6/Hl0KFJf9L+raQS\nhJkdDNwCfB34fHQGsRI43d1rzWwSUOPu8zKepwQh0od27oSTToKXXw7bY8fCe94DRx+dvdIfMULN\nO6Wo1BLEncA3gFHAF6IEUefuY6LHDdie2o49TwlCJA/WrYM//QkefDBMlTp6NJx9Npx1FpxxRtiW\n0tXTBFHw+w/N7Hxgi7s/S/aBAImygDKBSIFMnw4f+xjcdltobrrjDpg1C37wgzCa69q1SUcoSUji\nPoiTgQvM7G3AUGCUmf0MqDWzSe6+2cwmA1uyPXnJkiUH1qurq6murs5/xCL9yIABsGhRKF/4Ahxx\nRLhbWUpHTU0NNTU1vX6dRC9zNbPTSTcxXQdsc/drzewqYLQ6qUWSd+SRMGkSzJ4drjqaNCmU1PrE\nieqELnYl1Qdx4M1DgrgyuoppLHAHMB1d5ipSNF5+GV56KTQ9bd4cRm3NXA4b1jZpdLQ+YUK45FUK\nqyQTRHcpQYgUH/cw/EU8aWRLJJs3h/GSRo1qn0CmTGk7e92UKWFsJukbShAiUvRaW8P9FPGksXlz\nuIN6w4b0PNibNoVLazub8nTq1DCBkXRNCUJEykZLC2zZ0n6u7FQCSZXhw9snjfj2tGnhjKW/U4IQ\nkX7FPZyNZEsi69eH9fXrYeDAkCjiSSNzvaoq6U+TX0oQIiIZUv0jqeQRTxzxfQMGdJ5AZswIU6KW\nKiUIEZEs3GHv3pAospW6unAj4PLl8OKLUF/f/jUWLoQXXih87H1FCUJEypJ7mOe6owo+Vcl39nhF\nRZiWdPTo7KWzx6qqYFCJT62mBCEiZWPZMvjgB9MVfGNj9uPGjw/NQNOnp8u0aaEcdFC6gu/vl8wq\nQYhI2Whpgeeeg127wjAfqWV8vatlRUW4gqmysufLVCn1m/uUIEREIu6hL6GnyaU7yWbyZLj++tDR\nXax6miBKvGVNRKQ9szB3xYgR4U7trriHfoz4vRbxey5efRVefz0kDLPQbFVVFa5wmj27fOfI0BmE\niJQl91Ch19WFGfFSy9ratpV/ar2iIvsNd/H1MWNKMxmoiUlEyk7qCqZUBZ9Z2Wfui6/v2BHutB4z\nJgzbMWZMKBMmhE7szCQwcmTSnzZ/lCBEpKg1N8Pq1V1X8pn7BgxoW8FnW8+2b/To0u9c7itKECJS\n1G66Cb7yFZgzJ/fKfsyYMJS49I46qUWkqO3fD5deCjfckHQkkislCBEpiKFD4Wc/C0NapG5uS93U\nllov536AUqQmJhEpiJYWWLUqDI6XKuvWtV0fNqx90oivT52qu6J7Qn0QIlLSUsN3ZyaOv/8dli4N\n+8zCfQ3TpsERR8BPfpJ01KVBfRAiUnL2728/CVBm2bo1XJp6wgltJwOaNy/p6MufziBEJO/eeANu\nvbV95b9zZxiqItuUoqkyaVLpj6aaNJ1BiEjR+v3v4a674PLL4Zxz0pX/hAnFPYZRf6cEISJ51dgY\nbnpbvBg++cmko5HuUIIQkV7ZtSvMyLZuXVjG19etgy1bQjPSF7+YdKTSXeqDEJEOtbaGwe2yVfyp\nfU1NYc7mGTPC5ajx5YwZMGWK+hCSpstcRaRP1deHyn3IkOwVf2p97NjSHOG0P1EntYj0qf37Q8Vf\nW5t0JJIUJQgRyWrAANizJ4yfNGVK6EfILCNGJB2l5JOamESkQ488EmZS27QJNm4My1TZuDE0P8UT\nRmYiSW1XVqoZKknqgxCRgnIPk/LEE0ZH69A+cUycCOPGwfjxbZfjxmkeh76mBCEiRckddu9uf/ax\nZUsYe2nr1rbL7dtD01W25BFPIpn7hg5N+pMWLyUIESkLra1hCI5syaOzfRUVXSeRww6D445L+hMW\nnq5iEpGyYBaamIYPD9OGmoXKf+TIMMPcQQeFm/NSZefOUDZuDCO/rluXfq0RI2DUKKiqCsu5c+Gn\nP03us5VDzygRAAAOV0lEQVQanUGISJ+prw/9EvHKO7My72x7167QHDVkSLpST5XubFdVhYSiG/QC\nNTGJSKIaG8ONc5s353b88OGhw3rSpHSZPDl0XldWhsc7KsOGhSSiK6NyowQhIkWhpQUaGsLZRGrZ\nWenpMU1NnSeQrhJMZ4/Hj6moKP1EpAQhIv1Kc3P7xJGPhNTS0j55TJ0KZ58N554LCxcWfwIpmQRh\nZtOAnwITAAd+5O7/bWZjgduBGcAa4GJ335HxXCUIESmopqb2SeW11+CBB0JpaAhzXJx7bkga48cn\nHXF7pZQgJgGT3P05MxsJPAO8C/gosNXdrzOzLwFj3P2qjOcqQYhIQTU2wt697cuePWH5wgvwm9/A\nSy+FM4nFi+HKK+H970868rSSSRDtAjC7B/heVE5399ooidS4+7yMY5UgRKSdxsZ0hd1Z6ckxEC6X\nHTkyLDsq8cdPPhlOOSXZ7ySuJBOEmc0EHgGOANa5+5hovwHbU9ux45UgRMrQyy/D00/3vGKH7lXg\n3Xm8oiLZ76YvlNyNclHz0q+Bz7r7bov18ri7m1nWTLBkyZID69XV1VRXV+c3UBHJu2XL4M47w9Sk\ndXXhXoq6ujDkeFxlZZjLeubM9LzWqbmtR48ON9KlysiRxd95nC81NTXU1NT0+nUSOYMws8HAfcAf\n3P3b0b6VQLW7bzazycDDamIS6d8aGtonjXjJti+1f//+kDQyE0eqZNs/bly4l6PclEwTU9R8dCuw\nzd0/F9t/XbTvWjO7ChitTmoR6Yw77NuX/a7srVvT06KmyoYN4bLVbIYNC0njkUdg9uzCfo58K6UE\ncSrwKPAC4TJXgKuBpcAdwHR0matIWXMPv/C7Gn4jl8cGDmw/zEZnQ3Bke6yysjz6GjpSMgmiN5Qg\nRErb1VfDj38cKvampraPDRgQ+hbmzIFZs8Kv+a4q+srKMOSGdE4JQkSK3u7dsHp1GK+ptrbj5e7d\nYdTWiRPDGE2Zy/h6asRX6ZgShIiUjcbGMKFQKmF0lkzq60OySCWMK66At70t6U9QXEruMlcRkY5U\nVKQvYe3Myy/Dz38emq3eeANOPLH8OpiTpAQhIiWlsRGuuQZuvz1c0nrRRWGoi5NPDv0Y0neUIESk\npDQ2wuuvh/UdO+CJJ8KlritWhHGQjjhCHdd9RX0QIlKy9u6F55+HZ55Jl9deg3nzQrJYvBhOOgmO\nOirpSJOlTmoRKWmtraHC3727/RSk2aYl7Wi7ri68VmVluBR28mT4y1/C/RL9lRKEiCQidcNbVxV3\nV5X9nj3hbubU/Q3xeyC62s7cN2yYLn2NU4IQkS41NYWKOLPs3t35dmf7Wlvb3rjWWcXd2fbIkTBI\nvaJ5oQQhUmaam7tXYedyTHNzqIjjpbKy8+2u9g0Zol/rxU4JQqQAWlvDFTMNDellZsm2v7N99fXZ\nK/Wmpr6pwONl6FBV5v2REoT0O/HKuruVc0+PbWwMv5iHDg3t3PGSbV9H+zP3pSrzeKWuylz6ihKE\nlKQbbwx3wPakIk9V1r2tnLtz7JAhuhlLSo8ShJQkM/iP/4Dhw7tfYauyFsmNEoSUJDN4800YPz7p\nSETKlwbrk5L0trfBoYeGyxyPPhoWLUovZ85UG7xIknQGIYlrbYU1a+DZZ+G559LLPXvaJoxFi2D+\nfBg8OOmIRUqLmpik7GzZEsbZSSWMpUvDTGMPPph0ZCKlRQlCyt7SpfCpT4WliOSupwlC14CIiEhW\nShAiIpKVEoSUlMbGpCMQ6T+UIKRkzJ8f7qD+n/9JOhKR/kGd1FJUWlraDoaXuf7883DllfDww2EO\nYhHpmm6Ukz7jnq6QO6usu1rvyfNaWtoOq5EaWiO+fu654RJYEckvnUGUgD17YNu2wlXWmSOWdlRR\n52N90CDdPS3S13QfRBlbuBC2bw9DQeergo7v0wQwIuVFTUxlrKEBHnkEDjss6UhEpD/RVUwiIpKV\nEoSIiGSlJqYi1NAAq1fDq6/CK6/A1q1JRyQi/ZESREIaG+H110MCePXVdDJ49dUwgc6sWTBnDsye\nDd/6VtgWESkkXcVUAC0t4e7fFSvSSWDDBjj44HQSmD07vT59OgwcmHTUIlIudJlrEVu7NlyqumRJ\nOgnMmgUVFUlHJiL9gS5zLWLuMGYMfP7zSUciIpK7orqKyczOM7OVZvaqmX0p6Xh666WX4Oqr4fTT\nw7SZIiKlpGgShJkNBL4HnAfMBy4xs8OTjar73ngDrriihkWLwphBLS1w771wzz1JR9ZeTU1N0iHk\nRHH2LcXZd0ohxt4omgQBHA+sdvc17t4E/Ap4Z8IxdcvDD8OCBfD00zXccEPoe7juOjjyyKQjy65U\n/nErzr6lOPtOKcTYG8XUBzEVWB/b3gCckFAsPbJ+PVxwARxyCJxxRtLRiIj0TjGdQZTe5UkiImWs\naC5zNbMTgSXufl60fTXQ6u7Xxo4pjmBFREpMSd8HYWaDgFXAW4GNwFLgEnd/OdHARET6qaLpg3D3\nZjP7FPAAMBD4iZKDiEhyiuYMQkREiksxdVIf0NUNc2Y2z8yeNLN9ZnZlEjFGcXQV5wfN7Hkze8HM\n/mJmiVzwmkOc74zifNbMnjGzM4sxzthxx5lZs5ldWMj4Yu/f1fdZbWY7o+/zWTP7t2KLMRbns2b2\nopnVFDjEVAxdfZdfiH2Py6O/++gijHO8md1vZs9F3+dlhY4xiqOrOMeY2d3R//enzWxBpy/o7kVV\nCM1Lq4GZwGDgOeDwjGMOAo4F/gu4sojjPAmoitbPA54q0jhHxNYXEu5HKbo4Y8c9BNwHvKcY4wSq\ngd8WOrZuxjgaeAk4ONoeX4xxZhx/PvCnYowTWAL8n9R3CWwDBhVhnNcDX4nW53b1fRbjGUSXN8y5\n+5vuvgxoSiLASC5xPunuO6PNp4GDCxwj5Bbn3tjmSCCJGShyvVHy08BdwJuFDC4m1ziTnNU7lxg/\nAPza3TcAuHsx/81TPgDcVpDI2solzk3AqGh9FLDN3ZsLGCPkFufhwMMA7r4KmGlmB3X0gsWYILLd\nMDc1oVg60904Lwd+n9eIssspTjN7l5m9DPwB+EyBYovrMk4zm0r4B//9aFcSHWi5fJ8OnBydxv/e\nzOYXLLoglxhnA2PN7GEzW2ZmHy5YdGk5/x8ys+HAucCvCxBXplzi/DGwwMw2As8Dny1QbHG5xPk8\ncCGAmR0PzKCTH65FcxVTTKn0muccp5mdAXwMOCV/4XQopzjd/R7gHjN7C/AzwulnIeUS57eBq9zd\nzcxI5ld6LnH+DZjm7vVm9g/APcCc/IbVRi4xDgaOIVxWPhx40syecvdX8xpZW935v/4O4HF335Gv\nYDqRS5z/Cjzn7tVmdijwoJkd5e678xxbXC5xXgN8x8yeBZYDzwItHR1cjAniDWBabHsaIRMWm5zi\njDqmfwyc5+51BYotrlvfp7s/ZmaDzGycu2/Le3RpucS5GPhVyA2MB/7BzJrc/beFCRHIIc54peDu\nfzCzm8xsrLtvL5YYCb80t7p7A9BgZo8CRwGFTBDd+bf5fpJpXoLc4jwZ+DqAu79mZq8TfmQtK0iE\nQa7/Nj+W2o7i/HuHr1joDp8cOloGAa8ROloq6KTjitAxlFQndZdxAtMJnUYnFvP3CRxK+pLnY4DX\nijHOjONvBi4sxjiBibHv83hgTRHGOA/4E6Fjczjh1+T8YoszOq6K0Ok7rNB/7258n98Cvhb7+28A\nxhZhnFVARbT+CeCWzl6z6M4gvIMb5szsH6PHf2hmk4C/EjqDWs3ss4R/3HuKKU7gq8AY4PvRr94m\ndz++UDF2I873AJeaWROwh/BrraByjDNxOcb5XuAKM2sG6inw95lLjO6+0szuB14AWoEfu/uKYosz\nOvRdwAMeznYKLsc4vwHcbGbPE/p2/8ULd8bYnTjnA7dYGLboRULfaId0o5yIiGRVjFcxiYhIEVCC\nEBGRrJQgREQkKyUIERHJSglCRESyUoIQEZGslCCk7JjZuNgQ0ZvMbEO0XmdmL3Xztd5pZod34/iZ\nZrY+y/7nzOy4Tp6zvDtxiRSCEoSUHXff5u5Hu/vRwA+Ab0Xriwg3hXXHuwk3F+X63muAdWZ2Wmqf\nmc0DRrr7X7v53iKJUoKQ/sBiy4Fm9qNoUpcHzGwogJkdamZ/iEY2fdTM5prZyYRB4q43s7+Z2SFm\n9gkzWxqdEdxlZsOyvN9ttL17+v3AbWY2I3rtZ6JyUrtAzS4zs+/Gtu8zs9Oj9XPM7InouXeY2Yi+\n+XpEslOCkP5mNvA9dz8C2EEYZgTgR8Cn3f1Y4IvATe7+BPBb4Avufoy7/50wh8Lx7r4IeJnsQxXc\nCbzLzFL/vy4mJI0twNnuvpiQNP47h3gdcDMbD3wZeGv0/GeAz3f3w4t0R9GNxSSSZ6+7+wvR+jOE\nCVNGEEbjvDMaMwvCYGcp8WHFF5rZfxEGPRtJGPemDXevNbMXgbPMbAvQ7O4rzKwK+J6ZHUUYYjnX\nIcANOJHQ1PVEFGMF8ESOzxfpESUI6W/2x9ZbgKGEM+m6qJ8im/iAZbcAF7j7cjP7CGF60WxSzUy1\nwC+jfZ8DNrn7h81sILAvy/OaaXtmPzS2/qC7f6CD9xPpc2pikv7OPIyR/7qZvRfAgiOjx3eTnkoS\nwlnDZjMbDHyok9f9DfB24H2EqR+JXmdztH4pYcTNTGuARVEM0wjDhTvwFHBKNBkNZjbCzGZ365OK\ndJMShPQH3sF6fPuDwOVm9hxhGOQLov2/Ar4YdQwfAnyFML/444Q+iKzDIXuYi/wJYHN0ZRPATcBH\noveYSxhavU0c7v448DqwAvgOoRkMD3NGX0bo7H4+eu1Cz/on/YyG+xYRkax0BiEiIlkpQYiISFZK\nECIikpUShIiIZKUEISIiWSlBiIhIVkoQIiKSlRKEiIhk9f8BvhIJCGjEL24AAAAASUVORK5CYII=\n", "text": [ "" ] } ], "prompt_number": 96 } ], "metadata": {} } ] }