{ "cells": [ { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "# Week 9-11: Central limit theorem, Descriptive statistics\n", "\n", "\n", " #### [Back to main page](https://petrosyan.page/fall2020math3215)\n" ] }, { "cell_type": "code", "execution_count": 1, "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAlkAAAE/CAYAAAB1vdadAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjMsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+AADFEAAAgAElEQVR4nOzdd3hUVf7H8fc3CaFLR5TelVWxBARRsSCChWKjKyoga9fdnyK6rmVVXLsrioArKF2KoiL2Li0qgkhZBMWgAtJ7CDm/P+4EYwxkApOcKZ/X89xn7p05d+5nMpRvzj33XHPOISIiIiKRleQ7gIiIiEg8UpElIiIiUgRUZImIiIgUARVZIiIiIkVARZaIiIhIEVCRJSIiIlIEVGSJyD5mVsbM3jSzzWY23nOWZDPbZmZ1oi1bKE87M/sh1/ZSMzstQu99hZm9FVpPMTNnZvUi9N5/+LmKSNFRkSUSI8ysp5mlh/6D/MXM3jKzU0Ov3WNmY/K0z/nPNGfJNrOduba75XOYbkBloIpzrkcxfKz9cs7tdc6Vc86tCj0VNdny45xr6pz79EBtzKyRmRU4OaFzbrRzrmMkcpnZZ2bWN9d75/25ikgRSfEdQEQKZma3AoOAgcDbQCbQAegMfJbfPs65vUC5XO+RAfR2zn10gEPVBZY657IOImPKwexXCNGcLWJiKauIHJh6skSinJlVAO4DrnPOTXXObXfO7XHOve6c+78IHucBYDDQK9TTdYWZJZnZ3Wb2o5mtNbNRZnZYqH2j0GmsK81sFfBOPu/Zz8w+yrX9h1NfZjbGzJ4O9cptNbNZZlY/b9tIZMv1XF8zyzCzDWbW38xONrOFZrbJzJ46wM+njJm9bGYbzWwRcFKe1zPM7IzQeisz+8rMtpjZGjN7JNTsk9DrOb2JLUI/o09CP4cNwF15f24hF5rZSjP7zcyGmFlS6L3+ZWajcuXY11tmZg8DrYFhoeM9mc93UDH0Pawzsx/M7A4zs1zf38dm9kTo57PCzNrv72ckIn+kniyR6NcaKAVMK8qDOOfuDP3nXMs51xfAzAYAvYEzgN+AMcBTwJW5dj0dOAo42Ht09STolfsm9P73h44Z6Ww1Q8+lAQ2Bs4EpwEzgLKAk8I2ZTXLOfZ5PzvuA2kAD4DDgrQN8pv8AjzjnxptZeaBZrjz/c87l7mFsDpwCjAeqAalAr3zeszNwIlABeA9YAow6QAacc7ebWRtgpHNuVOh4ef/df5bgz1eD0PHfAX4GRodePyW0XgW4FniB4OcgIgVQT5ZI9KsC/ObpFFIv4FHn3Ern3FaC3qSeOb0oIf90zu1wzu08yGNMds6lO+f2AGOB44s42/3Oud3OuRkEp13HOOfWOecyCE69nrCf410G/Ms5t9E59yPwzAGy7QEam1kV59xW59ycAj7LKufcc6HxUvv7OQ4JHfsH4GngkMelmVkJgs81KJRzBfAE0CdXs++dc/8NnX4eDdQys6qHemyRRKAiSyT6rQeq5tMDURyOBH7Mtf0jQU9LtVzP/XSIx/g11/oOco0jK8BBZXPOrcm1uRPIu72/4x+R5/1+3E87CHrTmgFLzWyumZ13gLb55iygzY8En/9QVQeS+fPPsWau7bzfD4T/HYkkNBVZItFvFrAL6OLh2D8TDDjPUYeg92ddzhPOuQOdJtwOlMm1XSOKshXWr/zxNNl+p0Bwzi11znUnKGIeA6aYWSn2f0o1nJx5j/1zaL2gn/GB3nstsJc//xxXh5FHRAqgIkskyjnnNgN3A0PNrEtoAHYJM+toZv/O1TTJzErlWkpG4PDjgVtDg8/LAw8A451z2WHu/w1wnJkda2algX9GIFOkshXWJGBwaKB4HeD6/TU0sz5mVjWUZTNBoZNNUNQ4M2twEMe/LdexbwQmhp6fD7Q1s9pmVpHgKtTc1hCMt/qT0CnaycCDZlYudNHBLQTj20TkEKnIEokBzrnHgVuBuwh6an4i+E/+1VzNehCc7spZvo/AoUcQ/Gf+KbAC2ArcVIjc3wEPAh8BSwldXRchh5TtIPwT+AX4gWDQ+0sHaHsesNjMtgKPAt2cc5mhsWMPAXNCV+ulFeL4rxMUVF8TXAQxKvT8zND2QmAuMD3Pfk8CPULHezyf972WoAdwJfAxwbirA302EQmTRbY3XURERERAPVkiIiIiRUJFloiIiEgRUJElIiIiUgRUZImIiIgUARVZIiIiIkUgKu9dWLVqVVevXj3fMUREREQK9OWXX/7mnKuW9/moLLLq1atHenq67xgiIiIiBTKzfG+zpdOFIiIiIkVARZaIiIhIEVCRJSIiIlIEVGSJiIiIFAEVWSIiIiJFQEWWiIiISBEIq8gysw5mttTMlpvZoAO0a2Fme83sksLuKyIiIhJPCiyyzCwZGAp0BJoBPcys2X7aPQy8Xdh9RUREROJNOD1ZLYHlzrkVzrlMYALQOZ92NwBTgLUHsa+IiIhIXAmnyKoJ/JRrOyP03D5mVhPoCgwr7L4iIkUuKwt27oQ9e3wnEZEEEk6RZfk85/JsPwnc7pzbexD7Bg3NBphZupmlr1u3LoxYIiL7sX49PPccXHwx1KsHpUpBmTKQmgo1a0LHjvDoo/DTTwW+lYjIwQrn3oUZQO1c27WAn/O0SQMmmBlAVeA8M8sKc18AnHPDgeEAaWlp+RZiIiIHtHIl3HcfjB37x14rs6DQysyEn38Olpkz4bbb4IILgn2OP95fbhGJS+H0ZM0DGptZfTNLBboD03M3cM7Vd87Vc87VAyYD1zrnXg1nXxGRQ7ZrFwwaBE2bwqhRsHcvdOgAw4fD4sXBqcKdO2H3bvj++6AIu+wyKFECXn8dTjgB+vaFjRt9fxIRiSMFFlnOuSzgeoKrBhcDk5xzi8xsoJkNPJh9Dz22iEjIggVw0knw8MPB2KvevWHZMnjrLejfH446CkqWDNqmpECDBtCzJ0ycGJwuvPnm4DTi6NFwzDHw3nt+P4+IxA1zLvrOzKWlpbn09HTfMUQk2k2bFhRVO3ZAkyZBL1br1oV/n//9D664AmbNgqQkeOIJuOGG4DSjiEgBzOxL51xa3uc147uIxKannoKLLgoKrD594OuvD67AAmjcGD79FO66C7Kz4aabgiU7O7KZRSShqMgSkdjzyCPBaT6AIUOCU31lyhzaeyYnw/33w7hxwenF//wHrr1WhZaIHLRwri4UEYkeTz0VXBVoBs8/H4y7iqQePaByZejSJXj/lJSg4NKpQxEpJPVkiUhUqVG7LmaW73KpGdmhHqyrnMMGDNhv23CWGrXr5h/i3HODqw5LloShQ4M5tURECkk9WSISVdZkrCJ5xOo/Pd96+TzGPHYZSVmZ3HHRYF7qeB3Jh3qs/ge4AUW7dvDyy8FUD7fdBnXrBusiImFST5aIRL3qW9Yx4fmBlMzK5NkzruCRDtcWz4EvvRT+/e9g/YorgukiRETCpCJLRKJa8t4sxj//V2pu+pXPGrXk1m73Fu/4qL//Ha66Kpjw9JJLYMuW4ju2iMQ0nS4Ukah232uP0HbZLH6pUJ3u1wwjK6VE8QYwg2eegfT0oCdrwAAYP14D4SVuLVsG27b5TnHoypULps/zSUWWiEStU5fN4f9mDmWvJdFjwDB+rXi4nyClS8MrrwQzy0+cCO3bB71bInFo2zaoUMF3ikO3ebPvBDpdKCJRqtyubfz3xZtJco6HO17PZ01O9huoSRN47rlg/eabYdUqv3lEJOqpyBKRqPTopHtp8NsqvqpzDPdfeIvvOIFevaBrV9i6Fa6+GqLwtmQiEj1UZIlI1Dljyef0+3Qcu1JK0veqp9mTkuo7UsAMhg2DqlWDG0k//7zvRCISxVRkiUhUSQWGjhkEwIPn38h3NZv6DZRX9erw7LPB+h13wNq1fvOISNRSkSUiUWUQ0HTNChbXaMSj5/7Vd5z8XXIJdOgAmzYFE5WKSERMnTqazp1Ponnzw2jTphYPP3wbWVlZ+17ftGkDf/1rV449tiynn16X6dPHeUxbMBVZIhI9li9ncGj1ut5DyCxR0muc/TIL7mdYsmRwc+pPP/WdSCQu7Ny5g7vuepJ5835jypQ5fPHF+4wc+fttre655zpKlEhl9uw1PP74WO6++68sW7bIY+IDU5ElItHj//6PksDo1pfySdPWvtMcWKNGcPvtwfq118KePX7ziBSxtm3rMXLko5x//nEcf3wFbryxG7t374roMXr1+istWpxGamoqNWrUpFOnXnz55ecA7NixnbffnsItt9xP2bLlSEs7lbPP7sSrr74c0QyRpCJLRKLDhx/Cq6+yDbjzojt8pwnPoEFQvz58+y0MH+47jUiRmzFjEv/970w++mglS5cuYMqUUfm2S0//jBNOqLjfJT39s7CON2/eJzRu/BcAVq5cRlJSMvXr/z7D6NFHN+d//4venixNRioi/u3dC7feCsAQ8DfpaGGVLg2PPgoXXwz33gt9+sBhh/lOJVJkLr/8Rg4//EgAzjrrQhYvnp9vu7S0U/n6602HdKzJk19k4cJ0HnxwJAA7dmyjfPk/zpJarlwFtm/fekjHKUrqyRIR/0aPhvnzoVYtHvedpbC6doVTToF16+CRR3ynESlS1arV2LdeqlQZtm8vmvvvvPvuqzzyyCD++9+3qFy5KgBlypRj27Y/3jt027YtlC1bvkgyRIKKLBHxa8cOuOuuYH3IEHb6TVN4Zr8XV489BqtX+80jEgXmzfuU444rt99l3rz9Xyzy8cczGTy4P8OHv07Tpsfue75+/Sbs3ZvFDz/8b99zS5Z8s+90YjTS6UIR8WvoUPjlFzjxROjRA3r39p2o8E45JThlOGUK3H03vPCC70QiXrVocRoLFhS+l2vWrA/429968eyz02jevOUfXitTpizt21/Ek0/ezYMPjmTx4vm8995rTJr0RaRiR5x6skTEny1b4OGHg/UHHoCkGP4n6aGHICUFRo2CZct8pxGJSc88cz9bt26mX7/z9vV6XXVVx32v33vvs+zatZOTT67OzTf34L77nqNJE/VkiYj82VNPwfr10KYNnHtu8R8/pSRmFrG3ex4YALzUtClX5Hnt8Fp1+PWnHyN2LJHi9vHHP/xh+6ab7on4McaO/fCAr1esWJlhw16N+HGLioosEfFjw4bgyjwIerEiWOyELWs3ySMiN4bq4d9+4sq7TqVXdjZD7vuQZTUa7XttTf+aETuOiMSGsPrmzayDmS01s+VmNiif1zub2QIzm29m6WZ2aq7XfjCzhTmvRTK8iMSwxx4LThe2awdt2/pOExE/Vq3Ni226keyyGfzmU77jiIhnBRZZZpYMDAU6As2AHmbWLE+z94HmzrnjgauAkXleP9M5d7xzLi0CmUUk1m3aFNyWBuD++/1mibAh591IZnIJesx5lSa/LvcdR0Q8CqcnqyWw3Dm3wjmXCUwAOudu4Jzb5pxzoc2ygENEZH+GDoWtW+Hss6FVK99pImpVlVr7erPufEO9WSKJLJwiqybwU67tjNBzf2BmXc1sCfAmQW9WDge8Y2ZfmtmAQwkrInFg+3Z44olgffDgA7eNUft6s+ZOo9GaFb7jiIgn4RRZ+Y1G/VNPlXNumnPuKKALkLv/v41z7kSC043Xmdnp+R7EbEBoPFf6unXrwoglIjFpxIjgisJWreDMM32nKRI/VanJmFYXk+QcN7+rexqKJKpwiqwMoHau7VrAz/tr7Jz7BGhoZlVD2z+HHtcC0whOP+a333DnXJpzLq1atWphxheRmLJ79+9XFA4e7OeKwmLy+LkDAej7+SSqb9EvjiKJKJwpHOYBjc2sPrAa6A70zN3AzBoB3zvnnJmdCKQC682sLJDknNsaWm8P3BfRTyAiseOll4Lbzhx7LJx/vu80RWrJEY2Z3rw9nb55h+s+eJG7fAcSCVO5crB5s+8Uh65cOd8JwiiynHNZZnY98DaQDPzXObfIzAaGXh8GXAxcbmZ7gJ1At1DBdTgwLTTZXwowzjk3s4g+i4hEs+zs33ux7rgjtmd3D9MjHa6l0zfv8NcPR/OQ7zAiYWrSxHeC+BHWZKTOuRnAjDzPDcu1/jDwcD77rQCaH2JGEYkHM2YEt5upWxcuvdR3mmIxq1ELPm/UgjbL53G17zAiUuzi/1dJEYkOjz8ePN5wQ3CPvwTxyLnXAnArwJ49XrOISPFSkSUiRW/+fPjww2CQRL9+vtMUqzePa8fiGo2oCzBxou84IlKMVGSJSNHLmRfr6quhQgW/WYqZS0riifbXBBtPPQVOczWLJAoVWSJStH75BcaPD6ZruPFG32m8GHdyV9YDpKfDnDm+44hIMVGRJSJF69lng7FIXbpAgwa+03ixK7U0I3I2cu7ZKCJxT0WWiBSdnTvhueeC9Vtv9ZvFs+cgmLbilVeC3j0RiXsqskSk6EyYENxCJy0N2rTxncarVQCdOwe9es8/7zuOiBQDFVkiUnSefTZ4vO66uL6FTthyxqQ9/zxkZvrNIiJFTkWWiBSNefOCgd6VKkG3br7TRIe2beGYY+DXX2HyZN9pRKSIqcgSkaKR04t11VVQurTfLNHCLJiMFeDpp/1mEZEipyJLRCJv/fpgPBbAwIF+s0SbXr2C3r05c4LePhGJWyqyRCTyRo2CXbvg3HOhUSPfaaJL2bJB7x7AsGEHbisiMU1FlohEVnb278XDtdf6zRKtBgwIHidMgM2b/WYRkSKjIktEIuu992D5cqhTB84/33ea6NSkCZx5JuzYAWPH+k4jIkVERZaIRFbOgPdrroHkZL9ZollOb9bzz+t+hiJxSkWWiERORga8/jqUKBHcDFr2r2tXqFoVFiyAuXN9pxGRIqAiS0QiZ9SoYExWly5w+OG+00S3kiWhb99gffhwr1FEpGioyBKRyMjOhhdeCNb79fObJVb07x88agC8SFxSkSUikfHBB/DDD1C3LrRr5ztNbNAAeJG4piJLRCJj5Mjg8aqrIEn/tIRNA+BF4pb+JRSRQ/fbbzBtWnDbmCuv9J0mtmgAvEjcUpElIoduzBjIzIQOHaB2bd9pYkvuAfAjRniNIiKRpSJLRA6Nc7+fKtSA94OTM93FxImwfbvfLCISMWEVWWbWwcyWmtlyMxuUz+udzWyBmc03s3QzOzXcfUUkxs2ZA4sWQfXqcMEFvtPEpqOOglatYNs2mDLFdxoRiZACiywzSwaGAh2BZkAPM2uWp9n7QHPn3PHAVcDIQuwrIlGsRu26mNl+l5GtWwPwyNq1WMmSB2wbzpKwcsayvfii3xwiEjEpYbRpCSx3zq0AMLMJQGfgu5wGzrltudqXBVy4+4pIdFuTsYrkEavzfa3crm10//sJsHsHo+77mOQjGh3y8fb2r3nI7xGTunWDm2+Gjz6CFSugQQPfiUTkEIVzurAm8FOu7YzQc39gZl3NbAnwJkFvVtj7ikhsuiT9Dcrt3sHnjVqwNAIFVkKrUAEuuihYHz3abxYRiYhwiqz8+u//NJmLc26ac+4ooAtwf2H2BTCzAaHxXOnr1q0LI5aI+Hb5rFcAGHVKN89J4kTOKcPRo4MZ9EUkpoVTZGUAua/JrgX8vL/GzrlPgIZmVrUw+zrnhjvn0pxzadWqVQsjloj4VH/dj5y+bDY7UksxOU0D3iPizDODGfN//DE4bSgiMS2cImse0NjM6ptZKtAdmJ67gZk1stCIVTM7EUgF1oezr4jEpj6zJgMw7YTz2Fq6vOc0cSIpCa64IljXAHiRmFdgkeWcywKuB94GFgOTnHOLzGygmQ0MNbsY+NbM5hNcTdjNBfLdtyg+iIgUH8vOpneoyHr5lEs9p4kzOUXWlCm6abRIjAvn6kKcczOAGXmeG5Zr/WHg4XD3FZHYduryuTT4bRU/VTqCD45q4ztOfGnQANq2hY8/hkmToH9/34lE5CBpxncRKbTLv5gEwNhWF5OdlOw5TRzSnFkicUFFlogUSpndO7gk/Q0AXtKpwqJxySVQrhzMmgVLlvhOIyIHSUWWiBRKl6/fovzu7cxucCLLamhurCJRtixcdlmwrjmzRGKWiiwRKZTLvwjmxnqptXqxilTfvsHj2LGaM0skRqnIEpGw1dqwmrOWfMaulJJMatHJd5z41qYN1KsHP/0UDIIXkZijIktEwtZ71hSSnGP68e3ZVLai7zjxLSkJevUK1seM8ZtFRA6KiiwRCY9z9AndRkdzYxWT3r2Dx8mTYedOv1lEpNBUZIlIWE5e8RVN16zglwrVeadZW99xEsNRR0FaGmzZAq+/7juNiBSSiiwRCUvPOVMBmNCyC3uTw5rHWCKhT5/g8eWX/eYQkUJTkSUiBUrJ2kO3ea8BMKbVxZ7TJJju3SE5GWbOhHXrfKcRkUJQkSUiBTp30UdU3baRb49syje1/+I7TmKpXh3OPReysmDiRN9pRKQQVGSJSIFyThWOa3URmHlOk4B0ylAkJqnIEpEDKg90mv8OAONbdvUbJlF16gTly8PcubBsme80IhImFVkickAXAaX37OLjJq35qUpN33FiV0pJzOzglrJleXHrVgDua9o0rH1q1K7r+QOLiC4REpEDCs3UxNhWF3nNEfOydpM8YvVB7z5u8Wdc+Xg3eletw30PflHgads1/VUQi/imniwR2b/VqzkL2J2SypSTzvedJqF93LQ1GRVr0OC3VZzyfbrvOCISBhVZIrJ/48eTBLxxXDs2l6ngO01Cy05KZvzJQW9ir1mTPacRkXCoyBKR/QvdM2+s5saKCjmnbC9Nf4PUPbs9pxGRgqjIEpH8LVwI33zDBmDmMWf6TiPAt7WOZn6tZlTesYmOCz/wHUdECqAiS0TyN3YsAJOAzBIl/WaRfca2DnoVe8+e4jmJiBRERZaI/Fl29r4ia4znKPJHE1p2Ya8lcd7C96m0faPvOCJyACqyROTPPvkEMjKgXj2+8J1F/uCXijX44OhTKZmVySXpb/iOIyIHoCJLRP4sNOCdXr1wfpNIPsaGrjLUKUOR6KYiS0T+aNcueOWVYL1XL79ZJF/TTjyP7amlabN8HvXX/eg7jojsR1hFlpl1MLOlZrbczAbl83ovM1sQWr4ws+a5XvvBzBaa2Xwz0wx6ItHujTdgyxY46SQ4+mjfaSQf20uV5bXjOwDQY840z2lEZH8KLLLMLBkYCnQEmgE9zKxZnmYrgbbOueOA+4HheV4/0zl3vHMuLQKZRaQo5Zwq7N37wO3Eq3GhObN6zJkGTid1RaJROD1ZLYHlzrkVzrlMYALQOXcD59wXzrmcy1xmA7UiG1NEisX69TBjBiQlQffuvtPIAbzb7HTWlq/C0b8u58RVC33HEZF8hFNk1QR+yrWdEXpuf64G3sq17YB3zOxLMxtQ+IgiUmxeeQX27IFzzoEaNXynkQPYm5zCxBbB77s9Z0/1nEZE8hNOkZXfrd7z7Zs2szMJiqzbcz3dxjl3IsHpxuvM7PT97DvAzNLNLH3dunVhxBKRiMt1VaFEv/EndwWg27zXSMre6zmNiOQVTpGVAdTOtV0L+DlvIzM7DhgJdHbOrc953jn3c+hxLTCN4PTjnzjnhjvn0pxzadWqVQv/E4hIZKxYAZ9/DmXKQNeuvtNIGObWP4Hl1epxxOa1nLXkc99xRCSPcIqseUBjM6tvZqlAd2B67gZmVgeYCvRxzi3L9XxZMyufsw60B76NVHgRiaBx44LHLl2gXDm/WSQ8Zr8PgNcpQ5GoU2CR5ZzLAq4H3gYWA5Occ4vMbKCZDQw1uxuoAjybZ6qGw4HPzOwbYC7wpnNuZsQ/hYgcGud0VWGMGhc6Zdj167covXun5zQikltKOI2cczOAGXmeG5ZrvR/QL5/9VgDN8z4vIlHmq69g6VKoVi0Y9C4xY/nhDZhb/wRarvyaCxa8wystOhe8k4gUC834LiK/92L16AEpYf3uJVEkpzerpyYmFYkqKrJEEl1WFowfH6zrVGFMmtSiE1lJyXT49kOqbN3gO46IhKjIEkl0778Pa9ZAkyaQppsyxKK1h1XjvWanU2JvFpd8+brvOCISoiJLJNHlHvBu+U2LJ7FApwxFoo+KLJFEtm0bTA1d+q8JSGPaa8d3YHtqadosn0e9dat8xxERVGSJJLbXXoMdO6B1a2jQwHcaOQTbS5XlteM7ANBjrnqzRKKBiiyRRKa5seLK+FahU4aamFQkKqjIEklUa9bAO+8EUzZcdpnvNBIB7x59OmvLV+HoX5dzgu8wIqIiSyRhTZgA2dnQsSNUreo7jURAVkoJJqV1AkAj7ET8U5ElkqjGjg0edaowrowPXWXYA2DvXq9ZRBKdiiyRRLR0KcybB+XLw4UX+k4jETSnwYksr1aPIwE+/NB3HJGEpiJLJBHl9GJdcgmULu03i0SW2b7erH3fs4h4oSJLJNE4p6sK41zOxKRMmQI7d/oNI5LAVGSJJJpZs2DlSqhZE9q29Z1GisD/ajRkLsDWrfC6brMj4ouKLJFEk9OL1aMHJCf7zSJFZt+JQp0yFPFGRZZIIsnMhIkTg3WdKoxrEyAoomfMgPXrfccRSUgqskQSycyZsGEDHHMMHHec7zRShNYCtGsHWVnwyiu+44gkJBVZIokk94B3M79ZpOjl9FbqlKGIFyqyRBLF5s2/D4Lu2dNvFikeXbpAmTLw2Wfwww++04gkHBVZIoli6lTYtQvOOANq1/adRopDuXLQuXOwrt4skWKnIkskUWhurMSU832PGRPMkSYixUZFlkgiyMgIbrFSsiRcfLHvNFKc2reH6tVhyRJIT/edRiShqMgSSQTjxwe9GBdcABUr+k4jxSklBXr1CtZfeslvFpEEE1aRZWYdzGypmS03s0H5vN7LzBaEli/MrHm4+4pIMdCpwsR2+eXB4/jxwVxpIlIsCiyyzCwZGAp0BJoBPcysWZ5mK4G2zrnjgPuB4YXYV0SK0oIFwVKpEnTs6DuN+NC8ORx7bDAp6Vtv+U4jkjDC6clqCSx3zq1wzmUSTCTcOXcD59wXzrmNoc3ZQK1w9xWRwqtRuy5mFtbycPOgY3nYxo1YqVJh75ezSBww+703S6cMRYpNShhtagI/5drOAE4+QPurgZxflQq7r4iEYU3GKpJHrC6wnWVn03NQS9j4C+Nvm0pW90gAACAASURBVEZy45aFPtbe/jUPJqJEm5494fbbg7nSNmyAypV9JxKJe+H0ZOX3q2y+1wGb2ZkERdbtB7HvADNLN7P0devWhRFLRApy+v9mU3vjL6ysUpsvGqb5jiM+HXkknHMO7Nnz+/0rRaRIhVNkZQC5Zy6sBfyct5GZHQeMBDo759YXZl8A59xw51yacy6tWrVq4WQXkQL0nD0VgPEnd8Ul6WLihJdzynD0aL85RBJEOP/qzgMam1l9M0sFugPTczcwszrAVKCPc25ZYfYVkaJRevdOLvnyDQDGtbrIcxqJCl26QPnyMGcOLF3qO41I3CuwyHLOZQHXA28Di4FJzrlFZjbQzAaGmt0NVAGeNbP5ZpZ+oH2L4HOISB6d58+kws6tzK1/AkuOaOw7jkSDMmXg0kuD9Zdf9ptFJAGEdf7AOTfDOdfEOdfQOfdA6LlhzrlhofV+zrlKzrnjQ0vagfYVkaLX9/Ng3M3oUy71nESiSs4pw5dfhuxsv1lE4pwGaYjEodrrV3PWks/YnZLKxBaaNUVyOe00qFsXVq2CTz7xnUYkrqnIEolDvWdPJsk5Xjv+XDaV1W10JJekJOjTJ1jXnFkiRUpFlki8cY7Lv5gEwEunXOY5jESlnCLrlVdgxw6/WUTimIoskThzyvfpNF77A6sr1uCdv7T1HUd8SSm5/5n8mzZlFsC2bfQsW7bQdwHIu9SoXdf3pxWJSuHM+C4iMSSnF2tsq4vJTkr2nEa8ydp9wLsCjPloNK3HDqZvs7ZMumXcIR1qje4KIJIv9WSJxJHSu3dy2bxgKrqXdFWhHMCkFp3YnZJKu8WfUGtDwbdoEpHCU5ElEke6fj2Dw3ZtY47mxpICbCxbiVdP6ECSc1z+xSu+44jEJRVZInHkitCpwtEa8C5hGNWmOxD8uTHNmSUScSqyROJEnfUZnLnkc3allGRSi06+40gMeP/oU1lV+UgarvuRtstm+Y4jEndUZInEid6zgrmxXj2hg+bGkrBkJyXvm+aj7+cTPKcRiT8qskTiQa5xNZobSwoj59TyxV/O4LAdWzynEYkvKrJE4kCb5fNotC6YG+u9Zqf5jiMxZGW1unxwVBtK79lFt3mv+Y4jEldUZInEgStCN4Me0/oSzY0lhfZiaAD8laE/RyISGSqyRGJc+Z1b9/VAjNKpQjkI007syKbSh9Fy5df8ZfUS33FE4oaKLJEY12PONMpm7uTjJq35X42GvuNIDNqVWpoJLTsD0Fe9WSIRoyJLJMZd/WlwS5SRp/X0nERiWc6cWb1mT6FEVqbnNCLxQUWWSAw74ceFnLRqIRvKVGTqSef5jiMxLL1ecxbWPIrqW9dz/oL3fMcRiQsqskRi2NWfjgVgTOuL2V2ilOc0EtPMGNWmGwBXfTbecxiR+KAiSyRGldm9g55zpgHwgk4VSgSMbXUxu1NSOffbj6i9XjeNFjlUKrJEYtSl6a9z2K5tzG5wIotqHuU7jsSB38pXYdoJHUl22Vz92TjfcURinooskRjVL3SqcORpvTwnkXgyom1vAK78bALJe7M8pxGJbSqyRGJQM6D191+ypVQ53QxaIurjJq1ZUqMhNTf9qgHwIodIRZZIDOofehzfsgs7SpbxmkXijNm+3tEBn4zxHEYktqnIEok1u3bRJ7T6wuk6VSiR99Ipl7IrpSTtF31EvXWrfMcRiVlhFVlm1sHMlprZcjMblM/rR5nZLDPbbWZ/z/PaD2a20Mzmm1l6pIKLJKypU6kCfFnnWL6qe5zvNBKHNpSrzOS080lybt/YPxEpvAKLLDNLBoYCHQmGgvQws2Z5mm0AbgQe3c/bnOmcO945l3YoYUUEGDYM0LQNUrRGnB70l/b9fCIpWXs8pxGJTeH0ZLUEljvnVjjnMoEJQOfcDZxza51z8wD9TRQpSgsXwqefsgUY1+oi32kkjn3eqAWLjmhCjS3r6PTN277jiMSkcIqsmsBPubYzQs+FywHvmNmXZjZgf43MbICZpZtZ+rp16wrx9iIJ5NlnAXgJ2FaqnN8sEt/M9k3nMOBjDYAXORjhFFmWz3OuEMdo45w7keB043Vmdnp+jZxzw51zac65tGrVqhXi7UUSxObN8PLLADznOYokhjGtLmZniVK0W/wpDdeu9B1HJOaEU2RlALVzbdcCfg73AM65n0OPa4FpBKcfRaSwXn4Ztm+HM87gO99ZJCFsKluRiaF52K75+GXPaURiTzhF1jygsZnVN7NUoDswPZw3N7OyZlY+Zx1oD3x7sGFFEpZz+04Vcu21frNIQnnuzL5AMAN8md07/IYRiTEFFlnOuSzgeuBtYDEwyTm3yMwGmtlAADOrYWYZwK3AXWaWYWaHAYcDn5nZN8Bc4E3n3Myi+jAicevDD2HxYjjySOjSxXcaSSBf1mvO7AYnUmnHZnrOmeo7jkhMSQmnkXNuBjAjz3PDcq3/SnAaMa8tQPNDCSgi/N6LNWAAlCjhN4sknGfOuopWK77iug9eDGaDt/yG6opIXprxXSTaZWTAq69CSgr0719we5EIm3LS+fxSoTrHrl5C22WzfMcRiRkqskSi3fDhsHcvdO0anC4UKWZ7UlIZcXowncN1H7zoOY1I7FCRJRLNdu8OiizQgHfxavjpvdmTnELnr2dSe/1q33FEYoKKLJFoNn48rFkDxx0Hbdv6TiMJ7NeKhzPlxPNJdtkM/Pgl33FEYoKKLJFo5Rw88USwfsstGmws3g0960oArv50LKUyd3pOIxL9VGSJRKsPP4QFC+Dww6FHD99pRJjVMI0v6xxL1W0b6TYvrOkSRRKaiiyRaJXTi3XttVCypN8sIgBmPHP2VQDc9O6IoLdVRPZLRZZINFq2DN54IyiuBg70nUZkn4ktOvNLheoct3ox7RZ/6juOSFRTkSUSjZ56Knjs3RuqV/ebRSSXzBIlGXpmMDbr1reHFdBaJLGpyBKJNhs2wKhRwfrNN3uNIpKf58/ow/bU0rT/7mOOzdDtykX2R0WWSLQZMQJ27ID27eGYY3ynEfmTjWUr8WKb7gDc/M5wz2lEopeKLJFosmcPPPNMsH7LLX6ziBzA0+36sdeS6DH3VY7wHUYkSqnIEokm48YF9yo8+mg491zfaUT2a0X1ekw7sSOpe/dwg+8wIlFKRZZItMjOhocfDtZvv12Tj0rUe+KcawAYCLBtm9csItFIRZZItJg+HRYvhtq1oWdP32lECjSn4Ul83qgFlQBGjvQdRyTqqMgSiQbOwUMPBet//zuUKOE3j0iYHjk3dOPyRx8NbmguIvuoyBKJBh99BHPnQtWq0K+f7zQiYXvzuHYsBFi9Gl7SjaNFclORJRINcnqxbrwRypTxm0WkEFxSEg/mbAwZAllZPuOIRBUVWSK+ffklvPsulCsH113nO41IoU0CaNQIVqyAiRN9xxGJGiqyRHwbMiR4vOYaqFzZbxaRg5ANMGhQsPHgg8GVsiKiIkvEq0WLYMoUSE3V5KMS2/r0Ca6M/e47eO0132lEooKKLBGf7r03uLKwXz+oWdN3GpGDl5oKt90WrD/wQPDnWiTBqcgS8WXhQnjlleA/pzvu8J1G5NBdfTUcfngwznDmTN9pRLwLq8gysw5mttTMlpvZoHxeP8rMZpnZbjP7e2H2FUlY994bPA4YALVq+c0iEgmlSwfzvAHcfbd6syThFVhkmVkyMBToCDQDephZszzNNgA3Ao8exL4iiWfBgmAsVsmS6sWS+HLttVCjBqSnB3cxEElg4fRktQSWO+dWOOcygQlA59wNnHNrnXPzgD2F3VckIeX0Yl1zDRx5pN8sIpFUpgwMHhys/+MfutJQElo4RVZN4Kdc2xmh58JxKPuKxIwatetiZmEtx5vB1KnsBI54+umw98u9iES1/v2DU+ALF8Lkyb7TiHiTEkab/P5FD/dEe9j7mtkAYABAnTp1wnx7keiwJmMVySNWh9X23qFXwfy3Gd6uH+u63UvyQRxvb3/9riJRJKXkn4r//sBwYHG3bhzTrRuR7M86vFYdfv3pxwi+o0jRCKfIygBq59quBfwc5vuHva9zbjjB30nS0tI0WlLi0inL59F5/ttsTy3NvztodneJE1m7//RLxktZexj0j9M5+rdV9L7qKca2viRih1ujXzIkRoRzunAe0NjM6ptZKtAdCHc046HsKxJfnOOhKQ8A8ET7a1hTobrnQCJFJyulBPdfGEyw+8/pj5G6Z7fnRCLFr8AiyzmXBVwPvA0sBiY55xaZ2UAzGwhgZjXMLAO4FbjLzDLM7LD97VtUH0Ykml3wzbu0WT6PdeUq81j7gb7jiBS5cSdfxLdHNqXBb6v460ejfccRKXZhzZPlnJvhnGvinGvonHsg9Nww59yw0PqvzrlazrnDnHMVQ+tb9revSKJJ3pvFA1MfAuCBC25ma+nynhOJFL29ySkMvji40vDON56i4vZNnhOJFC/N+C5SDC7/4hX+8ssyVlStw/DTe/uOI1JsZhx7Nh82PYXKOzZxx4z/+I4jUqxUZIkUsdK7d/LP6cE8vf/ocjuZJUp6TiRSjMy4/dJ/AHD9B/+l7m8/FbCDSPxQkSVSxP72zjBqbfqVr+ocw6QWnXzHESl2X9U9jnEtu1IyK5P7Xv237zgixUZFlkgRqr1+NbfNfAaAv3W7F5ekv3KSmP7R9XZ2p6TSa85U0lbO9x1HpFjoX3yRIvTQlAcok7mLSWkX8mmTVr7jiHjzY9XaPH12PwCeHn8XptvtSAJQkSVSRE5dNofu815jZ4lS3H7JP3zHEfHugQtu4ucKh9Ny5ddcPusV33FEipyKLJEikJS9l8cn3g3AIx2u5acqmqFaZFupctx+yV1A0MtbYcdmz4lEipaKLJEicOVnEzhx1besqnwkj5x7re84IlFj/Mld+axRS6pvXc/d0x/zHUekSKnIEomwalt+46EpDwJw+yV3sbNkac+JRKKIGTf2/Bd7LYnrPhzFX1Yv8Z1IpMioyBKJsEdeuZfKOzbxTrO2vJKmKRtE8lpQ+y8837YPKdl7eWbsYA2Cl7ilIkskgs7+7hN6z57KzhKluL7Xg2DmO5JIVLq7y22sKV+V0/43h6s/G+c7jkiRUJElEiGlMncydMwdAPzrgptZUb2e30AiUWxT2Yrc1ON+AB6e/C+O2PSr50QikaciSyRCBr/5NI3W/cDCmkfxWPuBvuOIRL3JaRfyevNzqLBzK0+N1zQnEn9UZIlEwHHA399+DoBrew8hK6WE30AiscCMG3o+wNaSZbnoqxl0+WqG70QiEaUiS+RQZWbyEpC6dw/PnnEFsxq18J1IJGZkVK7J4IsHA/D0uLuouH2T50QikaMiS+RQ3XcfzYHl1epxx8V3+k4jEnOGtb2cLxqmceTmNfxnnP4OSfxQkSVyKObMgYceIhu4+son2F6qrO9EIjHHJSVx5VVPsj21ND3mvsplc1/zHUkkIlRkiRysnTvhiisgO5vHgM8bt/SdSCRmfV+9Pn/rdg8AQ8fewZEbf/EbSCQCVGSJHKzbboOlS+Hoo9F1USKHbuRpvXjz2LOptGMzL4y6VZOUSsxTkSVyMKZNg2eegRIl4OWX2e07j0g8MGPAFY+yrlxlzvnuE2744AXfiUQOiYoskcL68Ue46qpg/d//hpNO8ptHJI6sqVCday5/BIAhkx+gxcqvPScSOXgqskQKY88e6N4dNm2CTp3gppt8JxKJO9NP6MB/zrqK1L17GP/8QE3rIDFLRZZIYdx1F8yeDbVrw4sv6t6EIkXk9kvuYl6946m3PoMXRt0KzvmOJFJoKrJEwjVlSnB6MDkZxo+HypV9JxKJW5klStL9mmFsLFOBzvPf5tZ3nvcdSaTQwiqyzKyDmS01s+VmNiif183Mng69vsDMTsz12g9mttDM5ptZeiTDixSbBQuC6RoAHn4Y2rTxm0ckAfxYtTZXXvkkAA9NeYD2337kN5BIIRVYZJlZMjAU6Ag0A3qYWbM8zToCjUPLAOC5PK+f6Zw73jmXduiRRYrZb79B586wfTv06QO33uo7kUjCeOP49tx34a0ku2zGDf8rjX/93nckkbCF05PVEljunFvhnMsEJgCd87TpDLzkArOBimZ2RISzihS/PXvgssvghx+gRQt4/nmNwxIpZvdfcAvTTuhIxZ1bmDb0Sir4DiQSpnCKrJrAT7m2M0LPhdvGAe+Y2ZdmNmB/BzGzAWaWbmbp69atCyOWSBFzDq6/Hj78EGrUCObGKl3adyqRhOOSkuh71VMsqHk0R/36PZMAMjN9xxIpUDhFVn6/tue9zONAbdo4504kOKV4nZmdnt9BnHPDnXNpzrm0atWqhRFLpIjddx8MHw6lSsHUqVAz7+8WIlJctpcqy0XX/Ze15avQHqBfP11xKFEvnCIrA6ida7sW8HO4bZxzOY9rgWkEpx9FotuIEXDPPZCUBBMmQOvWvhOJJLwfqtXhwhtfZhvAyy/D4MG+I4kcUDhF1jygsZnVN7NUoDswPU+b6cDloasMWwGbnXO/mFlZMysPYGZlgfbAtxHMLxJ5r78OAwcG688+Gwx6F5Go8GW95lwCwVQqQ4bAf/7jO5LIfhVYZDnnsoDrgbeBxcAk59wiMxtoZqH/iZgBrACWAyOAa0PPHw58ZmbfAHOBN51zMyP8GUTyVaN2XcysUEs7M3Z26gTZ2dwH2MCBYe0nIsXn7ZSS9N27N9i48UauLuTf88IsNWrX9fthJaalhNPIOTeDoJDK/dywXOsOuC6f/VYAzQ8xo8hBWZOxiuQRq8Nuf+biz3jtmSsonbmL50/vzX29h5AcZgG1t7/Ga4kUm6zdjBmxmsrvjeDxifcwwoy9fZ9gzCmXRvxQa/R3Ww6BZnwXAU5fOotXn+lLmcxdvHBqD67v9ZCmahCJck+368+gi+8kyTleGHUr3edM8x1J5A9UZEnCa//tR7z+dB/KZu5k1CmXMbDPv3FJ+qshEgse7XAtd3f+P5JdNi+9cANXfzLWdySRffQ/iSS0y+a+xqvP9N1XYA244lEVWCIx5sELbuYfXW4jyTmef/k2/vZ23puOiPih/00kYQ38cBRjRl5H6t49PHHOAPpf8RjZScm+Y4nIQXjo/Ju4oecDADw8+V88MPUhLDvbcypJdCqyJOFYdjb/mvoQz4wLxnIMvugO/u/Su9WDJRLjnjuzL32u/g9ZScnc/tYzjB1xLaUyd/qOJQksrKsLReJFuV3beHnk9Vz4zbtkJSVzXa+HeOH0Xr5jiUiEjG91ERvKVWL88wO5LP116q7PoOv1L7L2MN1JRIqffnWXhFF/3Y989lAnLvzmXTaUqcj5N41RgSUSh94+5kxOG/QaP1Spxckrv2b2A+fRcsVXvmNJAlKRJQmh09czmfOv8zjm56V8d0RjWt/5Bu83y/c2miISBxbVPIrWg99kVsOTqLPhZz7690Xc8N5I3e9QipWKLIlrJffs4qlxdzL12aupvGMTrzc/h1MHTef76vV9RxORIrbusKqc9ffJPHX21cEFLhP/ycRhA6iwY7PvaJIgVGRJ3PoL8PmDF3Ldh6PITC7BLd3upet1L7KlzGG+o4lIMdmTksrfut/HZQOfZ0upclz81Qzm33M25yz62Hc0SQAqsiT+ZGbCPffwFXB8xncsr1aPU++Yzn/a9dMs7iIJaupJF9DyrreYW/8Eam/8hbee7MkzYwZRdtd239EkjqnIkvgydy6kpcG995IKDGvbhxb/mMlXdY/znUxEPFt+eANOu/1V7upyO5nJJRj48ct8c89ZdPp6psZqSZFQkSXx4ZdfoG9fOPlkWLgQGjbkDOD63kPYWrq853AiEi32Jqcw5PwbOfmuGXxV5xjqrc9g6rNX88bTfWi0ZoXveBJnVGRJbNuxAx56CJo0gdGjITUVbr8dFixAIy5EZH8W1mpG68FvcmOPf7Gp9GF0+PZDvrnnbB6b+E+qbl3vO57ECRVZEpt27oSnnoIGDWDwYNi2Dbp0ge++gyFDoEwZ3wlFJMrtTU7h2bOu5Oh/fcqLbbpRMiuTm94byf/uaM0/X3uU8ju3+o4oMU5FlsSWLVvgySehUSO4+WZYsyYYg/XuuzBtGjRs6DuhiMSYdYdVpX/fxznx7nd489izKb97O/944wlW3t6SByEYjiByEFRkSbGqUbsuZlbopb4ZT5ixpUIFuOUW+PlnvgYuBCw9HTvnnD/tIyJSGAtq/4XON75E29um8XGT1lTcuYU7AOrVg3794NtvPSeUWKN7F0qxWpOxiuQRq8Nqm7pnNxcseJcrPp9Eh28/JNllA/Bxk9Y8cc4A3jyuHS4pieT97L+3f80IpRaRRPJ545ac/X+TOfn7L7l1SCcu3rMHXnghWFq1gv79oVs3KFvWd1SJciqyJKpYdjatV3zJpfOm03POVKps3wRAZnIJxrXoytPt+vN13WM9pxSRRDCn4UlcArglS4IxoGPGwOzZwXLzzdC1K1x2GZxzTnDRjUgeKrLEuxJZmbRdNpsuX82g8/y3OWLz2n2vza/VjNFtujH+5K78Vr6Kx5QikrCaNIGhQ+Hf/4ZXXoERI+CLL+Cll4KlYkXo3Dkous46C8pr2hgJqMiS4uccx65ezNnffcrZiz/j9GWzKJu5c9/LK6vU5tUTOzK21cXMr3OMx6AiIrmULRvMx9e3LyxbFhRckybBggXBFDKjR0NKCrRpAx06QLt2cPzxwXOSkPTNS9HbsgXmzYM5c3gNaPW35lTPMw/Nt0c25bUTOjDtxPOYX/svuv2NiES3Jk3gzjuDZelSmDwZ3noLZs2Cjz8OljvuCAqzk08OCq82beCkk6BqVd/ppZioyJLIycqC5cuDK3ByL8uW7btlRSeArevJqFiD948+jfebncYHR53KrxUP9xpdROSgNW36e8G1cSO8/z7MnBkUWsuXwwcfBEuOI4+E5s1/X445JpjzT/P7xZ2wiiwz6wA8BSQDI51zQ/K8bqHXzwN2AH2dc1+Fs6/4VaN2XdZkrAqrbQpwOHAEUBdoANQPPTYIPZff0M9M4GtgNjAHSH/gc1ZUq6veKhGJP5UqwSWXBAsEc/l98QV8/nnQy/XNN/Dzz8Hy1lt/3LdWrWAOwJylTh2oWTNYjjwSSpcu/s8jh6TAIsvMkoGhwDlABjDPzKY7577L1awj0Di0nAw8B5wc5r7iy7p1HJaxioaDplNpx2Yq7thMpR2bqbx9I5W2B481Nq+jxua1HLF5DdW2bSjwLVdWqc2imk1ZVPMovj2yKYtqNmVJjUZkligJBNMqJFevV8QfTEQkShx+eDAgvmvXYDs7G77/Pii2cpYlS2DlSsjICJaPPsr/vSpVCgquI46AKlWCpXLlYMm9XqkSlCsXnKosVw5KlCi2jyt/FE5PVktguXNuBYCZTQA6A7kLpc7AS845B8w2s4pmdgRQL4x9xZeGDVkGMKRTWM33WhJrD6vKLxWqs7rSEayoVpcVVevyQ9XarKhWlx+q1GZnSf2mJSJxJKVksUxunAzUAZqUSKXBnkwaAbWAmqHlSCB148bgdGQhJ0XNBLYD20KPOetZpUpx9tatGphfhML5ydYEfsq1nUHQW1VQm5r7eT7vvuJLw4Ysnz+fjfWOZ2OZCmwsW4ENZSuyqUwFNpapyIayFfm1QnV+rVCdXypUZ135KmQn7W/qTxGROJS1O+wJlA/VKmBl/5r5Hs+ys6m6bQNHbvqVIzavpdL2TVTet2zct15l+0Yqbt9E2cydlNu9nbK7d5CavZdUoFKe98zatQuS9W96UTIXGpC83wZmlwLnOuf6hbb7AC2dczfkavMm8JBz7rPQ9vvAbQRDdQ64b673GAAMCG02BZYe4meLVlWB33yHkIOm7y+26fuLXfruYlu8f391nXPV8j4ZTk9WBlA713Yt4Ocw26SGsS8AzrnhwPAw8sQ0M0t3zqX5ziEHR99fbNP3F7v03cW2RP3+wrlB9DygsZnVN7NUoDswPU+b6cDlFmgFbHbO/RLmviIiIiJxp8CeLOdclpldD7xNMDbvv865RWY2MPT6MGAGwfQNywmmcLjyQPsWyScRERERiSJhXVLgnJtBUEjlfm5YrnUHXBfuvgku7k+Jxjl9f7FN31/s0ncX2xLy+ytw4LuIiIiIFF44Y7JEREREpJBUZHlkZn83M2dmultoDDGzR8xsiZktMLNpZlbRdyY5MDPrYGZLzWy5mQ3ynUfCZ2a1zexDM1tsZovM7CbfmaRwzCzZzL42szd8ZyluKrI8MbPaBLcbCu/GgRJN3gWOcc4dBywD7vCcRw4g1+29OgLNgB5m1sxvKimELOBvzrmjgVbAdfr+Ys5NwGLfIXxQkeXPEwQTtmpQXIxxzr3jnMsKbc4mmP9Note+W4M55zKBnNt7SQxwzv3inPsqtL6V4D/rmn5TSbjMrBZwPjDSdxYfVGR5YGadgNXOuW98Z5FDdhXwlu8QckD7u+2XxBgzqwecAMzxm0QK4UmCDoVs30F80F0hi4iZvQfUyOelO4HBQPviTSSFcaDvzzn3WqjNnQSnMsYWZzYptPzu7qse5BhjZuWAKcDNzrktvvNIwczsAmCtc+5LMzvDdx4fVGQVEedcu/yeN7NjgfrAN6E7u9cCvjKzls79f3t3jNJAFEVh+D+gYKOrsUshmEYkuABR3IKF67CysBW0sXADbiJdWhtBENzCtZhJaYzCyyTh//qBA8PAmcu8ufWxwoha4Kf7N5fkCpgA4/I/KOtumdVgWmNJdukK1lNVvQydR0sbAWdJToE94CDJY1VdDJxrZfxP1sCSvAGHVbXNizO3SpIT4BY4qqrPofNosSQ7dAcUxsA73bqvc7dPbIZ0b6MPwFdVXQ+dR//TT7JuqmoydJZV8pss6e/ugH3gNck0yf1vF2g4/SGF+XqvGfBswdooI+ASOO6ft2k/GZHWnpMsSZKkBpxkSZIkNWDJkiRJasCSGcsKRwAAAClJREFUJUmS1IAlS5IkqQFLliRJUgOWLEmSpAYsWZIkSQ1YsiRJkhr4Bh2Zs4BbpPJYAAAAAElFTkSuQmCC\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)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Box plots" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Below we display the blox plot of the following data:\n", "\n", "$$\\{-30,-1, -5, -0.5, 0.5, 0.6, 0, 2, 3, 4.6, 4, 7, 18, 35\\}.$$" ] }, { "cell_type": "code", "execution_count": 1, "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAA1IAAAEvCAYAAABPBfAxAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjMsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+AADFEAAAgAElEQVR4nO3df3hcZZ338e99ZiaTH81MfjRpJgXEVZTKz7ZA5ccKIg8UaK3rurqwoiLbVite4IU8i4vPNiAil/xwEUXagrKLK8sPH6DbXWBZWOTiWWk3LRXqglIBhcxJOk0ySTrJpJk59/NHEhpi2uZuMuec3Hm/rmsuOuecmfl+c9+Bfjjn3KO01gIAAAAAmDwn6AIAAAAAYKYhSAEAAACAIYIUAAAAABgiSAEAAACAIYIUAAAAABgiSAEAAACAoehB9rM2OgAAAPyggi4AMMEZKQAAAAAwRJACAAAAAEMEKQAAAAAwRJACAAAAAEMEKQAAAAAwRJACAAAAAEMEKQAAAAAwRJACAAAAAEMEKQAAAAAwRJACAAAAAEMEKQAAAAAwRJACAAAAAEMEKQAAAAAwRJACAAAAAEMEKQAAAAAwRJACAAAAAEMEKQAAAAAwRJACAAAAAEMEKQAAAAAwRJACAAAAAEMEKQAAAAAwRJACAAAAAEMEKQAAAAAwFA26AAAAAABTN3++LBCR6qDrmAX62trkFYIUAAAAYIdqEckGXcQsUCPCpX0AAAAAYIwgBQAAAACGCFIAAAAAYIggBQAAAACGCFIAAAAA3pFOl/9rV9eKU4KuI+xYtQ8AAADAO5qb8xcGXcPBuO7cm7XuOVXEqxSJZqLRBXc3Nm5/SERkYOCBZDb7pRu17jldJNodiy28taFh86b9vZfp8aM4IwUAAABgRqmsvGhdQ8NLH21uLi6qqvrKlwuFHVd2d3/hGBGRbHbNWhFnb339M6eVl3/860NDrdf19Hz1/ft7L9PjRxGkAAAAAMul02XPdHR88LJ0umJjOu1sb29v/HYut77edZMb0unINtetvndg4OeJ0WM7Oy88ddxrvzj82shW163/3tDQS2XBdSOSTN6xMxZbMDT8zNEiShcKvz1icPD5Cq27zq2svOz2ePys/rq6h7c6Tu3T+fxjKyZ6H9PjxwptkGppaQm6hJKwtS/4y9Z5ZGtfAACEgee9fV4yefulicSt53pe19k9PVfcXVHxqdsaGnYsERHV23v15/b32mLxrQsSiZv+urb2Z2drnTu6u/uznxy733Vr1qXTkdaJHq5bs64U/bS3p9am086vcrlbnxCJZhKJG38xMPCTI0WUl0x+983R45Sqf9Xzeo+a6D1Mjx8rtEHquuuuC7qEkrC1L/jL1nlka18AAIRBNHrsfVVVqzrnzPnaLqWqWx0n+auamnteicUWDEUiRzzled0L9vfaWOzYf5wz54pdFRWf6XGchmc8b9e7jk2lsqubm4snTfRIpbKrS9FPU5N7XVNTZmF5+Z9fHIk0/3sk8t69npetFHH6xh6nVMUekWLVRO9hevxYLDYBAAAAzAKOM3f3vmeRvFIVnfuexwZFCvsND47TkBn9s1LRvOcNNE5XXa5bfZ/WeyZcJVCpqq2p1J6L919XvVdX9/DW9vamj3d1nX9RWdmprSLenHcfla8SieQmfn1Nv8nxYxGkAAAAAEyJ6yY3aN130kT7lKpuTaV6Vu7vtalU3yVTr0BHPK/niIqKSx/u7/9xpLf3mvckEjf9XkTE8zqPdpzEaxO9qqLi0jdNjh8rtEFKax10CSVha1/wl63zyNa+AACw3YGC0nTL5e6p6++/+8OJxHeejcWOy3d3/9VpnpdZVlb2kavi8TMGlKp7KpfbcEU8vuzaXO72BZ7XdU5V1ZrPTPRepsePFdp7pNavXx90CSVha1/wl63zyNa+AADA9FEqoguFly/u7PzYc+3tja2Dg8/+TTR6/I1z5z77tIhITc2dLSJevLPzI7/M5x+9LRY7aW0yecfO0de7bnJDR8f737lv62DH77eOg/wf4MD+97BSysr/O21rX/CXrfPI1r4AAJOigi5gpps/X04RkWzQdcwCNW1tsiW0Z6QAAAAAIKwIUgAAAABgKLRBauPGjUGXUBK29gV/2TqPbO0LAADYJ7RBavHixUGXUBK29gV/2TqPbO0LAADYh8UmfGZrX/CXrfPI1r4AAJPCYhNTxGITvmGxCQAAAAA4FKH9Ql4AAAAAs8uuXcf+VaHw+idF8h90nPpNTU2Zaw50vOtW36d17kQRKQxviXU0Nw8u9aHU8AaplSt9+3JkX9naF/xl6zyytS8AADA5jlO3q6ys/keFwm/OECmWT+Y10ejx1zc2bn+o1LWNF9p7pAAAADCrcI/UFB3oHinP61OZzKJVxeKbF4no8mh0wfcLhR3X1tT85IzKyi90+1zqQXV0HHml1rmmyZyRikTet9HnIBXue6RsXb3L1r7gL1vnka19AQAQtEzm+Ms9L31WdXXLxXV1j55TLL7+CZFIdjIhynVr1qXTkdaJHq5bs86P+g+kUNhxVTodfcF1q+7v6lpxil+fG9pL+7Zt2xZ0CSVha1/wl63zyNa+AAAIUn//vbXF4h8ura5eu7y6+tq0iIjjND3reZmT8/lNc7q6Pn2vSP59lZVf/nRNzQ9fG//6VCq72veiJ6m8fPnNlZVf/J3jzN+bzV5yYT6/6a7e3m+uSCRueKvUnx3aIAUAAABg6nK5dacqFf9ddfXftY1u0zpfo1TiN9HocfmamjtX9vZe/Tel+OzhxSD2THiWSKmqranUnoun8v61tT97afTPDQ2tj7puYlk+/8iZicQNP53K+05GaINUKpUKuoSSsLUv+MvWeWRrXwAABEnrvlqRsq7R58WiG/G8zNmx2LF3RqPvKUSjX+ju7b16v6933eQGrftOmmifUtWtqVTPfleLSqX6LplS8caUFtG+3G8X2iCVTqeDLqEkbO0L/rJ1HtnaFwAAQYpE3vN6ofDKFX19LYdFo8f19vSsuVpk6PBo9Jg/uoxvIgcKStNtJORFRDxHRDtDQy+VOU5DMRJJFccfOzDwaHV//10nJBLf3eI4DcWurk9coHXfyfH4BTf6UWtoF5toaWkJuoSSsLUv+MvWeWRrXwAABKm+/l9/6TgNT/T1fWtjd/fFDzlO/W9FxKuqunxSQcpPu3eftiaTOeHlYvGt1Z7XuSKTOeHl3btPWzO633WTGzo63r9aRETrbHRw8LkrM5mFL3R0HLa5UNhxSTy+dE0yecsbftQa2uXPlVJykNpmJFv7gr9snUe29gUAmBSWP5+iAy1/PlZn5/mnDw4+s7a5efDc0W3t7Q03lZd/+p6JFpvAHwn38ucAAAAApl+x+Pv3KVX529Hnrptc73nZ0wcG7r0hk1n8Z0HWNpOE9h4pAAAAANPP87re5zg17wSpVKpnVZD1zFShDVKtra1Bl1AStvYFf9k6j2ztCwCAMGlqal8bdA024NI+AAAAADDEYhM+s7Uv+MvWeWRrXwCASWGxiSma7GITmDIWmwAAAACAQ0GQAgAAAABDoQ1Sa9faeQ+crX3BX7bOI1v7AgDAZu3tDTd1dBx5ZdB1+C2090gBAABgVuEeqSkK6h6p9vaGm5Sqap83782/n8r77Np17F8VCq9/UiT/Qcep39TUlLnGZP9Yrlt9n9a5E0WkMLwl1tHcPLh0KvWNEe57pJqbm4MuoSRs7Qv+snUe2doXAAA4OMep21VWdvKPHKfx4UPZP140evz1zc3ewuHHtIWoffVM9xtOF9d1gy6hJGztC/6ydR7Z2hcAAEFLp8ue6ej44GXpdMXGdNrZ3t7e+O1cbn296yY3pNORba5bfe/AwM8TIiJ79tze6Lq1d6TT0RfS6bKnd+360CVj3yubvWxBOl3xyPDr6r+ntRefjhrnzn3uqblzf/EfSpVPeFbtYPv9FtoghZmnrq5OlFJGD2lJGr9mso+6urqgfyQAAACh4Xlvn5dM3n5pInHruZ7XdXZPzxV3V1R86raGhh1LRET19l79Oc/rU72937jLcWpfbWx89U/nzLnm84XCzi90di49Q0SkUNgZ6++/785o9MjHGhtfPyUWW/CE1l3nTvR5rluzLp2OtE70cN2adaXut1DYcVU6HX3Bdavu7+paccp0v390ut9wuixatCjoEkrC1r5ERLq7u82/A6glWbLvDVLK3kutbZ1HtvYFAEAYRKPH3ldVtapTRKSv7/pWpeKdNTX3vCIiEokc8VSxmD61p2f1cSJDdfPmvf5DEZFE4vq3Bwbuf3BoaNuFIvJ8b+/VJ4joaEPD9nuVisvcuc8/6bpVL0/0ealUdrVvzY1TXr785srKL/7OcebvzWYvuTCf33RXb+83VyQSN7w1XZ8R2iC1devWoEsoCVv7gr9snUe29gUAQBg4ztzd+55F8kpVdO57HhsUKVQVCm/OFyk0ptOR1rEvVWpOq4hIsdjeKBLrUGrs1Xzl6dJWbq629mcvjf65oaH1UddNLMvnHzkzkbjhp9P1GaG9tG/VqlVBl1AStvYFf9k6j2ztCwCAmSISOcIVKXu7ubl40pjHolSqZ9Xw/nkZkaF5Wg+OeVU+NdF7Dd9/5bw40cN1kxt8aegdSovoab1cKbRBasMGn3+2PrG1L/jL1nlka18AAMwUNTU/fEkksqej4wMr9+7dGve8Tieb/cpR3d2XHCcikkjcsl1EipnMws8Vi25k9+6P/C+t+4+f6L1SqZ6V+1bNe/cjlepZOf74YtGNDA29VCbiOSLaGRp6qaxYdCOT3T9qYODR6s7OpWeM7s9klizXuu/kePyC56fxRxXeS/sAAAAA+Mtx6r1E4jtf6uu7/prdu5c8I6LLlCp/o6zsrO+JiESj7x+qqLjo8oGBh27o6DjsSqVqfqFU3b9Px2fv3n3ammLxzctHn2cyJ6yIRI78wbx5b9xxsP2um9zgOA2t8+btXKd1Njo4+NyVmczCPxERT6ny1+PxpWuSyVvemI46R4X2C3mVUiVbhCBItvYlcoi9tSRFWnrCU88MYWtvtvYFAJgUe1eJ8klQX8g7C4X7C3nb2tqCLqEkbO0L/rJ1HtnaFwAAsE9og5Stq3fZ2hf8Zes8srUvAABgHy7t85mtfYlwaZ+fbO3N1r4AAJPCpX1TxKV9vgn3pX0AAAAAEFYEKQAAAAAwFNogtW7duqBLKAlb+4K/bJ1HtvYFAADsE9p7pDDzcI8UAACYAu6RmiLukfJNTVubbAntF/La+pdgW/uCv2ydR7b2BQDAbOO6c2/WuudUEa9SJJqJRhfc3di4/aHhfdX3aZ07UUQKw0fHOpqbB5fu770GBh5IZrNfulHrntNFot2x2MJbGxo2b/KlkQMI7aV9AAAAAGamysqL1jU0vPTR5ubioqqqr3y5UNhxZXf3F44Z3R+NHn99c7O3cPix/xAlIpLNrlkr4uytr3/mtPLyj399aKj1up6er76/9F0cGEEKAAAAsJzn9amOjqNWp9OxZ9Pp6Au7dh13cTqtft3ff29tKT4vmbxjZyy2YGj4maNFlC4UfnuE6fsMDj5foXXXuZWVl90ej5/VX1f38FbHqX06n39sxXTXbCq0l/YtW7Ys6BJKwq++uERqYrb8XPj9AAAAJjKZ4y/3vF2nVVe3XByLndDb3f2ZH4tEspWVX+g+2Gtdt2ad1n2LJ9qnVPXWVCq7eqJ97e2ptZ7X8UkRXS5S/j+JxI2/GN1XKOy4Kp2OXqVU/I14/Jzv1dU9tmWi9xgY+MmRIspLJr/75r7PrH/V8zpOOVjdpRbaIPUv//IvQZdQErb2BX/ZOo9s7QsAgCD1999bWyz+4dLq6rXLq6uvTYuIOE7Ts56XOTmXW1/f2/u1H4qogogqJhLf+XpV1eWZsa/fX1A6mKYm9zrP6/xWNrt64dDQi6dEIu/dKyJSXr785srKL/7OcebvzWYvuTCf33RXb+83VyQSN7w1/j08L1sp4vSN3aZUxR6RYtWh1DSdQntp3/Lly4MuoSRs7Qv+snUe2doXAABByuXWnapU/HfV1X/XNrpN63yNUonfVFT8efe8ee0XpVJ7PhuNfuDRPXu+96np/GzHqffq6h7eqnWuqavr/ItERGprf/ZSPH5OLhZbMNTQ0PqoUlXb8vlHzpz49TX9It6cd2/NV4lEctNZ56EI7RmpTZsCX4ijJGztC/6ydR7Z2hcAAEHSuq9WpKxr9Hmx6EY8L3N2LHbsnY5T7+07Ll8Vicx/bfzrXTe5Qeu+kyZ6b6WqW1OpnpWTqCLieT37uUdKaRE94fL3FRWXvtnf/+NIb+8170kkbvq9iIjndR7tOIk/qtNvoT0jBQAAAGDqIpH3vK5138K+vpbDBgZ+nshkTmwRGTo8Gj3mNRGRbPavj3bdyocKhZ2fjcfP+5/xr0+lelbuW2Hv3Y+JQlQud09dJnPqBYODz1Z6XqfT2bn0DM/LLItGj3phYODR6s7OpWcMDb1UViy6kUxmyXKt+06Oxy94fqLa4/EzBpSqeyqX23DF4ODzFV1df7HI87rOKS9f8di0/6AMhfYLeW1ZFGA8v/oK4uc3E76Q15Z5ZUsf49naFwBgUvhC3ik60Bfytrc3Xe95mWXD3+n0vp8WCq/87dy5LywsK1uSHz0mkznt/GLx9Q83NbWvnUod/f331vb0XH6H1gNHi4gjEmuLRo++r7Fx+4P9/ffWZrNrNogM/omIeEqVv15Wdubf19f/23+Nvt51kxscp6F13ryd60TGf49UJBuLLbol4O+Rqmlrky2hDVKYGoKU3UEKAAALEaSm6EBBaqzOzvNPHxx8Zm1z8+C5Q0OvxEaXKR8+U7TjjKamt28qebEzW01bm2wJ7aV969evD7qEkrC1L/jL1nlka18AAIRJsfj79ylV+VsRkb6+645x3Tn/5LrV/7h37y8/X119zT1B1zdThPaMlK1nDri0bxzOSB0SW/oYz9a+AACTwhmpKZrsGan29qbrlKronDfvje/7UJaNatraZEtoV+0DAAAAMP2meg8UhoX20j4AAAAACKvQBqmNGzcGXUJJ2NoX/GXrPLK1LwAAYJ/QBqnFixcHXUJJ2NoX/GXrPLK1LwAAYJ/QBqn58+cHXUJJ2NoX/GXrPLK1LwAAYJ/QBikAAAD4o1gsyqZNm+Rb3/qWbNq0SYrFYtAlYZZx3er70mnn5XTaeXH4EX9i7P6BgQeSrlv7w3Ta2Z5Ol/1nJrNkmcn+Qz32QFi1DwAAYBYrFoty3nnnyebNmyWXy0lVVZUsWbJEnnzySYlEIkGXh1kkGj3++sbG7Q9NtC+bXbNWxNlbX//MabncDxbk84+s7+n56qvJ5B07J7Pf5L0mK7RnpFauXBl0CSVha1/wl63zyNa+ACDMHn/8cdm8ebPs2bNHtNayZ88e2bx5szz++ONBl4ZplE6XPdPR8cHL0umKjem0s729vfHbudz6etdNbkinI9tct/regYGfJ0REOjqOWpVOl/1HOh3Zlk6X/9vu3R85Z/R9enu/eXg6HdnS3X3ph0RE9uy5vTGdjr7Q1bXilFLVPjj4fIXWXedWVl52ezx+Vn9d3cNbHaf26Xz+sRWT2W/yXiZCG6TWr18fdAkl4WdfSilfH2E0E2o8FPx+AACmy4svvii5XO5d23K5nGzfvj2gilAqnvf2ecnk7ZcmEree63ldZ/f0XHF3RcWnbmto2LFERFRv79WfExFxnMY/JBI3X9zUlF0ciy36wd69z9+Sy/2gQUQkkbjhrWj0uFsGBu6/de/ezeV9ff/nO47T9Ehd3WNbxn6W69asS6cjrRM9XLdm3UT1FQo7rkqnoy+4btX9Y4PZwMBPjhRRXjL53TdHtylV/6rn9R41mf1jmRx7MKENUrau3uVnX1prXx9hNBNqPBT8fgAApsvChQulqqrqXduqqqrkxBNPDKgilEo0eux9VVWrOufM+doupapbHSf5q5qae16JxRYMRSJHPOV53QtERBoa/t8Tc+ZcsctxqnVDw3/9m0j89wMDDx8/+j6NjdsfVKr8zd27z3xI66GG+vqnbhv/WalUdnVzc/GkiR6pVHb1+OPLy5ffXF//xMcaGl7+02h0wQP5/Ka7enu/ebiIiOdlK0WcvrHHK1WxR6RYNZn9Y5kcezChDVLbtm0LuoSSsLUv+MvWeWRrXwAQZueff74sWbJE5syZI0opmTNnjixZskTOP//8oEvDNHOcubv3PYvklaro3Pc8NihSqBIRyWRO+kQ6XfHY6BkkkfxRntdTO/a9YrFFD4oMfiAa/eBPY7EFQ1Otrbb2Zy/F4+fkYrEFQw0NrY8qVbUtn3/kzOG6a/pFvDnvfkW+SiSSm8z+sUyOPZjQBikAAACUXiQSkSeffFLuv/9+uf766+X+++9noYlZrK/v281DQ9tuKC+/8Pqmpl2nNDcXTxIpf01Ev3OPxODgs5V79z5/reM0PlQo/PqrAwMPJMe/z/C9V6Mr8L374brJDQevROnRz6youPRNER3p7b3mPaN7Pa/zaMdJvDaZ/WOZHHswoV21L5VKBV1CSdjaF/xl6zyytS8ACLtIJCLLli2TZcsOaRVoWMTzOipEREciR3aJiGQyiz4pkn/X/UPd3Z++Vqk5O5qaOr7Z3j7vW9nsl6+rqPjMlWOPSaV6Jr2C1MDAo9X9/XedkEh8d4vjNBS7uj5xgdZ9J8fjF9woIhKPnzGgVN1TudyGK+LxZdfmcrcv8Lyuc6qq1nxmMvvHMjn2YEJ7RiqdTgddQknY2hf8Zes8srUvAABmimTy+7+LRA7/cS73vQfa2+f9V7G46wNKVb1z7f3u3Wd9zPO6P1JT86O1IiK1tQ98R+vcMZnMkuWH+plaZ6ODg89dmcksfKGj47DNhcKOS+LxpWuSyVveGD2mpubOFhEv3tn5kV/m84/eFoudtHbscuUH2u+6yQ0dHe9fPdn3mix1kBvwA7s7v6WlRVpaWoL6+JLxqy+llO+LKxzSZ7YkRVp6fKsniJ9LKfD7AQCwkD3L6wZk/nw5RUSyQdcxC9S0tcmW0AYpW/7CO55ffRGk7A5StvQxnq19AQAmhSA1RQQp39S0tcmW0F7aBwAAAABhRZACAAAAAEOhDVKtra1Bl1AStvYFf9k6j2ztCwAA2Ce0QQoAAAAAworFJnzGYhPjsNjEIbGlj/Fs7QsAMCksNjFFLDbhm5q2NtkS2i/kBQAAAGCkT0Rqgi5iFugTESFIAQAAABZoa5NXgq5hNgntPVJr164NuoSSsLUv+MvWeWRrXwAAwD6hvUcKU8M9UnbfIwUAgIW4RwozSmjPSDU3NwddQknY2hf8Zes8srUvAABgn9AGKdd1gy6hJPzqi7MuE7Pl58LvBwAAQLBCG6QAAAAAIKxCG6QWLVoUdAklYWtf8Jet88jWvgAAgH1YbALTZiYsNgEAAEKLxSYwo4T2jNSqVauCLqEkbO0L/rJ1HtnaFwAAsE9oz0jZejbB1r5EOCPlJ1t7s7UvAMCkcEYKM0poz0gBAAAAQFgRpAAAAADAUGiDVFtbW9AllIStfcFfts4jW/sCAAD2CW2Q2rp1a9AllIStfcFfts4jW/sCAAD2YbEJn9nalwiLTfjJ1t5s7QsAMCksNoEZJbRnpAAAAAAgrAhSAAAAAGAotEFq3bp1QZdQErb2BX/ZOo9s7QsAANgntPdIYebhHikAADAF3COFGSW0Z6SUsvN3yda+4C9b55GtfQEAAPuENkhhZlJKGT0O5TWTfdTW1gb80wAAAICtQntpn62XZdnaF/xl6zyytS8AwKRwWQJmlNCekVq2bFnQJZSErX3BX7bOI1v7AgAA9gntGSkAAADMKpyRwowS2jNSy5cvD7qEkrC1L/jL1nlka18AAMA+oT0jZeu9Erb2BX/ZOo9s7QsAMCmckcKMEtozUgAAAAAQVgQpAAAAADAU2kv7AAAAMKtwaR9mlNCekVq/fn3QJZSErX3BX7bOI1v7AgAA9gntGSlbbzq3tS/4y9Z5ZGtfAIBJ4YwUZpTQnpECAAAAgLAiSAEAAACAodAGqY0bNwZdQknY2hf8Zes8srUvAABgn9AGqcWLFwddQknY2hf8Zes8srUvAABgHxab8JmtfcFfts4jW/sCAEwKi01gRgntGSkAAAAACCuCFAAAAAAYCm2QWrlyZdAllIStfcFfts4jW/sCAAD2Ce09UgAAAJhVuEcKM0poz0jZunqXrX3BX7bOI1v7AgAA9gntGSlbV++ytS/4y9Z5ZGtfAIBJ4YwUZpTQnpECAAAAgLAKbZBKpVJBl1AStvYFf9k6j2ztCwAA2Ce0l/YBAABgVuHSPswooT0j1dLSEnQJJWFrX/CXrfPI1r4AAIB9QntGytabzm3tC/6ydR7Z2hcAYFI4I4UZJbRnpAAAAAAgrAhSAAAAAGAotEGqtbU16BJKwta+4C9b55GtfQEAAPuENkgBAAAAQFix2ITPbO0L/rJ1HtnaFwBgUlhsAjMKZ6QAAAAAwBBBCgAAAAAMRQ7yBZgH3FlqZ511VpAfXzK29gV/2TqPbO0LAHBQ1wVdAGAitPdIAQAAYFbhHinMKFzaBwAAAACGCFIAAAAAYIggBQAAAACGCFIAAAAAYIggBQAAAACGCFIAAAAAYIggBQAAAACGCFIAAAAAYIggBQAAAACGCFIAAAAAYIggBQAAAACGCFIAAAAAYIggBQAAAACGCFIAAAAAYIggBQAAAACGCFIAAAAAYIggBQAAAACGCFIAAAAAYIggBQAAAACGCFIAAAAAYIggBQAAAACGCFIAAAAAYIggBQAAAACGCFIAAAAAYIggBQAAAACGCFIAAAAAYIggBQAAAACGCFIAAAAAYIggBQAAAACGCFIAAAAAYIggBQAAAACGCFIAAAAAYIggBQAAAACGCFIAAAAAYIggBQAAAACGCFIAAAAAYIggBQAAAACGCFIAAAAAYIggBQAAAACGCFIAAAAAYIggBQAAAACGCFIAAAAAYIggBQAAAACGCFIAAAAAYIggBQAAAACGCFIAAAAAYIggBQAAAACGCFIAAAAAYIggBQAAAACGCFIAAAAAYIggBQAAAACGCFIAAAAAYPssigwAAAaISURBVIggBQAAAACGCFIAAAAAYIggBQAAAACGCFIAAAAAYIggBQAAAACGCFIAAAAAYIggBQAAAACGCFIAAAAAYIggBQAAAACGCFIAAAAAYIggBQAAAACGCFIAAAAAYIggBQAAAACGCFIAAAAAYIggBQAAAACGCFIAAAAAYIggBQAAAACGCFIAAAAAYIggBQAAAACGCFIAAAAAYIggBQAAAACGCFIAAAAAYIggBQAAAACGCFIAAAAAYIggBQAAAACGCFIAAAAAYIggBQAAAACGCFIAAAAAYIggBQAAAACGCFIAAAAAYIggBQAAAACGCFIAAAAAYIggBQAAAACGCFIAAAAAYIggBQAAAACGCFIAAAAAYIggBQAAAACGCFIAAAAAYIggBQAAAACGCFIAAAAAYIggBQAAAACGCFIAAAAAYIggBQAAAACGCFIAAAAAYIggBQAAAACGCFIAAAAAYIggBQAAAACGCFIAAAAAYEhprfe/U6knRGSuf+UcsrkisjvoIiAijEXYMB7hwViEC+MRLoxHeAQ5Fru11ksD+mzA2AGD1EyhlGrVWp8UdB1gLMKG8QgPxiJcGI9wYTzCg7EAJo9L+wAAAADAEEEKAAAAAAzZEqTWB10A3sFYhAvjER6MRbgwHuHCeIQHYwFMkhX3SAEAAACAn2w5IwUAAAAAvplRQUop9RdKqV8rpTyl1Enj9n1DKbVTKfUbpdR5Y7YvVkq9PLLv+0op5X/ls4NSaunIz3+nUuqaoOuxnVLqx0qpXUqpHWO21SmlnlJKvTbyz9ox+yb8HcH0UEodrpT6T6XUKyP/nrpiZDtj4jOlVLlSaotS6lcjY3HdyHbGIkBKqYhS6kWl1KaR54xHQJRSb4783Wi7Uqp1ZBvjARiaUUFKRHaIyCdF5LmxG5VSHxKRvxSRY0RkqYjcqZSKjOz+kYisEpGjRh58P0EJjPy8fygi54vIh0TkopFxQencK388n68Rkae11keJyNMjzw/2O4LpURCRq7TWC0TkwyLylZGfO2Piv0EROVtrfYKInCgiS5VSHxbGImhXiMgrY54zHsH6qNb6xDFLnTMegKEZFaS01q9orX8zwa4VIvLPWutBrfUbIrJTRE5RSqVEJKG1/qUevhnsH0XkEz6WPJucIiI7tdava633isg/y/C4oES01s+JSNe4zStE5B9G/vwPsm++T/g74kuhs4TW2tVabxv5c58M/4VxvjAmvtPD9ow8jY08tDAWgVFKHSYiF4rI3WM2Mx7hwngAhmZUkDqA+SLy1pjnb49smz/y5/HbMf32Nwbw1zyttSsy/Bd7EWkc2c74+EgpdaSILBSRzcKYBGLkMrLtIrJLRJ7SWjMWwfp7EfnfIuKN2cZ4BEeLyL8rpbYqpVaNbGM8AEPRoAsYTyn1HyLSNMGua7XWj+3vZRNs0wfYjunHzzrcGB+fKKXmiMjPReRKrXXvAW7LZExKSGtdFJETlVI1IvKIUurYAxzOWJSQUmqZiOzSWm9VSp01mZdMsI3xmF6na63TSqlGEXlKKfXqAY5lPID9CF2Q0lqfcwgve1tEDh/z/DARSY9sP2yC7Zh++xsD+KtDKZXSWrsjl7buGtnO+PhAKRWT4RD1T1rr/zuymTEJkNY6q5R6Vobv7WAsgnG6iHxcKXWBiJSLSEIp9VNhPAKjtU6P/HOXUuoRGb5Uj/EADNlyad9GEflLpVRcKfVeGV5UYsvIqek+pdSHR1br+5yI7O+sFqbmv0XkKKXUe5VSZTJ8Y+rGgGuajTaKyOdH/vx52TffJ/wdCaA+a438O+YeEXlFa33bmF2Mic+UUg0jZ6JEKVUhIueIyKvCWARCa/0NrfVhWusjZfi/Dc9orT8rjEcglFJVSqnq0T+LyLkyvJgX4wEYCt0ZqQNRSv2ZiNwhIg0i8q9Kqe1a6/O01r9WSj0oIv8jwytnfWXksg4RkS/L8OpmFSLy+MgD00xrXVBKXS4iT4pIRER+rLX+dcBlWU0pdb+InCUic5VSb4vIWhG5SUQeVEpdJiJ/EJG/EBE5yO8IpsfpInKJiLw8cm+OiMjfCmMShJSI/MPIymKOiDyotd6klPqlMBZhwu9GMObJ8OWuIsN/D/yZ1voJpdR/C+MBGFHDi9kBAAAAACbLlkv7AAAAAMA3BCkAAAAAMESQAgAAAABDBCkAAAAAMESQAgAAAABDBCkAAAAAMESQAgAAAABDBCkAAAAAMPT/AXIJ2+fqiPxxAAAAAElFTkSuQmCC\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "# nbi:hide_in\n", "import matplotlib.pyplot as plt\n", "import numpy as np\n", "\n", "# set up the figure\n", "fig, ax = plt.subplots(1,1, figsize=(12,5))\n", "\n", "# data\n", "data = np.array([-30, -5, 9, 10, 13, 20, 40, 500]).astype(float)\n", "# data = np.array([-30,-1, -5, -0.5, 0.5, 0.6, 0, 2, 3, 4.6, 4, 7, 18, 35]).astype(float)\n", "\n", "\n", "def percentile(data, p):\n", " \"\"\"\n", " Compute the percentiles the way we defined in class \n", " data : array of size N \n", " p : percentile\n", " \"\"\"\n", " data = np.sort(data, axis=0)\n", " rank = int(p * (data.shape[0] + 1) - 1) # the rank\n", " assert rank > 0, \"the rank does not exist\" \n", " alpha = p * (data.shape[0] + 1) - 1 - rank # the fractional part\n", " return data[rank] + alpha * (data[rank + 1] - data [rank])\n", "\n", "def box_plot(ax, data, width=0.4, showout = True, position = np.array([0.4])):\n", " \"\"\"\n", " ax : matplotlib ax\n", " data : the data \n", " width : box width\n", " showout : show the outliers \n", " position: the y axis of the box plot\n", " \"\"\"\n", " # compute the five number summary \n", " minim = np.min(data)\n", " maxim = np.max(data)\n", " q1 = percentile(data, 0.25)\n", " q2 = np.median(data)\n", " q3 = percentile(data, 0.75)\n", "\n", " # interquartile range\n", " iqr = q3 - q1\n", "\n", " # inner fences\n", " left_innerfence = q1 - 1.5 * iqr\n", " right_innerfence = q3 + 1.5 * iqr\n", "\n", " # compute outliers \n", " outliers = data[np.logical_or(data = right_innerfence)]\n", " \n", " # whiskers\n", " if showout==True:\n", " low_whisker = np.min(data[data >= left_innerfence])\n", " high_whisker = np.max(data[data <= right_innerfence])\n", " else:\n", " low_whisker = np.min(data)\n", " high_whisker = np.max(data)\n", "\n", "\n", "\n", " stats = [{'iqr': iqr,\n", " 'whishi': high_whisker,\n", " 'whislo': low_whisker,\n", " 'fliers': outliers,\n", " 'q1': q1,\n", " 'med': q2,\n", " 'q3': q3}]\n", "\n", " # add the box plot\n", " flierprops = dict(markerfacecolor='black', markersize=5)\n", " ax.bxp(stats, vert = False, widths=width, positions = position, \n", " flierprops=flierprops, showfliers=showout)\n", "\n", " # add Tukey's fences\n", " if showout==True:\n", " ax.vlines(q1-1.5*iqr, position-0.2,position+0.2, linestyle=\"dashed\", linewidth=1)\n", " ax.vlines(q3+1.5*iqr, position-0.2,position+0.2, linestyle=\"dashed\", linewidth=1)\n", "\n", " ax.vlines(q1-3*iqr, position-0.2,position+0.2, linestyle=\"dashed\", linewidth=1)\n", " ax.vlines(q3+3*iqr, position-0.2,position+0.2, linestyle=\"dashed\", linewidth=1)\n", "\n", " # \n", " ax.spines['top'].set_visible(False)\n", " ax.spines['right'].set_visible(False)\n", " ax.spines['left'].set_visible(False)\n", " ax.set_ylim(-0.1,position+0.3)\n", " ax.set_yticks([])\n", " plt.figtext(1,0.8,\n", " r\"$\\min={:.4}$\".format(minim)+\"\\n\"+\n", " r\"$q_1={:.4}$\".format(q1)+\"\\n\"+\n", " r\"med$={:.4}$\".format(q2)+\"\\n\"+\n", " r\"$q_3={:.4}$\".format(q3)+\"\\n\"+\n", " r\"max$={:.4}$\".format(maxim),\n", " ha=\"left\", va=\"top\",\n", " backgroundcolor=(0.1, 0.1, 1, 0.15),\n", " fontsize=\"large\")\n", " \n", "def disp_data(ax, data):\n", " ax.scatter(data, np.zeros(data.shape), zorder=2, s=10)\n", " ax.set_yticks([])\n", "# ax.set_xticks([])\n", " mean = np.mean(data)\n", " ax.scatter(mean, 0, zorder=2, s=20, color=\"red\")\n", " ax.set_ylim(-0.01,0.1)\n", " ax.axhline(y=0, color='k', zorder=1, linewidth=0.5)\n", "\n", " \n", " ax.spines['top'].set_visible(False)\n", " ax.spines['right'].set_visible(False)\n", " ax.spines['left'].set_visible(False)\n", " ax.spines['bottom'].set_visible(False)\n", "\n", " ax.set_ylim(-0.1,1.5)\n", " \n", "box_plot(ax, data, width=0.2, showout=True)\n", "\n", "plt.show();" ] } ], "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 }