{ "cells": [ { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "# 6: Central limit theorem\n", "\n", "\n", " #### [Back to main page](https://petrosyan.page/fall2020math3215)\n" ] }, { "cell_type": "code", "execution_count": 2, "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAlkAAAE/CAYAAAB1vdadAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjMsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+AADFEAAAgAElEQVR4nO3dd3gU1f7H8feXhNCLJYrSRVRQATU0wYp6wYYIijQ7iIJer3oV0atesfcG8gMEVJo0ERWwYUOKBC4IiiAiJSDNQhECJJzfH7PBiIFsyubsbj6v59lnZ2Zndj6bVfLNOTPnmHMOERERESlcJXwHEBEREYlHKrJEREREIkBFloiIiEgEqMgSERERiQAVWSIiIiIRoCJLREREJAJUZInIPmZW1szeN7MtZjbac5YEM9tuZjWiLVsoz3lmtjLb+lIzO6OQ3vsaM5saWk40M2dmtQrpvf/ycxWRyFGRJRIjzKyzmaWGfkH+bGZTzaxl6LWHzGzEfvtn/TLNeuw1s53Z1jvmcJqOwKHAYc65TkXwsQ7IOZfpnCvvnFsd2hQ12XLinDveOfflwfYxs2PNLNfBCZ1zrzvn2hRGLjObYWbXZnvv/X+uIhIhib4DiEjuzOwOoA/QE/gA2A20BtoCM3I6xjmXCZTP9h5pQFfn3GcHOVVNYKlzLiMfGRPzc1weRHO2QhNLWUXk4NSSJRLlzKwS8DDQyzk30Tn3h3Nuj3PuXefcvwvxPI8CfYEuoZaua8yshJk9YGarzGyjmQ03s4qh/Y8NdWNdZ2argQ9zeM8bzeyzbOt/6foysxFm9lKoVW6bmc0ys9r771sY2bJtu9bM0szsVzPrbmZNzWyRmf1uZi8e5OdT1szeNLPfzOxb4LT9Xk8zs7NDy83MbL6ZbTWzDWb2dGi3L0KvZ7UmNg79jL4I/Rx+Be7f/+cWcomZ/WRmm83sCTMrEXqvR8xseLYc+1rLzOxJoDkwMHS+F3L4DiqHvodNZrbSzO41M8v2/X1uZs+Hfj4rzOyCA/2MROSv1JIlEv2aA6WBtyN5EufcfaFfztWcc9cCmFkPoCtwNrAZGAG8CFyX7dAzgROA/M7R1ZmgVW5h6P37hc5Z2NmqhralAHWAVsAEYBpwLlAKWGhmY51zX+WQ82GgOnAMUBGYepDP9DLwtHNutJlVAOpny/ODcy57C2ND4HRgNJAMJAFdcnjPtsCpQCXgY+B7YPhBMuCcu8fMWgBDnHPDQ+fb/9/9AQT/fR0TOv+HwDrg9dDrp4eWDwNuAV4j+DmISC7UkiUS/Q4DNnvqQuoCPOOc+8k5t42gNalzVitKyIPOuR3OuZ35PMd451yqc24PMBJoFOFs/Zxzu5xzUwi6XUc45zY559IIul5POcD5rgQecc795pxbBbxykGx7gLpmdphzbptzbk4un2W1c+7V0PVSB/o5PhE690rgJaDA16WZWUmCz9UnlHMF8DzQLdtuPzrnhoa6n18HqpnZ4QU9t0hxoCJLJPr9AhyeQwtEUTgaWJVtfRVBS0tytm1rCniO9dmWd5DtOrJc5Cubc25DttWdwP7rBzr/Ufu936oD7AdBa1p9YKmZfW1mFx5k3xxz5rLPKoLPX1BHAAn8/edYNdv6/t8PhP8diRRrKrJEot8sIB24zMO51xFccJ6lBkHrz6asDc65g3UT/gGUzbZeJYqy5dV6/tpNdsAhEJxzS51zVxEUMc8CE8ysNAfuUg0n5/7nXhdazu1nfLD33ghk8vef49ow8ohILlRkiUQ559wW4AGgv5ldFroAu6SZtTGzp7LtWsLMSmd7lCqE048G7ghdfF4BeBQY7ZzbG+bxC4EGZnaymZUBHiyETIWVLa/GAn1DF4rXAHofaEcz62Zmh4eybCEodPYSFDXOzI7Jx/nvznbu24C3QtsXAGeZWXUzq0xwF2p2Gwiut/qbUBfteOAxMysfuungXwTXt4lIAanIEokBzrnngDuA+wlaatYQ/JKflG23TgTdXVmPHwvh1IMJfpl/CawAtgH/zEPu74DHgM+ApYTuriskBcqWDw8CPwMrCS56f+Mg+14ILDGzbcAzQEfn3O7QtWOPA3NCd+ul5OH87xIUVP8juAlieGj7tND6IuBrYPJ+x70AdAqd77kc3vcWghbAn4DPCa67OthnE5EwWeG2pouIiIgIqCVLREREJCJUZImIiIhEgIosERERkQhQkSUiIiISASqyRERERCIgKucubN26tZs2bZrvGCIiIiLhsJw2RmVL1ubNm31HEBERESmQqCyyRERERGKdiiwRERGRCFCRJSIiIhIBKrJEREREIkBFloiIiEgEqMgSERERiYCwiiwza21mS81suZn1Och+jc0s08w65PVYERERkXiSa5FlZglAf6ANUB/oZGb1D7Dfk8AHeT1WREREJN6E05LVBFjunFvhnNsNjAHa5rDfrcAEYGM+jhURERGJK+EUWVWBNdnW00Lb9jGzqkA7YGBejxURiTjn4I8/YO9e30lEpBgJp8jKaT4et9/6C8A9zrnMfBwb7GjWw8xSzSx106ZNYcQSETmAX36BgQOhQweoWRNKloTy5SExEY46Ci6+GJ59Ftasyf29RETyKZwJotOA6tnWqwHr9tsnBRhjZgCHAxeaWUaYxwLgnBsEDAJISUnJsRATETmo5cvh0Udh5EjYs+evr5UuDenpsH49vP9+8LjrLmjbFu6/H1JS/GQWkbgVTkvWXKCumdU2syTgKmBy9h2cc7Wdc7Wcc7WA8cAtzrlJ4RwrIlJgO3dC375Qrx4MHw6ZmfCPf8DgwfDtt7BrV7DPnj1BITZqFFxxBZQqBe+8A40bw/XXgyanF5FClGtLlnMuw8x6E9w1mAAMdc59a2Y9Q6/vfx1WrscWTnQREWDx4qBbcOlSMINrrw1apurU+fu+iYnB9jp1oFMn2LAh6DZ88UUYNgymToW33oIzzyzyjyEi8ceci76euZSUFJeamuo7hohEu3HjgqJqxw6oXx+GDIHmzfP+PsuWwY03wpdfQkICPPUU/OtfQdEmIpK7HP+x0IjvIhKbBgyAjh2DAuvqq2Hu3PwVWADHHQfTp0OfPkFX4513BtdrReEfoSISO1RkiUjseeIJ6NUrKIIefzy4Dqts2YK9Z2Ji8F5vvRXcjfjcc0HrloZ9EJF8UpElIrGlf3+4996gK2/w4KD1qTC79a68EiZPhjJlYOhQuOMOtWiJSL6oyBKR2DF6NPTuHSz/3/8FLU35UKV6TczswI82bThv5052A7z4IveVKHHw/XN4VKles9A+tojEpnDGyRIR8e/rr+G664LlJ5+E7t3z/VYb0laTMHjtQff5DOg67z3G/F9PHnWOH3r+HxNPuzj8c3TX5BYixZ1askQk+q1bB5ddFox31bMn3H13kZx24mkX8+8rHgBg6LB/ceLa74vkvCISH1RkiUh0y8gIrpP6+Wc444xgTKsi9OJ53RnZ9HLK79rBhP43UPmP34v0/CISu9RdKCLRrV8/+OorOPpoGD8ekpKK9vxm9Oz2FPV+XsapqxczYEQfOvd4VWNoSdxatgy2b/edouDKlw9GZ/FJRZaIRK8vv4RHHgkKmjffhCOO8BJjZ6kydOw5iPn/PZ8rU9/lvYbnM6pZey9ZRCJt+3aoVMl3ioLbssV3AnUXiki02rIFunQJxqm65x4491yvcX5KrsntnfoB8PKo+6i5eY3XPCIS/VRkiUh0uuceWLMGUlLg4Yd9pwHg9dOvZOKpF1Jp5zaGDvuXxs8SkYNSkSUi0eeLL4JxsEqWDCZuLlnSd6KAGTd3fZINFQ7nrGWzuPart3wnEpEopiJLRKJLevqfY2Ddey+cdJLfPPv5pcKh3NnxIQCeGteP5K2b/QYSkailIktEosvjjwe3N9WrB337+k6TozFNLuODE8/m0B2/8+xbD/mOIxI3Jk58nbZtT6Nhw4q0aFGNJ5+8m4yMjH2v//77r9x8cztOPrkcZ55Zk8mTR3lMmzsVWSISPVauhKeeCpYHDYJSpbzGOSAzenV5nB1Jpen89duc/f1XvhOJxIWdO3dw//0vMHfuZiZMmMPMmZ8wZMgz+15/6KFelCyZxOzZG3juuZE88MDNLFv2rcfEB6ciS0Six113Bd2FnTtDy5a+0xzUyuQaPNHmVgCeeeshSuzN9JxIJLLOOqsWQ4Y8w0UXNaBRo0rcdltHdu1KL9RzdOlyM40bn0FSUhJVqlTl0ku7MG9e8EfMjh1/8MEHE/jXv/pRrlx5UlJa0qrVpUya9GahZihMKrJEJDp8+ilMmABlywZzE8aA5y64iVWHVqVR2ndcN2OM7zgiETdlyliGDp3GZ5/9xNKl3zBhwvAc90tNncEpp1Q+4CM1dUZY55s79wvq1j0RgJ9+WkaJEgnUrv3nCKP16jXkhx+ityVLg5GKiH8ZGXDbbcHyvfdCtWp+84QpPakMfTrcx+hBt9Bv0pOMS7mErWUr+o4lEjFXX30bRx55NADnnnsJS5YsyHG/lJSW/O9/BZuCavz4YSxalMpjjw0BYMeO7VSo8NdRUsuXr8Qff2wr0HkiSS1ZIuLf66/D4sVQqxbceafvNHkyLuVSZtZJ4Yhtv9B3yku+44hEVHJylX3LpUuX5Y8/IjP/zkcfTeLpp/swdOhUDj30cADKli3P9u1b/7Lf9u1bKVeuQkQyFAYVWSLiV3o6PPRQsPzoo1CmjNc4eWbGHR3/C0DvT4ZS9dd1ngOJ+Dd37pc0aFD+gI+5c7884LGffz6Nvn27M2jQuxx//Mn7tteufRyZmRmsXPnDvm3ff79wX3diNFJ3oYj4NWAApKVBgwZw1VW+0+RLau1GjDvtYq6Y9x5933+RXt1i45oykUhp3PgMvvkm761cs2ZN5847uzBgwNs0bNjkL6+VLVuOCy64nBdeeIDHHhvCkiUL+Pjjdxg7dmZhxS50askSEX+2boXHHguWH38cSsTuP0kPtf03mVaC678aQ+1Nq3zHEYlJr7zSj23btnDjjRfua/W6/vo2+17/738HkJ6+k6ZNj+D22zvx8MOvctxx0duSZS4K595KSUlxqampvmOISKQ9+GAwL2HLlsFUOmZFclozI2Hw2kJ/39eG3s41s8bxRvMOXDNrPNH476tIbubPh0qVct8v2m3ZAqeeWmSny/Efr9j9s1FEYtuvv8LzzwfLjz9eZAVWJPW75A72JCTSZfZETvAdRkS8C6vIMrPWZrbUzJabWZ8cXm9rZt+Y2QIzSzWzltleW2lmi7JeK8zwIhLDXnwRtm2DCy6I+oFHw7UyuQavtexEgtvLw77DiIh3uRZZZpYA9AfaAPWBTmZWf7/dPgEaOucaAdcDQ/Z7/RznXCPnXEohZBaRWLdlS1BkAfznP397uUr1mphZxB6R9NhF/yQ9sRRXAHz3XUTPJSLRLZy7C5sAy51zKwDMbAzQFtj3r4dzLvstBOUAXYggIgfWv39QaJ11Vo6tWBvSVkfkmqksmd2rRuy91x1yFMNaduTmz94IRq5//fWInUtEols43YVVgTXZ1tNC2/7CzNqZ2ffA+wStWVkc8KGZzTOzHgUJKyJxYPt2eO65YDmHVqx48OwFN5MBMHJkMOm1iBRL4RRZObWt/62lyjn3tnPuBOAyoF+2l1o4504l6G7sZWZn5ngSsx6h67lSN23aFEYsEYlJAwfCL79A8+Zw7rm+00TEyuQajALIzISnn/YdR0Q8CafISgOqZ1uvBhxwSGPn3BdAHTM7PLS+LvS8EXiboPsxp+MGOedSnHMpycnJYcYXkZiSng7PPBMs339/XNxReCBPZC289hqsX+8zioh4Es41WXOBumZWG1gLXAV0zr6DmR0L/Oicc2Z2KpAE/GJm5YASzrltoeULQDfdiBRbI0fChg3QsCG0aZP7/jFsCUC7dvD228FQFU9qFHiJDeXLB5dMxrry5X0nCKPIcs5lmFlv4AMgARjqnPvWzHqGXh8ItAeuNrM9wE6gY6jgOhJ4O3Q3TyIwyjk3LUKfRUSimXN/Xot1111x3Yq1z733BkXWgAHQt298jPAoce+443wniB9hzV3onJsCTNlv28Bsy08Cf/szLXRHYsMCZhSReDBtWjCkQdWq0LGj7zRFo3FjOOcc+PTToNvwjjt8JxKRIqQR30WkaGRdi/XPf0LJkn6zFKV//St4fuklyMjwm0VEipSKLBGJvAULYPr04CKJ7t19pylaF10EdevCqlVB16GIFBsqskQk8p59Nnju3h0qV/abpaiVKAG33x4sZ83VKCLFgoosEYmstDQYMwYSEoKuwuLommvgkENg1iyYPdt3GhEpIiqyRCSyXn01uBapQweoWdN3Gj/KlYMeoQkv1JolUmyoyBKRyElPh0GDguXbbvObxbfevSExESZMCK7PEpG4pyJLRCJn3DjYvBkaNQqm0SnOqlWDK68MptoZODD3/UUk5qnIEpHI6d8/eO7du3gMPpqbXr2C59deg127/GYRkYhTkSUikZGaCnPmBBd8d+rkO010aN48mFJo06ag21BE4pqKLBGJjKxWrOuvh7Jl/WaJFmZw883B8oABfrOISMSFNa2OiEiebN4Mo0f/tagobhJLYTl0kZYD1gEVv/qKBmYsyufbH1mtBuvX6AJ6kWimIktECl/WNUcXXgh16vhO40fGLhIGr/3b5nTgjVH30/vTYdxyVjd6d30iX2+/oXvVAgYUkUhTd6GIFK7MzGBsLAgueJe/GXRWNwC6zp5AhZ3bPKcRkUhRkSUieValek3MLMdHm8REWLWKH4ESF154wP0O9oh331U9ns+Pa075XTvoOlsXwIvEK3UXikiebUhbnWNXGMBN/W+ABdN47fJ7KdEmfy1ZmcWgK+zVs6/mrGWzuOmzN3j17Gs0xIVIHFJLlogUmiq/b+Dibz5iT0Iir59+pe84Ue2dU1qzvmIyJ61byhk/zPEdR0QiQEWWiBSaa2aOJXFvJu81OJ8NlY7wHSeq7UlMYmjLYPywG74c6TmNiESCiiwRKRS2dy/XzxgNwGtndPacJjYMa3kVAO3nTaHSji2e04hIYVORJSKF4tzvZ1Bn0ypWHVqVD088y3ecmPBTck2mn9CCMnvS6TTnbd9xRKSQqcgSkUJxw5ejABjWshN7SyR4ThM7XmsZtPpdN2OM5yQiUthUZIlIgR2+7Rcu+980Mq0Ew1p09B0npkw6tTW/lq3MaasX0Wj1Yt9xRKQQqcgSkQLrNmscSZl7mHryuaw99GjfcWLKrpKlGdW0HQDXha5pE5H4oCJLRArGuX1dhUPO6OI5TGwaekZwl2HnOW9TevdOz2lEpLCEVWSZWWszW2pmy82sTw6vtzWzb8xsgZmlmlnLcI8Vkdh2+o+pnLD+R9ZVOpKpJ5/rO05M+qb6iaTWbMAhO7bQbv5U33FEpJDkWmSZWQLQH2gD1Ac6mVn9/Xb7BGjonGsEXA8MycOxIhLDrvnqLQDePP0KMhM0iUR+ZY2Zdb26DEXiRjgtWU2A5c65Fc653cAYoG32HZxz251zLrRaDnDhHisisavsrh1ckfouAK+ffoXnNLFtTJPL2JFUmnOWzuSYjSt9xxGRQhBOkVUVWJNtPS207S/MrJ2ZfQ+8T9CaFfaxIhKb2s2fQsX07cyqcxrLqhzrO05M21q2IhNOvQiA677ScA4i8SCcIiunWUvd3zY497Zz7gTgMqBfXo4FMLMeoeu5Ujdt2hRGLBHx7ZqZYwE0T2EhGRoaKf/qmeNIyMzwnEZECiqcIisNqJ5tvRqw7kA7O+e+AOqY2eF5OdY5N8g5l+KcS0lOTg4jloj4VOOXNM79/it2lizNuJRLfMeJC1/WbcoPR9Si6u/rabVkhu84IlJA4RRZc4G6ZlbbzJKAq4DJ2Xcws2PNzELLpwJJwC/hHCsisenqmeMAmHRKa7aUreQ5TZwwY2SzDgB0mT3BcxgRKahciyznXAbQG/gAWAKMdc59a2Y9zaxnaLf2wGIzW0BwN2FHF8jx2Eh8EBEpOgZcndVVqBHeC9XIZpcDcNn/plIu/Q/PaUSkIMK639o5NwWYst+2gdmWnwSeDPdYEYltZwDHbF7NmkOOYvoJLXzHiSs/JddkZp0UTv8xlbYLpjGqWXvfkUQknzTiu4jk2bWh5xHNO2gy6AjIas3qOktdhiKxTEWWiOTN9u1kjYiluwojY1zKJexOKEmrJV9S5fcNvuOISD6pyBKRvBk/nvLAV8c2ZvmRx/hOE5d+LX8oU08+lwS3l45z3/EdR0TySUWWiOTN8OGAWrEibUToWizdZSgSu1RkiUj4Vq6Ezz9nJ2hsrAib0qAVv5epyKmrF1N/7VLfcUQkH1RkiUj4Ro0CYBKwrUwFv1ni3K6SpRmfcjEAnedM9JxGRPJDRZaIhMc5ePNNAEZ4jlJcjAx1GXaa8za2d6/nNCKSVyqyRCQ88+bB999DcjIf+s5STMw4tgmrDq1KzV/XcsYPc3zHEZE8UpElIuEZEWq/6tQJTV1cNFyJEoxu2g5Ql6FILFKRJSK5y8iA0aOD5a5d/WYpZrK6DDukvkepPeme04hIXqjIEpHcffghbNwIxx8PKSm+0xQrS44+jnk1Tqbyzq1c9M3HvuOISB6oyBKR3GV1FXbrBmZ+sxRDo0LT7HSZrS5DkViiIktEDm7bNpg0KVju0sVvlmJqTJPLyLQStFk0nUO3/+o7joiESUWWiBzchAmwcyeccQbUquU7TbG0odIRfFz/TJIy93BF6ru+44hImFRkicjBZe8qFG9GqstQJOaoyBKRA0tLg+nTISkJrrjCd5pibdIpbdheqiyn/5jKMRtX+o4jImFQkSUiBzZqVDDS+yWXQOXKvtMUaztKlWXSKa0B6Dznbc9pRCQcKrJE5MDUVRhVRjbrAEDn2RM8JxGRcKjIEpGcLVwIixbBoYdCmza+0wgw/YQW/FzpCI7b+BNNfIcRkVypyBKRnIUmg6Zjx+CaLPEuMyGRtxq3BUDj7otEPxVZIvJ3mZnB9VigrsIoM6J50GV4FcCePV6ziMjBqcgSkb+bPh1+/hnq1IFmzXynkWwWVD+Rb486jmSADz7wHUdEDkJFloj8XdYF7126aBqdaGO2b5qdfd+TiEQlFVki8ld//AETQwNedtWVP9FoVNNQkfXOO7B1q98wInJAYRVZZtbazJaa2XIz65PD613M7JvQY6aZNcz22kozW2RmC8wstTDDi0gETJ4M27dD06ZQt67vNJKDNYdV5TOA9PRg2iMRiUq5FllmlgD0B9oA9YFOZlZ/v91+As5yzjUA+gGD9nv9HOdcI+dcSiFkFpFIyuqCUitWVNvXUaguQ5GoFU5LVhNguXNuhXNuNzAGaJt9B+fcTOfcb6HV2UC1wo0pIkVi48bgYuqEhGDoBola4wFKlYJPPw2mPxKRqBNOkVUVWJNtPS207UBuAKZmW3fAh2Y2z8x65D2iiBSZt94Khm9o3RqSk32nkYPYAsF0R879OdyGiESVcIqsnG4tcjnuaHYOQZF1T7bNLZxzpxJ0N/YyszMPcGwPM0s1s9RNmzaFEUtECp2m0YktWV266jIUiUrhFFlpQPVs69WAdfvvZGYNgCFAW+fcL1nbnXPrQs8bgbch59kgnHODnHMpzrmUZP0FLVL0li2Dr7+GChWCFhKJfm3aBNMeLVoE33zjO42I7CecImsuUNfMaptZEsFAw5Oz72BmNYCJQDfn3LJs28uZWYWsZeACYHFhhReRQjRyZPDcvj2ULes3i4QnKenPa+eypkESkaiRa5HlnMsAegMfAEuAsc65b82sp5n1DO32AHAYMGC/oRqOBGaY2ULga+B959y0Qv8UIlIwzumuwliV1bU7alRwPZ2IRI3EcHZyzk0Bpuy3bWC25RuBG3M4bgXQcP/tIhJlZs+GFSvg6KPh7LN9p5G8aNYMjjkm+P4++wxatfKdSERCNOK7iPzZitW5czB8g8QOsz9bH9VlKBJVVGSJxKEq1WtiZmE9kszYPGAAAA2feSasYyTKZBVZEybAjh1+s4jIPmF1F4pIbNmQtpqEwWvD2vfChR9y+CvXsajqCXz74MckhFFEZXY/2FB5UuTq1g2mQZozJ5gW6aqrfCcSEdSSJVLsdZ4dTAY9qunlQdeTxCZ1GYpEHRVZIsVYxR1buWThRwCMbnqZ5zRSIB07QmJiMC3Sxo2+04gIKrJEirV2/5tKmT3pfHZ8c9IOVRdgTEtODqZDyswMpkcSEe9UZIkUY11mTwBgZNP2npNIoVCXoUhUUZElUkxV/XUdZy+dSXpiKSaedqHvOFIYLr00mBZp7lxYutR3GpFiT0WWSDF11dx3KOEc7zU8jy1lK/mOI4WhTBno0CFYzpomSUS8UZElUkx1zuoqbKauwriS1WU4YkQwXZKIeKNxskSKoZPTvqNh2hJ+KVeZaSed4zuO5EdiqRwHhjVgFVD9p59oUaIEM/P59kdWq8H6NasKklCk2FORJVIMdZ79NgDjUi5hT2KS5zSSLxm7Djjg7Jjxj/LvDwbQ7axuzOn6RL7efoMGnBUpMHUXihQztncvneaEBiBVV2FcGtE8+F6vnPsuSXt2eU4jUnypyBIpZs5aNotqv6/np8OqM7NOiu84EgHfVj2BBdXqc+iO32m9+FPfcUSKLRVZIsXMvml0mmkanXiW1UqZNRaaiBQ9FVkixUjp3TtpP/99IDRXocSt0U0vY68ZF3/zMZX/+N13HJFiSUWWSDFy0TcfU2nnNubWbMjSo471HUci6OfKVZh+QktKZeymw7z3fMcRKZZUZIkUI12ydxVK3BsR6jLsNmu85yQixZOKLJFi4rBtv9Jm8XQySiTwVpO2vuNIEXj71Av5I6kMLZbPpfYmjXklUtRUZIkUEx3mvUvJzAw+rn8mGysm+44jReCP0uV4+9Q2wJ+tmCJSdFRkiRQTWb9kR+qC92JlZLNgLsOus8Zrmh2RIqYiS6QYqL1pFaf/mMr2UmV555TWvuNIEfqkXkvWVTqSYzetpNmKeb7jiBQrKrJEioHOoRHeJ53Smh2lynpOI0Vpb4kERjdtB2jMLJGiFlaRZWatzWypmS03sz45vN7FzL4JPWaaWcNwjxWRCHPuz5zShkMAAB8JSURBVAFIm2oaneIoa5qdjl9P1jQ7IkUo1yLLzBKA/kAboD7Qyczq77fbT8BZzrkGQD9gUB6OFZEISlm5kOM3rGB9xWQ+qdfSdxzxYFG1+vum2WmzaLrvOCLFRjgtWU2A5c65Fc653cAY4C/3fzvnZjrnfgutzgaqhXusiERWVlfhW03akpmQ6DmN+DIy1JrVVV2GIkUmnCKrKrAm23paaNuB3ABMzeexIlKIEjP2cNXXkwAY2UxdhcXZ6CbtyLQSXPTNxxy6/VffcUSKhXCKrJxmkM3xPmAzO4egyLonH8f2MLNUM0vdtGlTGLFEJDfnL/mCI7b9wpIqxzK/xsm+44hH6ysfySf1ziApcw9XpL7rO45IsRBOkZUGVM+2Xg1Yt/9OZtYAGAK0dc79kpdjAZxzg5xzKc65lORkDZQoUhg6Z59Gx3L6m0eKkzebh8bMUpehSJEIp8iaC9Q1s9pmlgRcBUzOvoOZ1QAmAt2cc8vycqyIREb59O20XTANCLqKRN45pTXbS5Wl+Y/zOHbDCt9xROJerkWWcy4D6A18ACwBxjrnvjWznmbWM7TbA8BhwAAzW2BmqQc7NgKfQ0T2c9n8qZTdnc6MY5uwMrmG7zgSBXaUKsvEUy8ENM2OSFEIa5ws59wU59xxzrk6zrlHQ9sGOucGhpZvdM4d4pxrFHqkHOxYEYm8LqG7Ckc20zQ68qcRoS7DLrMnaJodkQjTiO8icagKcO6SGexOKMn4lIt9x5Eo8tnxp5NWuQrHbF5Ni+VzfccRiWsqskTiUDcgwe1lysmt+K3cIb7jSBQJptkJWje7zRrnOY1IfFORJRJvnOO60OLwFh29RpHo9HqLKwG4cu5kyu7a4TmNSPxSkSUSb+bMoR6wocLhTDvpHN9pJAp9f1Rd5tQ+hYrp22k3f4rvOCJxS0WWSLwZNgwIplHJSCzpOYxEq6xWzmu/estzEpH4pSJLJJ7s2AGjRwMw/HR1FcqBvdW4LTtLluacpTOptWm17zgicUlFlkg8mTgRtm1jDvBd1eN9p5EotrVsRd4+pQ0AV+sCeJGIUJElEk9CXYXDPMeQ2JB1AfzVM8die/d6TiMSf1RkicSLlSth+nQoXZoxvrNITJh+QktWHVqVWr+kcfbSmb7jiMQdFVki8eL114Pndu3Y4jeJxAhXogRvnB60Zl0zUxfAixQ2FVki8WDvXhg+PFi+7rqD7iqS3euhIuvy+VOouGOr5zQi8UVFlkg8+PzzoLuwenU491zfaSSGrEyuwWfHN6fs7nSuSH3XdxyRuKIiSyQehC5455prICHBbxaJOVnDfVyrLkORQqUiSyTWbd0K48cHy9de6zWKxKaJp13E1tLlaf7jPI7/ebnvOCJxQ0WWSKwbOxZ27oQzz4Q6dXynkRi0o1RZxqVcAmgEeJHCpCJLJNYNGRI864J3KYBhLa8CgjGzSmbs9pxGJD6oyBKJZQsXwpw5UKkSXHml7zQSw2YfcxqLjz6eI7dt5uKFH/mOIxIXVGSJxLLBg4Pnrl2hbFm/WSS2mTHkzC4AdP9ipOcwIvFBRZZIrNqxA0aMCJa7d/ebReLCyGaXs7NkaS747nNq+Q4jEgdUZInEqnHjYMsWaNoUGjb0nUbiwG/lDmH8aRcBcKPnLCLxQEWWSKwaNCh47tHDbw6JK0PO7ArA9QB79njNIhLrVGSJxKLFi2HmTKhQATp29J1G4shXxzbmu6PqchTA++/7jiMS01RkicSi7Be8lyvnN4vEFzNeO6NzsJzVWioi+aIiSyTW7NwJb7wRLKurUCLgzeYd2AUwbRqsWuU7jkjMCqvIMrPWZrbUzJabWZ8cXj/BzGaZ2S4zu2u/11aa2SIzW2BmqYUVXKTYGj8efv8dGjeGRo18p5E49Gv5QxkP4By89prvOCIxK9ciy8wSgP5AG6A+0MnM6u+326/AbcAzB3ibc5xzjZxzKQUJKyLogncpEvs6Cl97DTIyfEYRiVnhtGQ1AZY751Y453YDY4C22Xdwzm10zs0FdCuKSCR99x3MmAHly8NVV/lOI3HsC4DjjoN163QBvEg+hVNkVQXWZFtPC20LlwM+NLN5ZnbAP73NrIeZpZpZ6qZNm/Lw9iLFyKuvBs9dugSFlkgk3XRT8Ny/v98cIjEqnCLLctjm8nCOFs65Uwm6G3uZ2Zk57eScG+ScS3HOpSQnJ+fh7UWKia1bYfjwYPmWW7xGkWLiuuugTBn46CNYtsx3GpGYE06RlQZUz7ZeDVgX7gmcc+tCzxuBtwm6H0Ukr958E7ZvhzPPhAYNfKeR4uCQQ4JWU4ABA/xmEYlB4RRZc4G6ZlbbzJKAq4DJ4by5mZUzswpZy8AFwOL8hhUptpz7s8umd2+/WaR46dUreB42LCjyRSRsuRZZzrkMoDfwAbAEGOuc+9bMeppZTwAzq2JmacAdwP1mlmZmFYEjgRlmthD4GnjfOTctUh9GJG59+iksWQJHHw2XXeY7jRQnjRrB6acH3dUjR/pOIxJTEsPZyTk3BZiy37aB2ZbXE3Qj7m8roJlrRQrqlVeC55tugpIl/WaR4qd372Aap/79g6FDLKdLdUVkfxrxXSTarV4N77wTFFcaG0t8aN8ejjwSFi2CL7/0nUYkZqjIEol2AwfC3r3QoQNUqeI7jRRHSUnQvXuwrOEcRMKmIkskmqWn/zkZtC54l6KUWAoz2/eo9sgjZAB7xo7l6Gzb8/uoUr2m708oEnFhXZMlIoWrSvWabEhbnet+3YA3gPnAaS1aRDqWyJ8ydpEweO2+1fXAO692p/38KfS66J88eNndBXr7Dd3zMqa1SGxSkSXiwYa01X/5BZYj57jt0Qth1TcMuOZZElqGP41Opn6BSQS83OpG2s+fQo/P3+TxC28lPamM70giUU3dhSJR6owf5pCy6hs2lT+UMU3a5n6ASITNqNuEuTUbkrz9V7rMnug7jkjUU5ElEqVu/2gQAAPPvkYtBhIdzHjx/OAC+H9+PDgYJFdEDkhFlkgUOnbDCi5Z+CHpiaV49ZxrfMcR2Wf8aReTVrkK9X/+gX98+5nvOCJRTUWWSBS67eMhlHCOkc0uZ2NFTZgu0SMjsSSvtLoB+LO1VURypiJLJMocuv1Xrp35FgAvntfdcxqRvxtyRmf+SCrD+d99wYlrv/cdRyRqqcgSiTI9Ph9B2d3pTDvpHL6rerzvOCJ/83u5ygxrEdzt+s+PBntOIxK9VGSJRJGkPbvo9ekwAF44X1PoSPR6+bwb2GtG5zlvc8TWTb7jiEQlFVkiUaTj3Hc4astGvqlaj4/rneE7jsgB/XhEbd5teAGlM3bRa/ow33FEopKKLJFo4Rz/+jC4kPiF87uDmedAIgf37D9uBuCWT4dTPn275zQi0UdFlkiUuOibj2mwdglrK1dhTJPLfMcRydXMYxvzZd2mHLJjCz0+H+E7jkjUUZElEg2c494pLwHw3AU3sbtkKc+BRMLzVJteQDCcQ9KeXZ7TiEQXFVkiUeCsZbNotmI+m8sfwpAzuviOIxK2qSedy4Jq9Tl6ywaunjXOdxyRqKIiSyQK9JnyMgAvt7qBP0qX85xGJA/MeKpNbwDumvYqCZkZngOJRA8VWSKepfy0gPO/+4JtpcrR/5zrfMcRybMJp13E8uRaHLtpJe3nve87jkjUUJEl4tk9U18Bgomgfy9X2XMakbzLTEjkmdbBnYZ3T31FE0eLhKjIEvGo3rpltPvfVNITSwXDNojEqDeaX8G6SkfSKO07Lv7mI99xRKKCiiwRj+59P7ijcFjLjmyodITnNCL5t7tkqX2tWf+Z/Jxas0RQkSXiTb11y7hq7iR2J5Tk6da9fMcRKbBBZ3ZlXaUjOW31Ii5Z+KHvOCLehVVkmVlrM1tqZsvNrE8Or59gZrPMbJeZ3ZWXY0WKq/+8+xwlnGPIGZ1ZfVg133FECiw9qcy+cbMefOcZbO9ez4lE/Mq1yDKzBKA/0AaoD3Qys/r77fYrcBvwTD6OFSl2TgauTH2X9MRSPBm6/V0kHgw+swtrK1ehUdp3tF0wzXccEa/CaclqAix3zq1wzu0GxgBts+/gnNvonJsL7MnrsSLF0UOh50FndWXtoUf7jCJSqHaVLL3vD4cHJj+n1iwp1sIpsqoCa7Ktp4W2haMgx4rEp/nzuRzYWbI0T+laLIlDr53RibTKVWiwdgnt/jfFdxwRb8IpsiyHbeHeNhL2sWbWw8xSzSx106ZNYb69SAx68EEAXj3nGtZXPtJzGJHCt6tkaZ648FYgaM0qsTfTcyIRP8IpstKA6tnWqwHrwnz/sI91zg1yzqU451KSk5PDfHuRGDNnDrz3Hn8AT//jFt9pRCJmaMtOrDysGietW0qX2RN8xxHxIpwiay5Q18xqm1kScBUwOcz3L8ixIvHFOfj3vwF4EdhU8XC/eUQiaHfJUjzYNvjv/b/vPE2pPemeE4kUvVyLLOdcBtAb+ABYAox1zn1rZj3NrCeAmVUxszTgDuB+M0szs4oHOjZSH0Ykqr37Lnz5JRx2GE/6ziJSBEY3bcfCavWo8es6bvl0uO84IkUurHGynHNTnHPHOefqOOceDW0b6JwbGFpe75yr5pyr6JyrHFreeqBjRYqdjAy4555g+YEH2Oo3jUiR2Fsigfsu7wvAve+/TKUdWzwnEilaGvFdpCgMHQrffw/HHAM9e/pOI1Jkpp10Dp8efzqH7vide6b29x1HpEipyBKJtO3b991RyOOPQ1KS3zwiRcmMe9vfB8Ctn7xGtV/Xeg4kUnRUZIlE2nPPwfr10LgxXHGF7zQiRS61diPGnXYxZfak8+jEJ3zHESkyKrJEImntWnjqqWD56afBcho6TiT+3dv+PtITS9FlzkSaL5/rO45IkVCRJRJJd98Nf/wBl10GZ53lO42INyuTa/DsP4LrEZ8b82COI1WLxBsVWSKR8uWXMGoUlC4ddBmKFHNPtulNWuUqNF61kKt9hxEpAiqyRCIhMxNuDaYV4e67oXZtv3lEosCOUmX3XQT/BMBWDWYi8U1Flkgk/N//wcKFUKPGn+NjiQijm7ZjVp3TqALwyCO+44hElIoskcK2eTPcf3+w/NxzULas3zwi0cSM26/qx16A55+HRYt8JxKJGBVZIoXt3/+G336DVq3g8st9pxGJOvNqNaQ/BDMh9OgBe/f6jiQSESqyRArTJ5/A8OFQqhQMGKAhG0QO4H6Ao4+G2bNh0CDfcUQiQkWWSGHZuRNuuilY/s9/4Ljj/OYRiWJbAV56KVjp0ycYsFckzqjIEiks//0v/PgjnHRS0GUoIgd3+eVw0UWwZQvcfrvvNCKFTkWWSGFYsACeeSboHhwyRPMTioTDDF55Jbg55K234J13fCcSKVQqskQKavduuOGGYGys3r2haVPfiURiR61awcTpEFwEv3mz1zgihUlFlkhBPfIIzJ8f/LJ49FHfaURiT+/ecPbZsHEj9OrlO41IoVGRJVIQc+bAY48F3R6vvw4VKvhOJBJ7SpSAoUOhfHkYOzZ4iMQBFVki+bVjB3TrFnQT3nknnHmm70Qisat27eC6RoBbbtHdhhIXEn0HEIlZ99wDP/wQ3E3Yr5/vNCKxJbEUlsM4ctOAf/zyC9OOOooLAZfPtz+yWg3Wr1lVkIQiBaYiSyQ/3nsvuCuqZEl4800oXdp3IpHYkrGLhMFr/7b5xt/XM/+/59N6+6/c3f4+nml9S77efkP3qgVNKFJg6i4UyavVq+Hqq4PlRx+FRo385hGJIz9XrsIN1z0PQL9JT9L0x3meE4nkn4oskbzYvRs6dgzmJrz44uBaLBEpVFManMfz5/egZGYGIwb3otKOLb4jieSLiiyRvOjbN5hrrXr1YI7CEvpfSCQS+l5+L6k1G1D7lzUMGX4nuPxenSXij35DiIRr0iR49llITAxGpz7sMN+JROLWnsQkOvd4lS1lKtDuf1PpM+Vl35FE8iysIsvMWpvZUjNbbmZ9cnjdzOyl0OvfmNmp2V5baWaLzGyBmaUWZniRIrNoUTBcA8ATT0Dz5n7ziBQDK46oxdU3vMxeMx5+5ynaLPrEdySRPMm1yDKzBKA/0AaoD3Qys/r77dYGqBt69ABe3e/1c5xzjZxzKQWPLBJ5VarXxMwwMw4346cGDWD7dkYCdtdd+17L70NEwvN+w/N56NK7KOEcbw7uTd31P/qOJBK2cIZwaAIsd86tADCzMUBb4Lts+7QF3nDOOWC2mVU2s6Occz8XemKRIrAhbTUJg9dSMmM3E57vTO1ls5hbqxE3/Xs8CUllCvz+mbq9XCRsj194G43WfMvl86cwYcANtOwzma1lK/qOJZKrcLoLqwJrsq2nhbaFu48DPjSzeWbW40AnMbMeZpZqZqmbNm0KI5ZIhDnHi6Pv56xls1hX6Uja3zKE9EIosEQkb1yJElx/3fMsPvp46v/8A+MGdqdkxm7fsURyFU6RlVPfxv63eRxsnxbOuVMJuhR7mVmOc4845wY551KccynJyclhxBKJrPvef4EeX4wkPbEU7Xu9xrpDjvIdSaTY2l66PG1vfZ31FZNptWQGg16/S3ccStQLp8hKA6pnW68GrAt3H+dc1vNG4G2C7keRqHYj8N93niHTStC1+yvMrX2K70gixd6qw6tz6W1vsL1UWbrNnsDDk57yHUnkoMIpsuYCdc2stpklAVcBk/fbZzJwdeguw2bAFufcz2ZWzswqAJhZOeACYHEh5hcpfJMmMTC0eGvnR5l06oVe44jIn+bXbECnmwaSUSKBvlNe4taPh/iOJHJAuRZZzrkMoDfwAbAEGOuc+9bMeppZz9BuU4AVwHJgMJA12dSRwAwzWwh8DbzvnJtWyJ9BpPBMnw6dOpEAPHzJHQw6+2rfiURkP1NPbkXPbkEr1vNvPcgNX4z0nEgkZ2FNEO2cm0JQSGXfNjDbsgN65XDcCqBhATOKFI3p04OpctLTeZWgyBKR6DS85VWU3/UHL4x5gFdH3MPOpNKMatbedyyRv9CI7yLwZ4G1cyfccEPwF4PGsxKJaq+0uoG+l99LCecYNvR2Os2e6DuSyF+oyBL55JO/FFgMGvS322dFJDo91aY3/S7+FwluL68PvY3un7/pO5LIPmF1F4rErbFjg+lydu/eV2Bp0meR2PLftnexI6kMj098jFdH9KFi+nZ036FEAxVZUny9/DL885/BWDu33govvKACSyRGPd2mF9tKl+OVUffx5PhHqAywd6/+nxav9F+fFD9790LfvnDbbUGB9fjj8OKL+sdYJMYNPOdarrn+RTJKJHAvQMeOsGOH71hSjOm3ihQv27ZB+/ZBYZWQAMOGQZ8+ushdJE6MbN6BS257ky0A48fD2WfD+vWeU0lxpSJLio/ly6F5c5g0CSpXhvffh2uv9Z1KRArZRyeexekAtWrB3LnQuDHMnu05lRRHKrKkeHjvPWjSBL79FurVg6+/hn/8w3cqEYmQ7wDmzAn+sEpLgzPPhJde0nyHUqRUZEl827kzuKj9kkvgt9/g0kuDv2jr1vWdTEQi7Ygj4LPPgusv9+wJbnTp2BG2bPGdTIoJFVkSc6pUr4mZ5fo4yYxvypaFV15hN3AnUGLyZKxSpVyPFZE4kZQU3NgydixUqADjxsHJJ8NHH/lOJsWAhnCQmLMhbTUJg9ce8PWSGbvpM+UV7p3yEkmZe1h2RG269hjA/JoNwv6rIrN71cIJKyLR4YoroFEj6No1uFzgggvg5pvhqaegfHnf6SROqSVL4krTH+cxt19rHnz3WZIy9zDozC40/s8HzK/ZwHc0EfGtbl346it49FEoWRJefRVOPDG4GUbXakkEqMiSuFDl9w0MHn4HXz7ZlpPWLeWHI2rR6q5x3NLtKf4oXc53PBGJFomJwTh5qalw6qmwejW0awcXXQQ//ug7ncQZFVkS00rv3kmf91/i+/tbct1Xb5FRIpEn2/TmlAc/5vPjT/cdT0SiVYMGQbfhK69ApUowdSrUrw+33w4bN/pOJ3FCRZbEpNK7d3Lbx4P5oe/pPDLpScrv2sGkRq05+eFPue/ye0lPKuM7oohEu4QE6NULli6Fa64J7kB88UWoUwceeEB3IUqB6cJ3iS3btnE7cHff0zlqS/DX5rwaJ9Onw/18Wq+l32wiEj0SS+X5TuGTgUeBS7Zvh379+L1fPwYALwEb9tv3yGo1WL9mVeFklbilIktiw4oVwYTOQ4fyPMCWjcyrcTL9Lr2D9xqcr2lxROSvMnYd9C7knHwHtANOXz6Xfm8/yVnLZtEXuDMxiRHN2vNyqxtYXK0eABt0B7KEQUWWRK9du4KR2ocPD6bACd398wXwXO9hKq5EJCJmHtuYVv8eT9Mf53HnhwO57H9TuWHGaG6YMZrZx5zKkDO6MNp3SIkJuiZLosvevTBrFtxyCxx1FHToEBRaJUsG10zMm8dZwHsNL1CBJSIRNafOaVx582BOfPhzBpx9DVvKVKDZivkMef1Ofga4+urgD8Ddu31HlSillizxb/du+PTTYKyad96Bn3/+87WGDYPiqkuXYIoMEZEi9kOVOtzW5TH6dLifDqnv0v3LkTT/cR68+WbwqFw5GAbikkugVSuoWNF3ZIkSKrKk6O3dC4sWwSefBI8vvoDt2/98vUYNaN8+KK4aNvSXU0Qkmx2lyvJGi4680aIjx3SvyrJ+/YLpehYtgmHDgkdiIrRoAa1bw3nnBaPMJ+pXbXGlb14ib8sWmDsX5swJxqWZORM2b/7rPiedFPwleNllcMop6goUkaj2A8D99wePJUtg4kSYNi243OHzz4PHvfdC2bLQtGlQeLVoASkpcPjhvuNLEVGRJYVr5UquSWlM9V82czLQEDghh93WAJ+EHtOBdYsXw+LF0K9fEYYVESkE9erBffcFj99+C1rop00LWul/+CG4HOLTT//c/+ijgxauhg2D5xNOgGOPDQoyiSthFVlm1hp4EUgAhjjnntjvdQu9fiGwA7jWOTc/nGMlznTpwuu//LWValdiEguqn8icY05lbq1GzK3diOVH1P5La1VCHk6hyZtFJGodckhww06HDsH6xo1B6/2MGUEr18KFsG5d8Jgy5a/HVq0azK+Y9TjvvKBlX2JWrkWWmSUA/YHzgTRgrplNds59l223NkDd0KMp8CrQNMxjJZ6cfTYfzpzJ4gtuYlHVE1hctR6Lqx7PnsQk38lERApPPgY7BTCgDkErf6PQ83HAMUDJtWth7Vr47LNg56efVpEV48JpyWoCLHfOrQAwszFAW4Jx27K0Bd5wzjlgtplVNrOjgFphHCvx5NFH+cdjj5FwxQO+k4iIRE4+BjvN8lPoMSnbtoTMDGr8upa6G36izqaV1Bl1H7e31CwWsS6cIqsqwSU0WdIIWqty26dqmMeKiIgUa5kJifyUXJOfkmsG66Pu4/ZmzTynkoIyFxpF+4A7mF0B/MM5d2NovRvQxDl3a7Z93gced87NCK1/AtxN0AJ60GOzvUcPoEdo9XhgaQE/W7Q6HNic614SrfT9xTZ9f7FL311si/fvb7NzrvX+G8NpyUoDqmdbrwasC3OfpDCOBcA5NwgYFEaemGZmqc65FN85JH/0/cU2fX+xS99dbCuu31840+rMBeqaWW0zSwKuAibvt89k4GoLNAO2OOd+DvNYERERkbiTa0uWcy7DzHoDHxDcaT/UOfetmfUMvT4QmEIwfMNygiEcrjvYsRH5JCIiIiJRJKxxspxzUwgKqezbBmZbdkCvcI8t5uK+SzTO6fuLbfr+Ype+u9hWLL+/XC98FxEREZG8C+eaLBERERHJIxVZHpnZXWbmzEyzhcYQM3vazL43s2/M7G0zq+w7kxycmbU2s6VmttzM+vjOI+Ezs+pm9qmZLTGzb83sn74zSd6YWYKZ/c/M3vOdpaipyPLEzKoTTDe02ncWybOPgJOccw2AZcC9nvPIQWSb3qsNUB/oZGb1/aaSPMgA7nTO1QOaAb30/cWcfwJLfIfwQUWWP88TDNiqi+JijHPuQ+dcRmh1NsH4bxK99k0N5pzbDWRN7yUxwDn3s3Nufmh5G8Eva80SHyPMrBpwETDEdxYfVGR5YGaXAmudcwt9Z5ECux6Y6juEHNSBpv2SGGNmtYBTgDl+k0gevEDQoLDXdxAfwhrCQfLOzD4GquTw0n1AX+CCok0keXGw7885905on/sIujJGFmU2yTPLYZtakGOMmZUHJgC3O+e2+s4juTOzi4GNzrl5Zna27zw+qMiKEOfceTltN7OTgdrAQjODoKtpvpk1cc6tL8KIchAH+v6ymNk1wMVAK6dxUKJdOFODSRQzs5IEBdZI59xE33kkbC2AS83sQqA0UNHMRjjnunrOVWQ0TpZnZrYSSHHOxfPEmXHFzFoDzwFnOec2+c4jB2dmiQQ3KLQC1hJM99VZs0/EBgv+Gn0d+NU5d7vvPJI/oZasu5xzF/vOUpR0TZZI3r0CVAA+MrMFZjYwtwPEn9BNClnTey0BxqrAiiktgG7AuaH/3xaEWkZEop5askREREQiQC1ZIiIiIhGgIktEREQkAlRkiYiIiESAiiwRERGRCFCRJSIiIhIBKrJEREREIkBFloiIiEgEqMgSERERiYD/B9AwVa7xzLDdAAAAAElFTkSuQmCC\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "# nbi:hide_in\n", "import numpy as np\n", "import matplotlib.pyplot as plt\n", "plt.rcParams[\"figure.figsize\"] = (10, 5)\n", "\n", "N=1000 #max number of samples\n", "num = 1000 # number of experiments\n", "a=2\n", "b=4\n", "# intsize the number of intervals\n", "plot_width = 2\n", "\n", "data =np.random.rand(N, num)*(b-a)+a\n", "mean = (b+a)/2\n", "sigma = np.sqrt((b-a)**2/12)\n", "\n", "def pdf_func(xdata, mu, sigma):\n", " val = np.exp(-np.power(xdata-mu,2)/(2*sigma**2))/(sigma *np.sqrt(2*np.pi))\n", " return val\n", "\n", "def epmf(x, inter, intsize):\n", " epmf_values = np.zeros(intsize-1)\n", " for i in range(intsize-1): \n", " length = inter[i+1]-inter[i]\n", " epmf_values[i] = np.sum((inter[i]<=x) & (x" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "# nbi:hide_in\n", "import numpy as np\n", "import matplotlib.pyplot as plt\n", "plt.rcParams[\"figure.figsize\"] = (10, 5)\n", "\n", "N=1000 #max number of samples\n", "num = 1000 # number of experiments\n", "# intsize the number of intervals\n", "theta=2\n", "\n", "data =np.random.exponential(theta, [N, num] )\n", "mean = theta\n", "sigma = theta\n", "\n", "def pdf_func(xdata, mu, sigma):\n", " val = np.exp(-np.power(xdata-mu,2)/(2*sigma**2))/(sigma *np.sqrt(2*np.pi))\n", " return val\n", "\n", "def epmf(x, inter, intsize):\n", " epmf_values = np.zeros(intsize-1)\n", " for i in range(intsize-1): \n", " length = inter[i+1]-inter[i]\n", " epmf_values[i] = np.sum((inter[i]<=x) & (x" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "# nbi:hide_in\n", "import matplotlib.pyplot as plt\n", "import numpy as np\n", "from scipy.special import comb, factorial\n", "from scipy.stats import poisson\n", "\n", "\n", "def pdf_func(xdata, mu, sigma):\n", " val = np.exp(-np.power(xdata-mu,2)/(2*sigma**2))/(sigma *np.sqrt(2*np.pi))\n", " return val\n", "\n", "def poiss_binom_pmf(n, p):\n", " lam = n*p\n", " q = 1-p\n", " N = 50\n", " brange_x = np.arange(0, np.minimum(n, 100), 1)\n", " prange_x = np.arange(0, np.minimum(n, 100), 1)\n", " \n", " mean = n*p\n", " sigma = np.sqrt(n*p*q)\n", " xvalues = np.linspace(mean-15,mean+15, 1000)\n", "\n", " def poiss_pmf(x, lam):\n", " pmf_val = np.exp(-lam) * np.power(lam, x) / factorial(x)\n", " return pmf_val\n", "\n", "# ppmf_values = np.array([poiss_pmf(x, lam) for x in prange_x])\n", " ppmf_values = np.array([poisson.pmf(x, lam) for x in prange_x])\n", " def binom_pmf(x):\n", " pmf_val = comb(n, x, exact=True) * p**x * q**(n-x)\n", " return pmf_val\n", " mean = n*p\n", "\n", " bpmf_values = np.array([binom_pmf(x) for x in brange_x])\n", "\n", " # plot setup\n", " plt.figure(figsize=(14,7)) \n", " plt.axhline(y=0, color='k')\n", " plt.gca().spines['top'].set_visible(False)\n", " plt.gca().spines['right'].set_visible(False)\n", " \n", " # Plotting poisson\n", " plt.bar(prange_x, ppmf_values, width=1, color=(0.1, 0.1, 1, 0.15), edgecolor=(0.1, 0.1, 1, 0.3), \n", " linewidth=1.3, label=\"Poisson\",zorder=-1)\n", " \n", " # Plotting binomial\n", " plt.bar(brange_x, bpmf_values, width=1, color=(0.3, 0.5, 0.3, 0.25), edgecolor=(0.1, 0.1, 1, 0.3),\n", " linewidth=1.3, label=\"Binomial\", zorder=-2)\n", " \n", " # Plotting normal for Binomial\n", " plt.plot(xvalues, pdf_func(xvalues, mean, sigma), linewidth=3, color=\"green\", \n", " label=r\"normal $\\approx$ Binomial\", zorder=3)\n", " \n", " # Plotting normal for Poisson\n", " plt.plot(xvalues, pdf_func(xvalues, lam, np.sqrt(lam)), linewidth=3, color=\"blue\", \n", " label=r\"normal $\\approx$ Poisson\", zorder=3)\n", " \n", " \n", " plt.figtext(0.81,0.7, \" n = {}\".format(n)+\"\\n p = {}\".format(p), ha=\"left\", va=\"top\",\n", " backgroundcolor=(0.1, 0.1, 1, 0.15), fontsize=\"large\")\n", " plt.legend()\n", " plt.plot();\n", "\n", "poiss_binom_pmf(50, 0.5)" ] } ], "metadata": { "kernelspec": { "display_name": "Python 3", "language": "python", "name": "python3" }, "language_info": { "codemirror_mode": { "name": "ipython", "version": 3 }, "file_extension": ".py", "mimetype": "text/x-python", "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", "version": "3.7.7" } }, "nbformat": 4, "nbformat_minor": 4 }