{ "cells": [ { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "# CS 236756 - Technion - Intro to Machine Learning\n", "---\n", "\n", "#### Tal Daniel\n", "\n", "## Tutorial 01 - Probability Refresher and Maximum Likelihood Estimator (MLE)" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "subslide" } }, "source": [ "### Agenda\n", "---\n", "* [Probability Basics](#-Probability-Basics)\n", "* [Bayes Rule](#-Bayes-Rule)\n", "* [Expectation & Variance](#-Mean-&-Variance)\n", "* [Correlation](#-Correlation)\n", "* [Maximum Likelihood Estimator](#-Maximum-Likelihood-Estimation-(MLE))\n", "* [Recommended Videos](#-Recommended-Videos)\n", "* [Credits](#-Credits)" ] }, { "cell_type": "code", "execution_count": 1, "metadata": { "slideshow": { "slide_type": "skip" } }, "outputs": [], "source": [ "# imports for the tutorial\n", "import numpy as np\n", "import pandas as pd\n", "import matplotlib.pyplot as plt\n", "%matplotlib notebook" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "## Probability Basics\n", "---\n", "We define the following:\n", "* **Experiment** - an experiment or trial is any procedure that can be infinitely repeated and has a well-defined set of possible outcomes, known as the sample space.\n", " * Example: toss a coin twice\n", "* **Sample Space ($\\Omega$)** - possible outcomes of an experiment\n", " * Example (coin toss): {HH, HT, TH, TT} (H = Heads, T = Tails)\n", "* **Event** - a subset of possible outcomes\n", " * Example (coin toss): A = {HH} , B= {HT, TH}" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "subslide" } }, "source": [ "* **Probability (of an event)** - a number assigned to an event\n", " * Example (coin toss): $Pr(A) = \\frac{1}{4}$\n", "* **Axioms**:\n", " 1. $0 \\leq Pr(A) \\leq 1$\n", " 2. $Pr(\\Omega) = 1 $\n", " 3. $Pr(A \\cup B) = Pr(A) + Pr(B) - P(A \\cap B)$ (if $A, B$ are independent $P(A \\cap B) = 0$)" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "subslide" } }, "source": [ "\n", "(image from tistats.com)" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "subslide" } }, "source": [ "#### Summary\n", "\n", "\n", "| Term | Usually donated by | Definition | Example |\n", "| --- | --- | --- | --- |\n", "| **Experiment** | |
any procedure that can be infinitely repeated and has a well-defined set of possible outcomes
|
toss a coin twice
|\n", "| **Sample** | $\\omega$ |
A single outcome of an experiment
|
A single outcome. for example: $H$
|\n", "| **Sample Space** | $\\Omega$ |
The set of all possible outcomes
|
The set of all possible outcomes, for example: {HH, HT, TH, TT}
|\n", "| **Event** | $A$ |
a subset of possible outcomes
|
A = {HH} , B= {HT, TH}
An empty set $\\emptyset$
The entire set (any outcome): $\\Omega$
|\n", "| **Event Space** | $\\mathcal{F}$ |
The space of all possible events
|
$\\{HH\\}, \\{HT\\}, \\{TH\\}, \\{TT\\}, \\{\\emptyset\\}, \\{\\Omega\\}$
|\n", "| **Probability** | P, Pr |
A function $P:\\mathcal{F}\\rightarrow\\left[0,1\\right]$ which assigns a probability to each event
|
$P(HT) = \\frac{1}{4}$
$P(\\emptyset)=0$
$P(\\Omega)=1$
|\n" ] }, { "cell_type": "code", "execution_count": 2, "metadata": { "slideshow": { "slide_type": "subslide" } }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "first toss result: H\n", "second toss result: T\n" ] } ], "source": [ "possible_outcomes = ['H', 'T']\n", "probabilities = [0.5, 0.5]\n", "# toss a coin twice\n", "first_toss = np.random.choice(possible_outcomes, p=probabilities)\n", "print(\"first toss result: \", first_toss)\n", "second_toss = np.random.choice(possible_outcomes, p=probabilities)\n", "print(\"second toss result: \", second_toss)" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "### Exercise 1 - Dice Probability\n", "---\n", "Find the proability of getting an even number **or** a number that is a multiple of 3." ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "subslide" } }, "source": [ "### Solution 1\n", "---\n", "* $S = \\{1, 2, 3, 4, 5, 6\\}$\n", "* $ P(A=even) = \\frac{3}{6} = 0.5$\n", "* $ P(B=\\textit{multiple of 3}) = \\frac{2}{6} = \\frac{1}{3}$\n", "* $ P(A \\cap B) = \\frac{1}{6}$, since $ A \\cap B = \\{6\\}$\n", "* $\\rightarrow P(A \\cup B) = \\frac{3}{6} + \\frac{2}{6} - \\frac{1}{6} = \\frac{4}{6} = \\frac{2}{3} $\n", "\n", "Why is it different than just giving the probability of getting a 6? Because we ask about 2 different events (even number and a multiple of 3), and we ask for the union of 2 events." ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "## Joint Probability\n", "---\n", "* **Joint Probability - Pr(A,B)** - the probability the both event A and event B happen ($Pr(A,B) \\geq 0$)\n", "* **Marginal Distributions** -\n", " * $\\sum_i Pr(A_i, B_j) = Pr(B_j)$\n", " * $\\sum_j Pr(A_i, B_j) = Pr(A_i)$\n", "* **Law of Total Probability** - suppose the events $B_1, B_2, ..., B_k$ are mutually exclusive (intersection of all events is zero) and form a partition of the sample space (i.e. one of them must occur), then for any event $Pr(A)$: $$Pr(A) = \\sum_{j=1}^k Pr(A|B_j)P(B_j)$$" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "subslide" } }, "source": [ "* **Conditioning** - if A and B are events and $Pr(B) > 0$, the **conditional probability of A given B** is: $Pr(A|B)$\n", " * $Pr(A|B) = \\frac{P(A,B)}{P(B)}$\n", " * $\\rightarrow Pr(A,B) = Pr(A|B)Pr(B)$\n", " * **the chain rule** - in the general case: $$Pr(\\bigcap_i A_i) = \\prod_{i=n}^1 Pr(A_i|A_{i-1}, ..., A_1)$$\n", " * $Pr(A,B,C) = P(A|B,C)P(B,C) = P(A|B,C)P(B|C)P(C)$\n", " * Is that the only option? No! $Pr(A,B,C) = P(C|A,B)P(A,B) = P(C|A,B)P(B|A)P(A)$" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "subslide" } }, "source": [ "* **Independence** - \n", " * Two events A and B are independent iff\n", " * $Pr(A|B) = Pr(A)$\n", " * $Pr(A,B) = Pr(A)Pr(B)$\n", " * For a set of events $\\{A_i\\}$, independence:\n", " * $Pr(\\bigcap_i A_i) = \\prod_i Pr(A_i)$" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "subslide" } }, "source": [ "* **Conditional Independence** - Events A and B are conditionally independent given C:\n", " * $Pr(A,B|C) = Pr(A|C)Pr(B|C)$\n", " * $Pr(A|B,C) = Pr(A|C)$" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "## Continuous Random Variables\n", "---\n", "Assume $X$ is a continuous random variable. We define the following:\n", "* **Cumulative Distribution Function (CDF)** - $F(x) = P(X \\leq x)$\n", " * The CDF is monotonically non-decreasing\n", " * $ P(a < X \\leq b) = F(b) - F(a)$\n", " * $\\lim_{a\\to\\infty} F_x(a) = 1$\n", " * $\\lim_{a\\to 0} F_x(a) = 0$" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "subslide" } }, "source": [ "* **Probability Density Function (PDF)** - $f(x) = \\frac{d}{dx}F(x)$\n", " * $P(a < X \\leq b) = \\int_a^b f(x)dx $\n", "* All we have seen for **discrete** variables hold for **continuous** by replacing the sum with an integral\n", "* Note that unlike the discrete case the PDF can be larger than one, i.e. $p(x) > 0$ but is not upper bounded. In order for the integral of the PDF to be smaller than 1, the PDF can be larger than one, but not for long intervals (recall that the CDF is just the area under the PDF)." ] }, { "cell_type": "code", "execution_count": 3, "metadata": { "slideshow": { "slide_type": "skip" } }, "outputs": [], "source": [ "def plot_normal_pdf_cdf(mu=0, sigma=1):\n", " x = np.linspace(-10, 10, 1000)\n", " x_pdf = (1 / np.sqrt(2 * np.pi * sigma ** 2)) * np.exp(- (x - mu) ** 2 / (2 * sigma ** 2))\n", " x_cdf = np.cumsum(x_pdf) / (len(x) / (np.max(x) - np.min(x))) # normalization\n", " fig = plt.figure(figsize=(8,5))\n", " ax = fig.add_subplot(1,1,1)\n", " ax.plot(x, x_pdf, label='PDF')\n", " ax.plot(x, x_cdf, label='CDF')\n", " ax.grid()\n", " ax.legend()\n", " ax.set_xlabel('x')\n", " ax.set_title('PDF and CDF of Normal Distribution')" ] }, { "cell_type": "code", "execution_count": 4, "metadata": { "slideshow": { "slide_type": "subslide" } }, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAeMAAAFNCAYAAAAzYQemAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8QZhcZAAAgAElEQVR4nO3dd5wb5Z3H8c9vu+1dd3td1h0bsDHuNs3EBEILYODihIRQciE+7iBXklwScncJHLlLuQvhCCQESIAAwUBoBkwNNt1t3Y0L7ltc123X3q7n/hitLWu11ZJG0n7fr9da0syj0e/RrPereWY0Y845RERExD9pfhcgIiLS0SmMRUREfKYwFhER8ZnCWERExGcKYxEREZ8pjEVERHymMJYOz8zuNLMn/a6jNcws38zeN7NyM/uV3/W0lZnNMLPik1zG62Z2U5TqmW5mG0IebzOzi6Kx7ODy1prZjGgtT1KXwliiJviHrNLMKsxst5k9ama5wXkLzKwqGCKHzazQzH5oZtkhz7/TzGqDz2/4+b5/PTpWV1czu9fMdgRr2hR83Ds4v6Hf5WZ20Mw+NrNbzSwtZBmPmVlNWN++0o5yZgP7gK7Oue9GqPUxM3NmNjVk2ilmlhQnFAjWfiT4/pSZ2V/D3yfn3GXOucdbuaxTmmvjnPvAOXfqydYdfL3HzOynYcsf45xbEI3lS2pTGEu0XemcywUmAlOAfw+Zd7tzLg/oD3wXuA6YZ2YW0uYZ51xuyM8v41Z5BGaWBfwVGANcCnQFzgHKgKkhTa8M9m0I8HPgB8Afwhb3y7C+PdOOkoYAn7rmz9azH/hpM/NbzcwyorGcNhoX/B06FXgMuN/MfhLtF/GpbyIRKYwlJpxzJcDrwBkR5h0Jbi1cBZwNfLGtyw9uVW8Obo1+ambXhMy72cw+NLP/NbMDZrbVzC4LmT/MzN4LPvdtoHczL3UjMBi4xjn3qXMu4Jzb45y72zk3L0LfDjnn5gJfAW4ys0b9b0XfzjGzJWZ2KHh7TnD6Y8BNwPeDW45NDac+DpxpZp9rYvkDzGyume0PbuV/K2TenWb2FzN70swOAzcHpz0XnFZuZqvNbJSZ3WFme8ysyMwuDlnGN8xsXbDtFjP7u7a+BwDOuX3OuSeAvwfuMLNeweUvMLNbgvdPCa7LQ2a2z8yeCU5/P7iYlQ2jEBYcIjezH5jZLuBRizxsPiX4O3UgOLqTE1zmzWb2Ydh76YI1zAau5/i6eSU4/9iwt5llB0dUSoM/91pwZCiktu8G39OdZvaN9rxvkpwUxhITZjYIuBxY3lQb59wOYCkwvR0vsTn4vG7AXcCTZtY/ZP40YANe0P4S+EPIFvifgcLgvLvxAq4pFwFvOOcq2lKcc24xUEwb+2ZmPYHXgPuAXsA9wGtm1ss5dzPwFMe3sN9pYjFHgf8G/quJ+U8HaxsAfAn4bzO7MGT+TOAvQPfg6wFcCTwB9MBbp2/i/f0YCPwn8PuQ5+8BrsAbRfgG8Gszm9ia/jfhZSCDE0ciGtwNvBWsqwD4DYBz7vzg/HFhoxD9gJ54Iwyzm3i964FLgBHAKE4c3YnIOfcQJ66bKyM0+zfgLGA8MC7Yn9Bl98P7fR4IfBN4wMx6tPTakhoUxhJtL5nZQeBD4D28UGhOKd4fxwZfNm+/a8PPgEhPcs4955wrDW6pPgN8xol/rLc75x52ztXjbSn2B/LNbDDe8Pl/OOeqnXPvA680U18vYGcLfWht374X0q99TTzni8BnzrknnHN1zrmngfV4YdgWvwcGh44IwLEPSecBP3DOVTnnVgCPADeENPvEOfdS8L2tDE77wDn3pnOuDngO6AP83DlXC8wBhppZdwDn3GvOuc3O8x5eWLbnAxfB5dXi7SfvGWF2LV6wDgj258MIbUIFgJ8E131lE23ud84VOef2432g+Wp7aw9zPfCfwZGVvXgfIkPf99rg/NrgqEsF3lC9dAAKY4m2q51z3Z1zQ5xz/9DMH7wGA/H2cTZ4Nvj8hp/SSE8ysxvNbEVDuOENh4cON+9quOOcOxq8m4u3NXjAOXckpO32Zuorwwvy9gjv2/+G9KupofEBEerZHlxWqznnqvG2Gu8GQvfJDwD2O+fKm1l+UYRF7g65XwnsC37QaXgM3vuLmV1mZguDw+AH8UZImtsV0Cwzy8QL//0RZn8fr3+LzTty+W9bWNxe51xVC21C+78d7z2LhvB1G77ssuCHnQZHCb6nkvoUxuKb4FbaJOCDNj5vCPAwcDvQyznXHVjDiaHTlJ1ADzPrEjJtcDPt3wEuCWvfmhqn4AVcS1tq4UrxtvRCDQZK2rgcgEfxhj2vCZlWCvQ0s7xmlt/uI6+D+0CfB/4XyA+um3m0bt00ZSZQBywOn+Gc2+Wc+5ZzbgDwd8BvrfkjqFvTt0Eh9wfjvWcAR4DODTPMrF8blx2+bkOXLR2cwljizsw6Bw8uehnvD2yjA6Fa0AXvD9/e4PK+QYQDxSJxzm3H2099l5llmdl5ND8E/ATeltLzZnaamaWZWS8z+5GZXR7e2LyvQV2BN3T7pHNudZt65r0Xo8zsa2aWYd7XekYDr7ZxOQS3su7EO7K7YVoR8DHwMzPLMbMz8fZPPhVxIW2XBWTjrZu64DD5xc0/JTIz62lm1wMPAL9wzpVFaDPLzAqCDw/g/V40bLHvBoa346VvM7OC4P77HwEN+5tXAmPMbHzwoK47w57X0us9Dfy7mfUx72txPwaS4vvtEnsKY4mn+82sHO+P1r14W1CXOucCbVmIc+5T4FfAJ8FljQU+asMivoZ3gNd+4CfAn5p5rWq8g7jWA28Dh/E+QPQGFoU0fSXYtyK8A3XuwTt4qU2CgXMF3le/yvCGYa9wzjW1j7klT9N4n/dXgaF4W2Uv4u1Dfbudyz9BcPj7H4Fn8cLxa8DcNi5mpZlVAJuAW4B/cc79uIm2U4BFwfZzgX9yzm0NzrsTeDy4K+PLbXj9P+Pt594S/PkpgHNuI97Bau/gHaMQPurxB2B08PVeirDcn+J9EFwFrAaWEaWvoEnys+a/rigiIiKxpi1jERERnymMRUREfKYwFhER8ZnCWERExGcKYxEREZ/5dtWS3r17u6FDh0ZteUeOHKFLlzadlyFhqS+JKVX6kir9APUlUaVKX2LRj8LCwn3OuT7h030L46FDh7J06dKoLW/BggXMmDEjasvzk/qSmFKlL6nSD1BfElWq9CUW/TCziKff1TC1iIiIzxTGIiIiPlMYi4iI+My3fcaR1NbWUlxcTFVVS1c4a6xbt26sW7cuBlVFT05ODgUFBWRmZvpdioiIJJCECuPi4mLy8vIYOnQoZm274lp5eTl5eXktN/SJc46ysjKKi4sZNmyY3+WIiEgCSahh6qqqKnr16tXmIE4GZkavXr3atdUvIiKprcUwNrM/mtkeM1vTxHwzs/vMbJOZrTKziSdTUCoGcYNU7puIiLRfa7aMHwMubWb+ZcDI4M9s4HcnX5Z/0tPTGT9+PGeccQazZs3i6NGjJ0wfM2YM48aN45577iEQ8C7Du2DBArp168b48eMZP348F110kZ9dEBGRJNNiGDvn3se7CHtTZgJ/cp6FQHcz6x+tAuOtU6dOrFixgjVr1pCVlcWDDz54wvS1a9fy9ttvM2/ePO66665jz5s+fTorVqxgxYoVvPPOO36VLyIiSSgaB3ANBIpCHhcHp+2MwrJ9NX36dFatWtVoet++fXnooYeYMmUKd955Z/wLE5HWC9RDdTnUVHi3tZUQqIP6WgjUQn0d1NcE79eCc8EnuuD94OOG+y3Nb6f+pRuhcFu7n59IUqUv/Us3Qt05kJEV89eKRhhH2hEa8TfSzGbjDWWTn5/PggULTpjfrVs3ysvL21VEfX19u58brry8nLq6Ol555RUuuuiiY8sNXX6fPn0IBAJs2bKFo0eP8sEHH3DmmWcCcPXVV/Ov//qvEZddVVXVqN/hKioqWmyTLNSXxJMq/QA4Un6Ixa89SeejRXQ+WkR2dRnZ1fvJqtlPdvUBMuqOkB5IjoMmTwXY6HcV0ZEqfTkV+GDBudRnxP4829EI42JgUMjjAqA0UkPn3EPAQwCTJ0924ef8XLdu3bGvJ931ylo+LT3c6iLq6+tJT09vts3oAV35yZVjmm1TWVnJ9OnTAW/L+LbbbiMry/tUFP7VKeccubm5dO7cmenTp/Pqq6+2WGdOTg4TJkxotk2qnNcV1JdElNT9qK2CbR/C9g9h+8cEigtJc3XH53fuBXn9ocdQyD0LOnWH7DzIyvVus3MhszOkZUJ6RvA2+NNw3xr23hk0HHR57ODLhmmh0y3C/Lb7+JNPOOfss9v13ESTKn35+JNPmP75yyAt9l88ikYYzwVuN7M5wDTgkHMuaYeoG/YNt2TLli2kp6fTt2/fhD/ZiEhSCwRg87uwag5seANqyiEtAwZMoGTgFQyafAn0OQ16j/LCNknVZPeCrgP8LiMqUqUvNdm94hLE0IowNrOngRlAbzMrBn4CZAI45x4E5gGXA5uAo8A3olFYS1uw4eJ50o+9e/dy6623cvvtt+vrSiKxUlcNy/4Ei34PZZ9Bp54w5moYPROGnANZXdi8YAGDxs/wu1KRk9ZiGDvnvtrCfAfcFrWKElRlZSXjx4+ntraWjIwMbrjhBr7zne/4XZZI6nEO1jwPf70LDu6AgZPg2ke8EI7DgTQifkio02EmgoqKiojT6+vrm3zOjBkzkncfnEgiOVQCr/wTbHob+o2Fr78Ap1zod1UiMacwFpHEsPEteOEW7+tFl/0SpnwrbvvrRPymMBYRfzkHH90L79zlbQ3Pegx6jfC7KpG4UhiLiH+cgzfugEW/gzHXwswHIKuz31WJxJ3CWET84RzM+x4seQSm/T1c+rN2f0dXJNlph4yI+OPdn3pBfM4/Koilw1MYi0j8LX8SPvhfmHgjfOE/FcTS4SmMI9i1axfXXXcdI0aMYPTo0Vx++eVs3LiRTp06MWHCBE4//XSmTp3K448/fuw5jz32GH369Dl2GcUbb7zRxx6IJLAdC72vLw2/AL54j4JYBO0zbsQ5xzXXXMNNN93EnDlzAFixYgW7d+9mxIgRLF++HPBOh3nttdcSCAT4xje8k4595Stf4f777/etdpGEV3kAnr8FuhXAlx/3zgUtItoyDjd//nwyMzO59dZbj00bP348gwYNOqHd8OHDueeee7jvvvviXaJIcnIO5v4jlO+Ev/kj5HTzuyKRhKEwDrNmzRomTZrUqrYTJ05k/fr1xx4/88wzx4apH3300ViVKJKcVj0D6+bChT+Ggtb9HxPpKBJ3mPr1H8Ku1a1u3qm+zrskWnP6jYXLfn6ShR3n3ImXbdYwtUgTjpTBmz+Cgqlw9rf9rkYk4WjLOMyYMWMoLCxsVdvly5dz+umnx7gikRTw9o+h6hBcea9OcSkSQeJuGbdxC7YySpdQ/PznP8+PfvQjHn74Yb71rW8BsGTJEo4ePXpCu23btvG9732Pb39bn/JFmrX9E1jxJJz3L5DftkujinQUiRvGPjEzXnzxRf75n/+Zn//85+Tk5DB06FDuvfdeNm/ezIQJE6iqqiIvL49vf/vbx46kFpEInIO3/wPy+sP53/e7GpGEpTCOYMCAATz77LONpldWVjb5nJtvvpmbb745hlWJJKF1r0DxErjqNzrntEgztPNGRGKjvhb+ehf0OQ3Gfc3vakQSmraMRSQ2Vs6Bsk1w3dMtf9NBpIPTlrGIRF+gHj78NfQfB6de5nc1Igkv4cI4/Lu7qSSV+yZygnVzYf9mOO87Ove0SCskVBjn5ORQVlaWkqHlnKOsrIycnBy/SxGJLefgg3ug10g4/Uq/qxFJCgm1I6egoIDi4mL27t3b5udWVVUlfNDl5ORQUFDgdxkisbX5Xdi1CmY+AGnpflcjkhQSKowzMzMZNmxYu567YMECJkyYEOWKRKTNFv0euvSFsV/2uxKRpJFQw9QikuQObIPP3oJJN0NGlt/ViCQNhbGIRM+SP4ClwWSdmU6kLRTGIhIdtZWw/Ak47YvQdYDf1YgkFYWxiETH2heh8gBM/ZbflYgkHYWxiETH8qeg53AYOt3vSkSSjsJYRE7egW2w/UPvHNQ6yYdImymMReTkrXwGMBh3nd+ViCQlhbGInBznYOWfYdh06D7I72pEkpLCWEROzo5PvGHq8df7XYlI0lIYi8jJWf0cZHbWeahFToLCWETaL1AP616BUZdAVhe/qxFJWgpjEWm/7R/Bkb0w+mq/KxFJagpjEWm/tS96Q9QjL/a7EpGkpjAWkfY5YYi6s9/ViCQ1hbGItI+GqEWiRmEsIu2z9iUNUYtEicJYRNouEID1r3pBrCFqkZOmMBaRtitdDhW7vcslishJa1UYm9mlZrbBzDaZ2Q8jzO9mZq+Y2UozW2tmurK4SCrb+AZYGpxykd+ViKSEFsPYzNKBB4DLgNHAV81sdFiz24BPnXPjgBnAr8wsK8q1ikii2PgGDJoGnXv6XYlISmjNlvFUYJNzbotzrgaYA8wMa+OAPDMzIBfYD9RFtVIRSQyHSmDXKhh1qd+ViKSM1oTxQKAo5HFxcFqo+4HTgVJgNfBPzrlAVCoUkcTy2ZvercJYJGrMOdd8A7NZwCXOuVuCj28Apjrnvh3S5kvAucB3gBHA28A459zhsGXNBmYD5OfnT5ozZ07UOlJRUUFubm7Ulucn9SUxpUpfTrYfY1fdTeejRSya9nswi2JlbZcq6wTUl0QUi35ccMEFhc65yY1mOOea/QHOBt4MeXwHcEdYm9eA6SGP38UL7CaXO2nSJBdN8+fPj+ry/KS+JKZU6ctJ9aP6iHN393Vu3vejVs/JSJV14pz6kohi0Q9gqYuQia0Zpl4CjDSzYcGDsq4D5oa12QFcCGBm+cCpwJa2fmIQkQS37QOoq9KJPkSiLKOlBs65OjO7HXgTSAf+6Jxba2a3Buc/CNwNPGZmqwEDfuCc2xfDukXED5vfhYwcGHKu35WIpJQWwxjAOTcPmBc27cGQ+6WAPiqLpLrN73pBnJnjdyUiKUVn4BKR1jlUAvs2wogL/K5EJOUojEWkdbbM925HfN7fOkRSkMJYRFpn87uQmw99w0/AJyInS2EsIi0LBGDLAhg+w/fvFoukIoWxiLRs92o4WqYhapEYURiLSMs2v+vdDp/hZxUiKUthLCIt27LA21ec18/vSkRSksJYRJpXVwM7FsGw8/2uRCRlKYxFpHmly6CuUmfdEokhhbGING/bh96twlgkZhTGItK8bR96+4u79PK7EpGUpTAWkabV10LRIhh6nt+ViKQ0hbGINK10OdQeVRiLxJjCWESatu0D71b7i0ViSmEsIk3b9hH0OR269Pa7EpGUpjAWkcjqa2HHQhiqrWKRWFMYi0hkO1dC7RHtLxaJA4WxiESm/cUicaMwFpHIdiyEXiMht6/flYikPIWxiDQWCHjfLx48ze9KRDoEhbGINFb2GVQegEFn+V2JSIegMBaRxooWebeDtGUsEg8KYxFpbMci6NQTeo/0uxKRDkFhLCKNFS3ytorN/K5EpENQGIvIiY6UefuMB031uxKRDkNhLCInKl7s3Q7WwVsi8aIwFpET7VgIaZkwYILflYh0GApjETlR0SLoPw4yO/ldiUiHoTAWkePqaqBkmYaoReJMYSwix+1cCfXVOnhLJM4UxiJy3LGTfWjLWCSeFMYiclzRQugxFPLy/a5EpENRGIvIcUVLoEBD1CLxpjAWEc/hUqjYBQWT/a5EpMNRGIuIp6TQux0w0d86RDoghbGIeEoKIS0D+o31uxKRDkdhLCKekmWQfwZk5vhdiUiHozAWEQgEoHQ5DNQQtYgfFMYiAmWboPowDJzkdyUiHZLCWESgdJl3qzAW8UWrwtjMLjWzDWa2ycx+2ESbGWa2wszWmtl70S1TRGKqpBAyu0DvUX5XItIhZbTUwMzSgQeALwDFwBIzm+uc+zSkTXfgt8ClzrkdZtY3VgWLSAyUFHqXTExL97sSkQ6pNVvGU4FNzrktzrkaYA4wM6zN14AXnHM7AJxze6JbpojETF0N7Fqtg7dEfNSaMB4IFIU8Lg5OCzUK6GFmC8ys0MxujFaBIhJju9dAfY32F4v4yJxzzTcwmwVc4py7Jfj4BmCqc+7bIW3uByYDFwKdgE+ALzrnNoYtazYwGyA/P3/SnDlzotaRiooKcnNzo7Y8P6kviSlV+hLejwEl8xj12e/55KyHqc5Jrj1MqbJOQH1JRLHoxwUXXFDonGt0ztkW9xnjbQkPCnlcAJRGaLPPOXcEOGJm7wPjgBPC2Dn3EPAQwOTJk92MGTNa3YGWLFiwgGguz0/qS2JKlb406sdLz0KXPpx9ySww862u9kiVdQLqSyKKZz9aM0y9BBhpZsPMLAu4Dpgb1uZlYLqZZZhZZ2AasC66pYpITJQUekPUSRbEIqmkxS1j51ydmd0OvAmkA390zq01s1uD8x90zq0zszeAVUAAeMQ5tyaWhYtIFFSXw94NMOZavysR6dBaM0yNc24eMC9s2oNhj/8H+J/olSYiMVe6AnA6eEvEZzoDl0hH1nDZRH2tScRXCmORjqykEHoMhc49/a5EpENTGIt0ZKXLNUQtkgAUxiIdVfluOFSkMBZJAApjkY6q4UpNA7S/WMRvCmORjqpkGVg69D/T70pEOjyFsUhHVVIIfUdDVhe/KxHp8BTGIh2Rc94w9cAJflciIiiMRTqmA1uh8oAO3hJJEApjkY6oJHjwlsJYJCEojEU6opJCyOgEfU73uxIRQWEs0jGVLIP+4yC9VaenF5EYUxiLdDAWqIOdKzVELZJAFMYiHUyXIzugrlIXhxBJIApjkQ4mr3yTd0dhLJIwFMYiHUxe+Ubo1AN6DPO7FBEJUhiLdDBdD3/mnY/azO9SRCRIYSzSkdQc8fYZ6+AtkYSiMBbpSHauwggojEUSjMJYpCMpKfRudfCWSEJRGIt0JKXLqMruA7l9/a5EREIojEU6kpJCDncd6XcVIhJGYSzSURwpgwPbKM9TGIskGoWxSEdRuhxAYSySgBTGIh1FSSFglOeN8LsSEQmjMBbpKEoKoc9p1Gd09rsSEQmjMBbpCJyD0mX6SpNIglIYi3QEh4rgyF6FsUiCUhiLdATHTvahM2+JJCKFsUhHULIM0rOg7xi/KxGRCBTGIh1ByTLodyZkZPldiYhEoDAWSXWBeti5QkPUIglMYSyS6vZthJoKHbwlksAUxiKpruHgrQEKY5FEpTAWSXUlhZDdDXqd4nclItIEhbFIqiteCgMnQJr+u4skKv3vFElltZWwey0MnOx3JSLSDIWxSCrbuQpcvY6kFklwCmORVFay1LtVGIskNIWxSCorKYRugyAv3+9KRKQZCmORVFZSqO8XiySBVoWxmV1qZhvMbJOZ/bCZdlPMrN7MvhS9EkWkXY7sgwPbNEQtkgRaDGMzSwceAC4DRgNfNbPRTbT7BfBmtIsUkXYoWebd6khqkYTXmi3jqcAm59wW51wNMAeYGaHdt4HngT1RrE9E2qukECwN+o/zuxIRaUFGK9oMBIpCHhcD00IbmNlA4Brg88CUqFUnIu1XshT6nA7ZuRTtP8ofPtzK25/uZu/hSgYve48rzxzAzecMpVvnTL8rFenwzDnXfAOzWcAlzrlbgo9vAKY6574d0uY54FfOuYVm9hjwqnPuLxGWNRuYDZCfnz9pzpw5UetIRUUFubm5UVuen9SXxJRUfXGOcz+6gX29p/FQ7q08ta6GQADG9U2nW3odu6vTWVsWoFu28Q/jsjm1Z7rfFbdLUq2TFqgviScW/bjgggsKnXON9x0555r9Ac4G3gx5fAdwR1ibrcC24E8F3lD11c0td9KkSS6a5s+fH9Xl+Ul9SUxJ1Zeyzc79pKub/9Qv3ZAfvOquf3ihKz141Dl3vB+rig66Gf8z353yo9fcu+t3+1hs+yXVOmmB+pJ4YtEPYKmLkImt2We8BBhpZsPMLAu4DpgbFujDnHNDnXNDgb8A/+Cce6ntnxlEJCqKvSs1/WJ1F64aN4BHvzGF/t06ndBkbEE3XrrtXE7tl8etTxSyouigH5WKCK04gMs5VwfcjneU9DrgWefcWjO71cxujXWBItJ2ZRs/4ajLpufQcfzqy+PITI/8X71bp0z+9LfT6J2bzW1PLePAkZo4Vyoi0MrvGTvn5jnnRjnnRjjn/is47UHn3IMR2t7sIuwvFpH4qK0PsHvdR6y34fzm61OaDOIGPbtk8dvrJ7KnvIq7XlkbpypFJJTOwCWSYh6ev54RdZvpMfJsenbJatVzxg3qzt9/bgQvrSjlw8/2xbhCEQmnMBZJIaUHK3lnwXyyrZZh46a36bn/cMEpDO3VmX9/aTXVdfUxqlBEIlEYi6SQe97eyJm2yXvQxtNg5mSmc9fMM9hWdpQ5i4tafoKIRI3CWCRFrN91mOeXFfM3fUshNx+6D27zMs4f2ZuzhvfkN+9u4mhNXQyqFJFIFMYiKeK+v35GbnYGo+vXQ8EUMGvzMsyMf73kNPZVVPPoR9uiX6SIRKQwFkkBW/cd4fU1u/i7SXmkH9wGg6a1+JymTBrSgwtO7cMfP9xKVa32HYvEg8JYJAU89P5mMtPTuKEgeJ2WQVNPanl/97kRlB2p4YVlJVGoTkRaojAWSXJ7DlfxfGEJsyYV0G3fMkjLhP7jT2qZ04b1ZOzAbjzywRYCgebPXy8iJ09hLJLkHv9kG3WBALPPHw5FS2DAeMjMOallmhnfOn84W/Yd4a/rdVVUkVhTGIsksZq6AM8sKeLC0/MZ0i0TSpdBwckNUTe4/Ix+DOzeiT98uCUqyxORpimMRZLYG2t3sa+ihq+fNQR2rYa6qpPeX9wgIz2N688azMIt+9m8tyIqyxSRyBTGIknsyU+2M6RXZ6af0huKF3sToxTGALMmDSIjzXh60Y6oLVNEGlMYiySpDbvKWbxtP1+bOpi0NIOiRdBtEHQdELXX6JOXzSVj+vGXZcX6mpNIDCmMRZLUnCU7yEpPY9bkQd6EoiXeyT6i7GvTBnPwaC1vrNkV9WWLiEdhLJKE6uoDvLJyJ58/ra93ZaZDxXC4+KQldhEAABr8SURBVKRO9tGUs4f3Ymivzjy9WEPVIrGiMBZJQgu37GdfRTUzxweHpIuiv7+4QVqa8aVJBSzaup+i/UejvnwRURiLJKWXV5SQl53BBaf19SYUL4GMTtBvbExe7+oJAwF4abnOyCUSCwpjkSRTVVvPG2t2cfGYfuRkpnsTd3wCAydCemZMXrOgR2fOGt6TF5aX4JzOyCUSbQpjkSSzYMNeyqvrjg9RV5fDzpUw5JyYvu61EwrYuu8IK4oOxvR1RDoihbFIkpm7soTeuVmcM6KXN6FoEbhAzMP4srH9yM5I08UjRGJAYSySRMqranln3R6+OLY/GenB/77bP4a0jJgcSR0qLyeTS8b045VVpdTUBWL6WiIdjcJYJIm8tXY3NXUBrho/8PjE7R97V2nK6hLz179m4kAOHq1l/gZdPEIkmhTGIknk5ZWlFPToxMTB3b0JtZVQUhjzIeoG00/pTe/cbF7UULVIVCmMRZLEvopqPtq0j6vGDcDMvIklhVBfA0POjUsNGelpXHFmf+Zv2EN5VW1cXlOkI1AYiySJeat3Uh9wzAwfosZgcGz3F4e6clx/qusCvLNud9xeUyTVKYxFksTLK0o5rV8ep/bLOz5x24eQfwZ06hG3OiYM6sGAbjm8unJn3F5TJNUpjEWSQNH+oxRuP8CV40KuyFRX450GM077ixukpRlXjBvA+5/t5eDRmri+tkiqUhiLJIFXVpUCcFVoGO9cCXWVcQ9jgCvPHEBtvePNtbqSk0g0KIxFksDcFaVMGtKDQT07H5+4/SPv1ocwPmNgV4b26swrGqoWiQqFsUiC27CrnPW7yk/cKgbY+j70HgW5feNek5lxxZkD+HjzPvaWV8f99UVSjcJYJMHNXVlCeppx+dj+xyfW1XgXhxg+w6+yuHLcAAIO3lijrWORk6UwFklgzjleXlHKuaf0pk9e9vEZxUug9igM+5xvtZ3aL49R+bkaqhaJAoWxSAJbtuMgxQcqGw9Rb1kAlgZDz/OlrgZXnjmAxdv2s/NQpa91iCQ7hbFIAntlZSlZGWlcMib/xBlb34MBE6BTd38KC7oi+CHhtVXaOhY5GQpjkQRVVx/g1VWlXHR6X/JyMo/PqDoMxUt9HaJuMKx3F84Y2JVXVpb6XYpIUlMYiySojzeXsa+ipvEQ9faPwdX7evBWqCvPHMDK4kNsLzvidykiSUthLJKg5q4sJS87gxmnhn11aet7kJET8+sXt1bDULW2jkXaT2EskoCqaut5Y80uLj2jHzmZ6SfO3PKeF8SZOf4UF2Zg905MHtKDuQpjkXZTGIskoPnr91BRXXfiFZoADu+EPWsTZoi6wVXjB7BxdwXrdx32uxSRpKQwFklAc1eW0js3m7NH9DpxxqZ3vNuRF8e/qGZcPrY/6WnG3BXaOhZpj1aFsZldamYbzGyTmf0wwvzrzWxV8OdjMxsX/VJFOobDVbX8df0erjjTC7gTfPYW5A2A/DH+FNeE3rnZnDOiF6+sKsU553c5IkmnxTA2s3TgAeAyYDTwVTMbHdZsK/A559yZwN3AQ9EuVKSjeHPNLmrqAswcH3YUdX0tbJ4PI78AZpGf7KOrxg2gaH8ly4sO+l2KSNJpzZbxVGCTc26Lc64GmAPMDG3gnPvYOXcg+HAhUBDdMkU6jrkrSxncszPjB4Wd0GPHQqgpT7gh6gaXnNGPrIw0DVWLtENrwnggUBTyuDg4rSnfBF4/maJEOqq95dV8tGkfV47rj4Vv/W56G9IyYbj/J/uIpGtOJhec2ofXVu+kPqChapG2sJb275jZLOAS59wtwcc3AFOdc9+O0PYC4LfAec65sgjzZwOzAfLz8yfNmTPn5HsQVFFRQW5ubtSW5yf1JTHFoy9vbavlz+tr+K/zOjEw98TPypOX/CO1md1YOf7uk3qNWPZj8a46fruimu9PyWF0r/SWn3CS9PuVmFKlL7HoxwUXXFDonJvcaIZzrtkf4GzgzZDHdwB3RGh3JrAZGNXSMp1zTJo0yUXT/Pnzo7o8P6kviSkefbnyNx+4L973fuMZB3Y495Ouzn30m5N+jVj2o7Kmzo3+j9fdvz63ImavEUq/X4kpVfoSi34AS12ETGzNMPUSYKSZDTOzLOA6YG5oAzMbDLwA3OCc29jeTwwiHdmmPRWsKj7E1eHfLQbvKGrwDt5KYDmZ6Vw8ph+vr9lFdV293+WIJI0Ww9g5VwfcDrwJrAOedc6tNbNbzezWYLMfA72A35rZCjNbGrOKRVLUS8tLSDMan4saYP2r0HM49B4V/8La6KpxAyivquO9DXv9LkUkaWS0ppFzbh4wL2zagyH3bwFuiW5pIh1HIOB4aUUJ557Sm75dw05zWXkQtr4PZ9+WkF9pCnfeyN706JzJ3JWlXDymn9/liCQFnYFLJAEU7jhA8YFKrp3YxBB1oA5OuzL+hbVDZnoal4/tzzvrdlNeVet3OSJJQWEskgBeXF5Cp8x0Lh4dYUty3SuQ2w8GTop/Ye30N5MKqKoNMG/1Tr9LEUkKCmMRn1XX1fPaqp1cMiafLtlhe45qK73zUZ/2RUhLnv+uEwZ1Z3ifLjy3tNjvUkSSQvL87xZJUfPX7+FQZS1XT4gwRL35Xag9CqcnxxB1AzNj1qRBLN1+gC17K/wuRyThKYxFfPbMkiLyu2Zz3im9G89c+yJ06gFDz4t/YSfp2okDSTN4fpm2jkVaojAW8dHOQ5W8t3EvX5pUQEZ62H/HmiOw/jUYfTWkZ/pT4EnI75rD+aP68HxhiU6PKdIChbGIj/6ytJiAgy9PHtR45obXvSHqsbPiX1iUzJo0iF2Hq/ho0z6/SxFJaApjEZ8EAo5nC4s4Z0QvhvTq0rjB6r9A14Ew+Oz4FxclF57el26dMnl2aVHLjUU6MIWxiE8+2VJG0f5KvjIlwlbx0f3eUdRnXJtUR1GHy8lM55oJA3lr7W72VVT7XY5Iwkre/+UiSW7OkiK6dcrkkkhnqfr0ZQjUwhlfin9hUfb1swZTUx/Q1rFIMxTGIj7YU17FG2t2cs2EgeRkRrjU4PInoc9p0H9c/IuLslP65nHW8J48tXCHDuQSaYLCWMQHf160g9p6x41nD2k8c/daKFkKE29KinNRt8YNZw2l5GAl723c43cpIglJYSwSZzV1AZ5atIMZp/ZheJ8IFy5f9gSkZ8GZX4l/cTFy8Zh8+uRl88Qn2/0uRSQhKYxF4uz1NTvZW17NTecMbTyztgpWzYHTroAuveJeW6xkpqfx1amDWbBxL9vLjvhdjkjCURiLxNljH29jWO8ufG5kn8Yz178KlQdg4o3xLyzGrp82mIw0448fbvW7FJGEozAWiaMVRQdZvuMgN549hLS0CPuDFz8EPYbCsM/FvbZYy++aw8zxA3lmaRH7j9T4XY5IQlEYi8TR7xZsomtOBl+aVNB4ZnEhFC2Cabcm9XeLmzP7/OFU1Qa071gkTGr+jxdJQJv2lPPm2t3cdM5Q8nIinGt64W8huytM+Hr8i4uTUfl5fP60vjz+yTaqauv9LkckYSiMReLkwfe2kJOZxs2RDtw6VAKfvgQTboDsvLjXFk+zzx/O/iM1PKeTgIgcozAWiYOSg5W8tLyE66YMplduduMGi38PLgDTZse/uDibNqwn4wd158H3tlBdp61jEVAYi8TFfe98RpoZs88f3njmkTJY/AiMucY7eCvFmRnf+cIoSg5W8uwSbR2LgMJYJOY2763gucIirj9rMAO6d2rc4JPfeJdK/NwP4l+cT6aP7M2UoT24f/4m7TsWQWEsEnP3vL2RnMx0brvglMYzj5TBoofgjL+BPqfGvzifmBnfvfhUdh+u5smFOrJaRGEsEkOriw/x2qqdfPO8YfSOtK/44/8LbhV/P/7F+eys4b2YPrI3v3l3Ewf0vWPp4BTGIjESCDh+MncNvXOz+FakfcX7t8DC38G4r3aoreJQ//7F0ZRX1fLrdzb6XYqIrxTGIjHywvISlu04yA8uPY2ukb5X/NZ/QFomXPjj+BeXIE7tl8fXzxrCkwu3s37XYb/LEfGNwlgkBg5X1fLz19czYXB3/mZihLNtbX3fOw/19H+Brv3jX2AC+ZeLRpGXk8mPX15LQNc7lg5KYSwSA//92jr2H6nmP686o/E5qGsr4dXvQPfBcPbt/hSYQHp0yeKOy05j8db9PLV4h9/liPhCYSwSZfM37GHOkiJmnz+CsQXdIjT4byj7DK68DzIjfNWpA/rKlEFMH9mbn89bR/GBo36XIxJ3CmORKDp4tIYfPr+KUfm5/MsXRjZuULQEPrkfJt0MIy6Ie32Jysz42bVjAfjecyup13C1dDAKY5EoCQQc3312JWUVNdzz5fFkZ6Sf2KDyADz/TcgbAF+4258iE1hBj87cedUYFm7Zz706ulo6GIWxSJTcP38Tf12/hx9fOZozBoYNTwcC8MJsOFwKsx6DnK6+1JjoZk0exKxJBfzm3U3M37DH73JE4kZhLBIFb3+6m1+/s5FrJwzkhrOGNG6w4Gfw2Vtw6c9g0JT4F5hE/nPmGZzWL49/fHo5G3aV+12OSFwojEVO0qItZdz+52WcWdCd/7pmLGZhR08veQTe/yWMvx6m3OJPkUmkU1Y6j9w0mU6Z6dz0x8WUHqz0uySRmFMYi5yEVcUHueXxpRT06MSjN0+hU1bYfuLVf4HXvgejLoMr/w/Cg1oiKujRmce+MZUj1XXc8IdF7DpU5XdJIjGlMBZpp4837eOrDy2ka6dMnvjmNHp2yTqxwZJH4PlbYMi5MOtRSI9wFi5p0ugBXfnDzVPYfbiaWb//mB1l+sqTpC6FsUg7vLCsmJsfXcLAHp14/u/POfHSiIF6eOcueO27MOpSuP45fZ+4naYO68lTt0zjcGUd1/7uYxZv3e93SSIxoTAWaYOq2nrueGE133l2JRMGd+fZvzubft1yjjco3w1/mgkf3gMTb4KvPAlZnf0rOAWMG9Sdv9x6Nnk5GXzt4YX84cOtOm2mpByFsUgrfVpWz2X/9wFPL97BrZ8bwVO3TKN75+DQdCAAhY/DA1OheCnM/C1cdR+kZ/hbdIoYmZ/HS7edy4xT+3D3q59y3cML2bK3wu+yRKJGfylEWrBxdzm/fnsjr6+pYkivzjz5zWmcN7K3N9M52PSOd4rL0mXe/uErft1hL4kYS906ZfLwjZN5rrCYn776KZfc+z6fG5jOmEnV9MmLcK1okSSiMBaJoD7geH/jXp5atIO/rt9Nl6wMZo7I5Bc3n09OZjpUHYY1f4HCx2DnSug2GK4OXptYR0zHjJnx5cmDmDGqD79+5zOeWbKD837xLleNG8CNZw+NfC5wkSTQqjA2s0uB/wPSgUeccz8Pm2/B+ZcDR4GbnXPLolyrSEwdqa5jybb9vLNuN+98uoddh6vonZvFbTNO4ZvnDWP9B3PJWfUEbHwLtsyH2qPQdzRcca/3HeKMrJZfRKKib9ccfnbtWMZl72V1TR9eXF7Cc4XFjOjThUvG9OPzp/VlbEG3xqckFUlQLYaxmaUDDwBfAIqBJWY21zn3aUizy4CRwZ9pwO+CtyIJp7Y+QMmBSrbuO8LWfUfYtLeC5TsOsmHXYQLO0T+riqsHVXHJ+ErGdtpPxu7H4MFlnF1e6i2g2yAY/zVvK3jgJG0J+6hflzSu++JYfnDZaby8opQ31uzk9+9v4bcLNpOVkca4gm6cMbAbp/TNZUSfXIb26kLv3Cwy0nW4jCSW1mwZTwU2Oee2AJjZHGAmEBrGM4E/OeccsNDMuptZf+fczqhXHMH+PSUc2vQxy2uKm23nlXciIxDWBiDCkZphz42wKIiwLAtblgv598SJx6cdKi5hacXGyC8SMs0RYfmR6o8wLfy5xxcbqWOu0d0TnhtWV2jjw7t2sbhs+bHWFlZ/2Is3fk3XeFropIY6nPNCtq4+QG3wp67eBe/XU1NdRV31UeprKrG6KnKoIYca+lgNp6UfZXbmEXrlldO5/hDp9dVQgvcD0HMEDD2XTUe7csrFt0Df0xXACaZrTiY3nDWEG84awoEjNSzaup/C7ftZsu0AcxYXUVlbf6ytGfTsnEWfvGx65WbRJSuD3OwMOmen0yU7g86ZGWRmGBlpRnpaGpnpRnqakZmWRnqakZFux86yZiHL9B5b2OMT59PkfGPNnjpqP90d9ffGD6nSlzV76ji3PkBmHD68tSaMBwJFIY+LabzVG6nNQOCEMDaz2cBsgPz8fBYsWNDGciM7tH0FM4t/4b1qCpgMsMvvKqJor98FHBfAqM3Mot6yCKRnQXo29Zm51Gb1pCJzCAcyu1GT1YPKTv2o7NSfqpx8AunewUEV2RUUr9sD65L7AgYVFRVR+7/nt6b6kgOc2wXOHQOB0dkcqHLsPBJgz1HHoWrHoeoAh2qOsnvfEarqHFX1UF3nqKyHukCjxcXPsqU+vniUpUhfTu3xHp0zY//huzVhHKmK8E2Z1rTBOfcQ8BDA5MmT3YwZM1rx8i2rODyOV+Z14fTTR59YVKStF2v8CSe8XeSNHms0I9LyLfytsGP/tLCs4/dXrVrNmWeODT7VTnxe2LQTP5+HzA5fPoaL1C7SqgtdfsTZxz7Pn1CHi/C8ZYXLmDhp0rEaItffxLRj7SPXFjIRM8jJTCc7PY2M9Aj7CdMzIbMTaelZZLdzq3bBggVE63fWT6nSD4hNX+rqA9QFHHUBR329oy5w/LE38tLwp+34qMzxR6GPw+Y3Mb3B0qVLmTx5clT74pdU6cvSpUu55MIZpKclRhgXA4NCHhcApe1oEzO5XXuQ128kp4w7N14vGVNbSg8yeNR4v8uIik6biskvGOF3GSKtlpGehh/Hfe37LL3xpTeTVKr0Zd9n6XEJYmjdST+WACPNbJiZZQHXAXPD2swFbjTPWcCheO0vFhERSXYtbhk75+rM7HbgTbyvNv3RObfWzG4Nzn8QmIf3taZNeF9t+kbsShYREUktrfqesXNuHl7ghk57MOS+A26LbmkiIiIdg75sJyIi4jOFsYiIiM8UxiIiIj5TGIuIiPhMYSwiIuIzhbGIiIjPFMYiIiI+s0hXMorLC5vtBbZHcZG9gX1RXJ6f1JfElCp9SZV+gPqSqFKlL7HoxxDnXJ/wib6FcbSZ2VLnXPKfmRz1JVGlSl9SpR+gviSqVOlLPPuhYWoRERGfKYxFRER8lkph/JDfBUSR+pKYUqUvqdIPUF8SVar0JW79SJl9xiIiIskqlbaMRUREklJShbGZzTKztWYWMLPJYfPuMLNNZrbBzC5p4vk9zextM/sseNsjPpU3z8yeMbMVwZ9tZraiiXbbzGx1sN3SeNfZGmZ2p5mVhPTn8ibaXRpcV5vM7IfxrrM1zOx/zGy9ma0ysxfNrHsT7RJyvbT0HpvnvuD8VWY20Y86W2Jmg8xsvpmtC/7//6cIbWaY2aGQ37sf+1Fra7T0+5IM68XMTg15r1eY2WEz++ewNgm7Tszsj2a2x8zWhExrVT7E7G+Xcy5pfoDTgVOBBcDkkOmjgZVANjAM2AykR3j+L4EfBu//EPiF332KUOOvgB83MW8b0NvvGluo/07gey20SQ+uo+FAVnDdjfa79gh1XgxkBO//oqnfl0RcL615j4HLgdcBA84CFvlddxN96Q9MDN7PAzZG6MsM4FW/a21lf5r9fUmW9RJSbzqwC+/7s0mxToDzgYnAmpBpLeZDLP92JdWWsXNunXNuQ4RZM4E5zrlq59xWYBMwtYl2jwfvPw5cHZtK28fMDPgy8LTftcTYVGCTc26Lc64GmIO3bhKKc+4t51xd8OFCoMDPetqoNe/xTOBPzrMQ6G5m/eNdaEucczudc8uC98uBdcBAf6uKqaRYLyEuBDY756J5EqeYcs69D+wPm9yafIjZ366kCuNmDASKQh4XE/k/a75zbid4/8GBvnGorS2mA7udc581Md8Bb5lZoZnNjmNdbXV7cHjtj00M9bR2fSWSv8XbWokkEddLa97jpFsPZjYUmAAsijD7bDNbaWavm9mYuBbWNi39viTbermOpjcgkmWdQOvyIWbrJiMaC4kmM3sH6Bdh1r85515u6mkRpiXUYeKt7NdXaX6r+FznXKmZ9QXeNrP1wU94cdVcX4DfAXfjvf934w27/234IiI815f11Zr1Ymb/BtQBTzWxmIRYL2Fa8x4nzHpoDTPLBZ4H/tk5dzhs9jK8YdKK4HEKLwEj411jK7X0+5I068XMsoCrgDsizE6mddJaMVs3CRfGzrmL2vG0YmBQyOMCoDRCu91m1t85tzM47LOnPTW2R0v9MrMM4FpgUjPLKA3e7jGzF/GGTOL+R7+168jMHgZejTCrtesr5lqxXm4CrgAudMGdRhGWkRDrJUxr3uOEWQ8tMbNMvCB+yjn3Qvj80HB2zs0zs9+aWW/nXMKdH7kVvy9Js16Ay4Blzrnd4TOSaZ0EtSYfYrZuUmWYei5wnZllm9kwvE9fi5tod1Pw/k1AU1vafrgIWO+cK44008y6mFlew328g4vWRGrrp7B9W9cQucYlwEgzGxb8ZH0d3rpJKGZ2KfAD4Crn3NEm2iTqemnNezwXuDF49O5ZwKGGYbpEEjyW4g/AOufcPU206Rdsh5lNxfvbVha/Klunlb8vSbFegpoczUuWdRKiNfkQu79dfh7R1tYfvD/uxUA1sBt4M2Tev+Ed5bYBuCxk+iMEj7wGegF/BT4L3vb0u08hdT4G3Bo2bQAwL3h/ON6ReyuBtXjDqL7XHaEfTwCrgVXBX9L+4X0JPr4c76jYzQncl014+4dWBH8eTKb1Euk9Bm5t+D3DG3J7IDh/NSHfUEikH+A8vKHAVSHr4vKwvtwefP9X4h1sd47fdTfRl4i/L0m6XjrjhWu3kGlJsU7wPkDsBGqDmfLNpvIhXn+7dAYuERERn6XKMLWIiEjSUhiLiIj4TGEsIiLiM4WxiIiIzxTGIiIiPlMYi4iI+ExhLCIi4jOFsUgHYWZTghfwyAmeCWqtmZ3hd10igk76IdKRmNlPgRygE1DsnPuZzyWJCApjkQ4leD7dJUAV3ukJ630uSUTQMLVIR9MTyAXy8LaQRSQBaMtYpAMxs7nAHGAY3kU8bve5JBEhAa9nLCKxYWY3AnXOuT+bWTrwsZl93jn3rt+1iXR02jIWERHxmfYZi4iI+ExhLCIi4jOFsYiIiM8UxiIiIj5TGIuIiPhMYSwiIuIzhbGIiIjPFMYiIiI++3+YeMZURgQN8AAAAABJRU5ErkJggg==\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "plot_normal_pdf_cdf(mu=0, sigma=1)" ] }, { "cell_type": "code", "execution_count": 5, "metadata": { "slideshow": { "slide_type": "subslide" } }, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAeMAAAFNCAYAAAAzYQemAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8QZhcZAAAgAElEQVR4nO3deZxcdZ3v/9enq/cte0I2SIIgq0QScQNtl58CV4erdxwBHyKMyuXewd94Ha+oMz+XUefqLF5/XnSYOCKuoDOOihjEjVYUUUADBMISQkI6CUln6+7q6u7avvePc6pT3anldKeqTtXp9/Px6FRXnVPf8/3WSfe7v9/zPeeYcw4REREJT1PYFRAREZnrFMYiIiIhUxiLiIiETGEsIiISMoWxiIhIyBTGIiIiIVMYy5xnZh8zs2+EXY8gzGyZmf3KzEbM7J/Crs9MmVmfmQ2cYBl3mtk7KlSfi8zsibznO83stZUo2y/vUTPrq1R5El0KY6kY/xfZmJnFzWy/mX3FzLr9Zf1mNu6HyLCZPWhmHzSztrz3f8zMUv77c18fCK9Fk/XqNbPPmdmzfp22+88X+8tz7R4xs6Nmdq+ZXWdmTXll3GJmyWlte+ssqnMtcBDodc79VYG63mJmzswuyHvteWbWEBcU8Os+6n8+h8zs59M/J+fcJc65rwYs63ml1nHO3eOce/6J1tvf3i1m9slp5Z/tnOuvRPkSbQpjqbQ3Oue6gfOBFwF/k7fseudcD7Ac+CvgcmCzmVneOt92znXnff19zWpegJm1Aj8HzgYuBnqBlwGHgAvyVn2j37ZTgE8DNwBfnlbc309r27dnUaVTgMdc6av1HAY+WWJ5YGbWXIlyZug8///Q84FbgBvN7KOV3khIbRMpSGEsVeGc2wPcCZxTYNmo31v4E+ClwH+aafl+r/ppvzf6mJm9KW/Z1Wb2azP7RzM7YmbPmNklecvXmtkv/ff+FFhcYlNXAScDb3LOPeacyzrnDjjnPuGc21ygbUPOuduBtwLvMLPj2h+gbS8zs/vNbMh/fJn/+i3AO4AP+D3HYsOpXwVeYGavLFL+CjO73cwO+738d+ct+5iZ/buZfcPMhoGr/df+zX9txMweMbPTzexDZnbAzHab2evyyrjGzLb56+4ws/86088AwDl30Dn3deC/AR8ys0V++f1m9i7/++f5+3LIzA6a2bf913/lF/NQbhTC/CFyM7vBzJ4DvmKFh81f5P+fOuKP7rT7ZV5tZr+e9lk6vw7XAm/j2L75ob98ctjbzNr8EZW9/tfnzB8ZyqvbX/mf6T4zu2Y2n5s0JoWxVIWZrQYuBf5YbB3n3LPAA8BFs9jE0/775gEfB75hZsvzlr8YeAIvaP8e+HJeD/xbwIP+sk/gBVwxrwV+7JyLz6RyzrnfAwPMsG1mthD4EfB5YBHwWeBHZrbIOXc18E2O9bB/VqSYBPB3wKeKLL/Vr9sK4E+BvzOz1+Qtvwz4d2C+vz2ANwJfBxbg7dO78H5/rAT+FviXvPcfAN6AN4pwDfC/zez8IO0v4gdAM1NHInI+AfzEr9cq4P8AOOde4S8/b9ooxEnAQrwRhmuLbO9twOuBU4HTmTq6U5BzbhNT980bC6z218BLgPXAeX578ss+Ce//80rgncAXzGxBuW1LNCiMpdK+b2ZHgV8Dv8QLhVL24v1yzPkz84675r5WFHqTc+7fnHN7/Z7qt4GnmPrLepdz7kvOuQxeT3E5sMzMTsYbPv//nHMTzrlfAT8sUb9FwL4ybQjatvfntetgkff8J+Ap59zXnXNp59ytwON4YTgT/wKcnD8iAJN/JF0I3OCcG3fObQH+FXh73mq/dc593/9sx/zX7nHO3eWcSwP/BiwBPu2cSwG3AWvMbD6Ac+5HzrmnneeXeGE5mz+48MtL4R0nX1hgcQovWFf47fl1gXXyZYGP+vt+rMg6NzrndjvnDuP9QXPFbOs+zduAv/VHVgbx/ojM/9xT/vKUP+oSxxuqlzlAYSyV9p+dc/Odc6c45/57iV94OSvxjnHmfMd/f+5rb6E3mdlVZrYlF254w+H5w83P5b5xziX8b7vxeoNHnHOjeevuKlG/Q3hBPhvT2/aPee0qNjS+okB9dvllBeacm8DrNX4CyD8mvwI47JwbKVH+7gJF7s/7fgw46P+hk3sO3ueLmV1iZvf5w+BH8UZISh0KKMnMWvDC/3CBxR/Aa9/vzZu5/Odliht0zo2XWSe//bvwPrNKmL5vp5d9yP9jJyeB/5lK9CmMJTR+L20DcM8M33cK8CXgemCRc24+sJWpoVPMPmCBmXXlvXZyifV/Brx+2vpB6vgivIAr11Obbi9eTy/fycCeGZYD8BW8Yc835b22F1hoZj0lyp/1zGv/GOh3gX8Elvn7ZjPB9k0xlwFp4PfTFzjnnnPOvds5twL4r8AXrfQM6iBtW533/cl4nxnAKNCZW2BmJ82w7On7Nr9smeMUxlJzZtbpTy76Ad4v2OMmQpXRhfeLb9Av7xoKTBQrxDm3C+849cfNrNXMLqT0EPDX8XpK3zWzM8ysycwWmdmHzezS6SubdxrUG/CGbr/hnHtkRi3zPovTzexKM2s277Ses4A7ZlgOfi/rY3gzu3Ov7QbuBf6XmbWb2Qvwjk9+s2AhM9cKtOHtm7Q/TP660m8pzMwWmtnbgC8An3HOHSqwzlvMbJX/9Aje/4tcj30/sG4Wm/4LM1vlH7//MJA73vwQcLaZrfcndX1s2vvKbe9W4G/MbIl5p8V9BGiI89ul+hTGUks3mtkI3i+tz+H1oC52zmVnUohz7jHgn4Df+mWdC/xmBkVciTfB6zDwUeBrJbY1gTeJ63Hgp8Aw3h8Qi4Hf5a36Q79tu/Em6nwWb/LSjPiB8wa8U78O4Q3DvsE5V+wYczm3cvwx7yuANXi9su/hHUP96SzLn8If/v5/ge/gheOVwO0zLOYhM4sD24F3Af/DOfeRIuu+CPidv/7twF86557xl30M+Kp/KOPPZrD9b+Ed597hf30SwDn3JN5ktZ/hzVGYPurxZeAsf3vfL1DuJ/H+EHwYeAT4AxU6BU0an5U+XVFERESqTT1jERGRkCmMRUREQqYwFhERCZnCWEREJGQKYxERkZCFdteSxYsXuzVr1lSsvNHRUbq6ZnRdhrqlttSnqLQlKu0AtaVeRaUt1WjHgw8+eNA5t2T666GF8Zo1a3jggQcqVl5/fz99fX0VKy9Makt9ikpbotIOUFvqVVTaUo12mFnBy+9qmFpERCRkCmMREZGQKYxFRERCFtoxYxERmbtSqRQDAwOMj5e7o2V45s2bx7Zt22b13vb2dlatWkVLS0ug9RXGIiJScwMDA/T09LBmzRrMTuQOm9UzMjJCT09P+RWncc5x6NAhBgYGWLt2baD3aJhaRERqbnx8nEWLFtVtEJ8IM2PRokUz6vUrjEVEJBRRDOKcmbZNYSwiInNSLBZj/fr1nHPOObzlLW8hkUhMef2CCy7gvPPO47Of/SzZrHfb9f7+fubNm8f69etZv349r33taytSFx0zFhGROamjo4MtW7YA8La3vY2bbrqJ973vfZOvj4yMMDY2xpVXXsnQ0BAf//jHAbjooou44447KloX9YxF5pADw+P8Zk+K8VQm7KqI1JWLLrqI7du3H/f60qVL2bRpEzfeeCPOuaptX2EsMod87IeP8qVHknzttzvDropI3Uin09x5552ce+65BZevW7eObDbLgQMHALjnnnsmh6k/9alPVaQOGqYWmSPSmSy/eNz7ZfLjrc9x7StODblGIp6P//BRHts7XNEyz1rRy0ffeHbJdcbGxli/fj3g9Yzf+c53Fl03v1dcjWFqhbHIHLHj4CjjqSw9LbBt3wiZrCPWFN3ZrCLl5B8zLmXHjh3EYjGWLl0664uAlKMwFpkjtu3zeh4vW9nMXTvT7Do0yrol3SHXSoSyPdgwDQ4Oct1113H99ddX9VQshbHIHLHn6BgAZy+KcdfONLuPjCmMRQrIDV9PTEzQ2trK29/+dt73vvdVdZsKY5E5Yv/QOD3tzazq8eZtDhxJhFwjkXDF4/GCr2cy3tkGhS6H2dfXV5V7NWs2tcgcsW9onOXz2pnfZrTEjN2Hx8Kukoj4FMYic8Rzw+OcNK+DJjOW9rQzODIRdpVExKcwFpkj9g2Ns7y3HYBF3a0cGlUYi9QLhbHIHJDKZDkYn+CkeX4Yd7VyKJ4MuVYikqMwFpkDDoxM4BzHwri7jUNx9YxF6oXCWGQOyB0fXtrTBnjD1AdHk1W91q6IBKcwFpkDjia8Ien5na0ALO5qI5nOEp9Ih1ktkdA999xzXH755Zx66qmcddZZXHrppTz55JN0dHRw4YUXcuaZZ3LBBRfw1a9+dfI9t9xyC0uWLJm8PvVVV111wvXQecYic8DRRAqA+Z0tjOD1jAEOxpP0tLeEWDOR8DjneNOb3sQ73vEObrvtNgC2bNnC/v37OfXUU/n1r39NT08PO3bs4M1vfjPZbJZrrrkGgLe+9a3ceOONFauLesYic8Bkz7jDC95F3d5wtY4by1x2991309LSwnXXXTf52vr161m9evWU9datW8dnP/tZPv/5z1etLgpjkTng6JjXM56XC+OuYz1jkblq69atbNiwIdC6559/Po8//vjk829/+9uTw9Rf+cpXTrguGqYWmQOOJlL0tDfTHPP+/s4NUx8eVRhLHbjzg/DcI5Ut86Rz4ZJPV6y46ZMdNUwtIjN2NJFkgT95C471kEfGU2FVSSR0Z599Ng8++GCgdf/4xz9y5plnVq0u6hmLzAFHEinmdx6bqNXREqO5yRhWGEs9qGAPdiZe/epX8+EPf5gvfelLvPvd7wbg/vvvJ5GYehOVnTt38v73v5/3vOc9VauLwlhkDjg6lpo8rQnAzOjtaGF4TKc2ydxlZnzve9/jve99L5/+9Kdpb29nzZo1fO5zn+Ppp5/mwgsvJJlM0tPTw3ve857JmdTVoDAWmQOGEklOWdg55bXe9mb1jGXOW7FiBd/5zneOe31sbKzgLRQBrr76aq6++uqK1qPsMWMzu9nMDpjZ1iLL32ZmD/tf95rZeRWtoYicsCOJFAs6p55P7PWMFcYi9SDIBK5bgItLLH8GeKVz7gXAJ4BNFaiXiFRINusYHk9NTtrK6W1vYXhcw9Qi9aBsGDvnfgUcLrH8XufcEf/pfcCqCtVNRCogkcrgHMddaau3o1k9Y5E6YUEuFG9ma4A7nHPnlFnv/cAZzrl3FVl+LXAtwLJlyzbkLj9WCfF4nO7u7oqVFya1pT41aluOjGf5H/1jXH12K32rWybbcfPWCR4ezPC5V3WWL6RONeo+KWSutWXevHmceuqpmFmNajVzmUyGWCw2q/c653j66acZGhqa8vqrXvWqB51zGwu+odwXsAbYWmadVwHbgEVBytywYYOrpLvvvrui5YVJbalPjdqWp/aPuFNuuMP9YMse59yxdnzqR4+5M/7mzhBrduIadZ8UMtfasmPHDjc4OOiy2Wz1KzRLw8PDs3pfNpt1g4ODbseOHcctAx5wBTKxIrOpzewFwL8ClzjnDlWiTBGpjNydmXrapv6497Y3M5bKkExnaW3W9X+ktlatWsXAwACDg4NhV6Wo8fFx2tvbZ/Xe9vZ2Vq0KftT2hMPYzE4G/gN4u3PuyRMtT0QqK+5P0uqaHsZ5V+HK3ThCpFZaWlpYu3Zt2NUoqb+/nxe+8IU12VbZMDazW4E+YLGZDQAfBVoAnHM3AR8BFgFf9Mf+067QeLiIhCLXM+4+rmfshfHweFphLBKysmHsnLuizPJ3AQUnbIlI+CaHqdun/rjnnmtGtUj4dKBIJOLi/lW2pg9T556PJnWusUjYFMYiEZfrGXe1TT1Fo6vVD+OJTM3rJCJTKYxFIi4+kaG1uYm25mlh7IdzQj1jkdApjEUiLj6ROm7yFhwbps71nEUkPApjkYiLj6dLhvGowlgkdApjkYiLT2QKhnFnizdMrWPGIuFTGItEXLFh6qYmo7M1pp6xSB1QGItEXHwiTXd74UsKdLU169QmkTqgMBaJuNEiw9TgXZVLw9Qi4VMYi0TcyHj6uAt+5GiYWqQ+KIxFIi4+kTruUpg5GqYWqQ8KY5EIS2eyjKeyRYepu1pjGqYWqQMKY5EIywVt0TBua9YwtUgdUBiLRNjIhHeTiJITuDRMLRI6hbFIhE32jIscM+5s1WxqkXqgMBaJsNx1pztbYwWXd7fFGE2mcc7VsloiMo3CWCTCxpJer7eztfgxY+dgLKXesUiYFMYiEZa7PWKxnnGn7twkUhcUxiIRluvxdpQYpgZI6LixSKgUxiIRdmyYunAYd7WqZyxSDxTGIhGWyIVxS/FjxqB7GouETWEsEmG5Yer21sI/6rkwzoW2iIRDYSwSYYlkmliT0RorEsb+8LWGqUXCpTAWibCxZJbOlhhmVnC5hqlF6oPCWCTCxlLpojOp4dgELg1Ti4RLYSwSYYlkpmQY544l66IfIuFSGItEWCKZoaOleBi3xpqINdnkxUFEJBwKY5EIG0tmip5jDGBmdLbEGEtma1grEZmubBib2c1mdsDMthZZbmb2eTPbbmYPm9n5la+miMzGWCpT9LrUOe2tMcZS6hmLhClIz/gW4OISyy8BTvO/rgX++cSrJSKVUO6YMXhX59IELpFwlf6TGXDO/crM1pRY5TLga867B9t9ZjbfzJY75/ZVqI4iMktjyXTJY8YAHS2xyctmnjDnIJsBl/G+r7KmTBJS41XfTi2oLfWnKZP0/h8XOTWwksqGcQArgd15zwf81xTGIiFLlDlmDF7POPBs6ok4PH4H7PoN7H8MRgdh7AikEl4IU9v7Ir8C4J6abrJq1Jb68wqAlz8L7fOqvq1KhHGhPxkK/kSa2bV4Q9ksW7aM/v7+CmzeE4/HK1pemNSW+tSIbRkZm+DQgX309x+efG16O8ZHxxgepnTbXJYVe3/Muh3foDkzSqq5m3j3OpKta0h1nUsm1o6zJqAJZ03+99XvTUwkk7S1tlZ9O7WgttSfiWSSwd/8jmys+m2pRBgPAKvznq8C9hZa0Tm3CdgEsHHjRtfX11eBzXv6+/upZHlhUlvqUyO2JfWTzZy+7hT6+s6YfG16O76x6wH2Hh2jr++iwoVk0vD9/wZPfQfW9UHfh2hZ/WIW1GDorpxG3CfFqC31p5btqMSpTbcDV/mzql8CDOl4sUj4kuks6awre8y47DD1L/4WHvkOvPpv4O3fh5NfUpNjaCJzSdmesZndCvQBi81sAPgo0ALgnLsJ2AxcCmwHEsA11aqsiASXm5TVUebUpo6WWPGLfuz+Pfzm/4cNV8Mr/meFaygiOUFmU19RZrkD/qJiNRKRikj45w6Xm8DVUezUJufgzhugZwW87lPVqKKI+HQFLpGIyvWMg8ymHi80TP3ML2HvH+CVH4C27mpUUUR8CmORiMr1doOcZ5zKOFKZaZfE/N2/QNdSWH9ltaooIj6FsUhE5SZllbsCV275lKHq0YPw1E/gvLdCc1vV6igiHoWxSEQlAg9Te1NHpgxVP/Z9yKbhvJJTRkSkQhTGIhE1OZu6pcxsav+exlN6xk/8GBaug6VnVa1+InKMwlgkosaCzqb2w3ry9KbkKDzzKzj9Yp1PLFIjCmORiJqcwBVgNjXkDVPv+i1kJuB5r61q/UTkGIWxSESNBQzj4yZw7f4dWBOsvqCq9RORYxTGIhE1eZ5xgFObID+M74Nl50BbT1XrJyLHKIxFIiqRytAaa6I5VvrHPDdMPZbMeDeFGHjQu/60iNSMwlgkosaSmbJD1HBsmHoslYH9WyE1CqtfXO3qiUgehbFIRCWS6bJX3wLonJxNnYGB+70XdbxYpKYqcT9jEalDiWSm7GlNkNczTqbh8CPQsQDmrS7zLhGpJPWMRSJqPBVsmLolZsSazBumPvAYLD1b5xeL1JjCWCSigvaMzYzOlhiJiRQc2AbLdNUtkVpTGItEVCKZoT3AMWPwhqrbR/dAMq5LYIqEQGEsElFjAXvG4IXxwvh278mys6tYKxEpRGEsElFjqczkHZnK6WiJsTjxtPdk6ZlVrJWIFKIwFomoRMDzjMG78Mey8Wdg3sm68pZICBTGIhE1lkyXvRRmTkdrjCWpPbD4eVWulYgUojAWiSDnHImApzaBdxvFkzJ7vXsYi0jNKYxFImgincW58ndsylkcG6XHxWHB2irXTEQKURiLRFDQOzblrMzu875Rz1gkFApjkQhKpPwwDjibeoXCWCRUCmORCBpLpgFoDzhMvSy9h6wzWLCmirUSkWIUxiIRlJjhMPXi5B72sZCktVazWiJShMJYJIImjxkH7BkvmBhgV3aZd7MIEak5hbFIBOWOGQedTd07todn3dLJEBeR2goUxmZ2sZk9YWbbzeyDBZbPM7MfmtlDZvaomV1T+aqKSFC5UA0UxukJOpKH2OMWk/CPNYtIbZUNYzOLAV8ALgHOAq4ws+m3dfkL4DHn3HlAH/BPZjr4JBKWY8eMA8ymHt4DwD4WaZhaJCRBesYXANudczucc0ngNuCyaes4oMfMDOgGDgP6E1skJLnZ1IF6xkMDAOxxizVMLRKSIGG8Etid93zAfy3fjcCZwF7gEeAvnXPZitRQRGZsLDWDCVx+GO91iyZ71CJSW0GuCGAFXnPTnr8e2AK8GjgV+KmZ3eOcG55SkNm1wLUAy5Yto7+/f8YVLiYej1e0vDCpLfWpkdry2PYkAL+79x6abOqP8PR2nLLzHtYCz7mF3P/Hh8juDXahkHrQSPukHLWl/tSyHUF+6gaA1XnPV+H1gPNdA3zaOeeA7Wb2DHAG8Pv8lZxzm4BNABs3bnR9fX2zrPbx+vv7qWR5YVJb6lMjteW3iW207dzJq1/1quOWHdeO279Lev8SJsZbOfX0M+l74fSBr/rVSPukHLWl/tSyHUGGqe8HTjOztf6krMuB26et8yzwGgAzWwY8H9hRyYqKSHCJZCbwOcYMDZDtXTn5PhGpvbI9Y+dc2syuB+4CYsDNzrlHzew6f/lNwCeAW8zsEbxh7RuccwerWG8RKWEslQl8XWqG9mALTwPQqU0iIQn00+qc2wxsnvbaTXnf7wVeV9mqichsjSUD3svYORgaoGndqyffJyK1pytwiURQIpmmI8h1qceOQGqU2PxVNDeZzjMWCYnCWCSCEkF7xiPPeY+9y+lojemYsUhIFMYiEeQdMw4QxnE/jLtPorM1pmFqkZAojEUiaCzobOpcz7jnJDpbmydvMCEitaUwFomgRDJDR5DrUueFcXtLbPIymiJSWwpjkQgaS2XoaA3w4x3fD6090NpFp44Zi4RGYSwSQYlkOth5xiPPQc8yAIWxSIgUxiIRk806xlPZYKc2xfdD90kAdLRoApdIWBTGIhEznp7BHZtG9kHPSZPr6zxjkXAojEUiJjfUXPY8Y+dgZP9kGHe0NmuYWiQkCmORiMkNNZcdpp4YhvQYdB87ZqzZ1CLhUBiLREyud1t2AtfIfu8xb5g6kcrg3QlVRGpJYSwSMbnjvmWPGY/s8x79nnF7SwznYCKdrWb1RKQAhbFIxORug1j2mHE81zNeDhwLbx03Fqk9hbFIxAQ+Zjx59a1jx4xB9zQWCYPCWCRijh0zDtAzbu6Atl7Am00NuqexSBgUxiIRkztmXHaYOnf1LTMAOv2etM41Fqk9hbFIxIwFnU0d3z85ectbX8eMRcKiMBaJmETQY8ajB6FryeTTXE9aw9QitacwFomYsWQaM2hvKfPjPToIXYsnn+Z60uoZi9SewlgkYkaTGTpbYph/LLigbAbGDk/tGbdoNrVIWBTGIhGTSGbobCtzvHjsCLhs4WFqTeASqTmFsUjEjCXT5U9rGh30HqcMU2sCl0hYFMYiETOazASYvJUL40LD1ApjkVpTGItEzFgyQ1e5YepcGHce6xk3NRntLU2Ma5hapOYUxiIRkwg0TH3Ie8zrGYM3o1oTuERqT2EsEjGJwMPUBp0Lp7zc0RLTMLVICBTGIhGTCDpM3bkImqaGdkdrTBf9EAlBoDA2s4vN7Akz225mHyyyTp+ZbTGzR83sl5WtpogElUhmyl+XenTwuCFq8GZUq2csUntl/nwGM4sBXwD+H2AAuN/MbnfOPZa3znzgi8DFzrlnzWxptSosIqUlkunJmz4UNXpwymlNOR0t6hmLhCFIz/gCYLtzbodzLgncBlw2bZ0rgf9wzj0L4Jw7UNlqikgQzjnGUgEu+pEoHMadrTESKU3gEqm1IGG8Etid93zAfy3f6cACM+s3swfN7KpKVVBEghtPZXEuwL2Miw5TN2uYWiQEZYepgUIXuHUFytkAvAboAH5rZvc5556cUpDZtcC1AMuWLaO/v3/GFS4mHo9XtLwwqS31qRHaMjzh/WgO7NxBv9tdcJ3R4aMwPsQzB+LsmtaeocMTHB3J1H07cxphnwSlttSfWrYjSBgPAKvznq8C9hZY56BzbhQYNbNfAecBU8LYObcJ2ASwceNG19fXN8tqH6+/v59KlhcmtaU+NUJbdh9OwN13s/6cM+nbsKrgOvfe9V0A1p7zItZu7Juy7BdDW9l6ZG/dtzOnEfZJUGpL/allO4IMU98PnGZma82sFbgcuH3aOj8ALjKzZjPrBF4MbKtsVUWknNwQc6lh6tbkkPdNZ4EJXJpNLRKKsj1j51zazK4H7gJiwM3OuUfN7Dp/+U3OuW1m9mPgYSAL/Ktzbms1Ky4ixxv1r55V6tSmlpQfxgWOGXe0xEims2SyjlhTiVswikhFBRmmxjm3Gdg87bWbpj3/B+AfKlc1EZmp3GlJXa3Ff7Qne8ZFzjMG7/SonvaWyldQRArSFbhEIiTIMPWxnnGhYWovxHWusUhtKYxFIiQRYJi6NXkUmlqgfd5xyzp1G0WRUCiMRSIkEWCYuiU17PWK7fhjwseGqRXGIrWkMBaJkFyIlp3AVWCIOv99Y7qnsUhNKYxFIiQx4Q1Tlz616WjByVve+3TMWCQMCmORCEmkMrTGmmiJFf/R9nrGxcL42GxqEakdhbFIhIwFuH1ia7J4GLe3aJhaJAwKY5EIGZ1Il75JRHKUWHYCOmqR2HoAABgMSURBVBcVXKwJXCLhUBiLREgilSkdxqMHvceyw9QKY5FaUhiLRMhYMjM5CaugsmGcm8ClY8YitaQwFomQ0Yl06WPGo4PeY5Ewbm1uoiVmjKpnLFJTCmORCBlLZegKFMaFzzMG6GprZnRCPWORWlIYi0SIN4GrxDB1IjdMXSKMW5uJK4xFakphLBIhZU9tGj1IpqkVWruKrtLVFiMxoWFqkVpSGItESCLAMHWydX7JMrramifviywitaEwFomQxERm8jaIBY0eJNVy/N2a8mmYWqT2FMYiEZHOZElmsmXOMx4sH8YaphapOYWxSEQk/EtYlrvoR7K1XBirZyxSawpjkYjI9WaLzqZ2DhLBhql1owiR2lIYi0REfCIFQHd7kTCeGIZMMlDPeFTD1CI1pTAWiYi4H6DdbUWGqf1LYaZaekuW090WI5nJkkxnK1o/ESlOYSwSEbmrZnW3tRRZIRfGpU9tyg1za6hapHYUxiIRMTLuhWdX0Z6xdynMZGu5nrEXxprEJVI7CmORiDjWMy5yzDiR6xmXP2bslafjxiK1ojAWiYh4uTD2e8blwrjT71nrKlwitaMwFomIXBh3FQ3jQ9DaQzbWWrKc7smescJYpFYUxiIRMTqRpiVmtDUX+bEeHSx5t6acrlaFsUitKYxFIiI+kaarrRkzK7xC0DDODVPrmLFIzQQKYzO72MyeMLPtZvbBEuu9yMwyZvanlauiiAQRn0hP9moLShyCriVly5mcwKVjxiI1UzaMzSwGfAG4BDgLuMLMziqy3meAuypdSREpb3QiTU+xq29B4J5xt2ZTi9RckJ7xBcB259wO51wSuA24rMB67wG+CxyoYP1EJKDcMHVB2azXM+4sH8ZtzU00mY4Zi9RSkDBeCezOez7gvzbJzFYCbwJuqlzVRGQm4hOZ4mE8fhSy6UDD1GamOzeJ1FiJMa1JhWaDuGnPPwfc4JzLFJ08ApjZtcC1AMuWLaO/vz9gNcuLx+MVLS9Makt9qve2DB5O0NrTVLCOnaMDXAA89uwB4l0nl21HCxme3jVAf/9gVepaKfW+T2ZCbak/tWxHkDAeAFbnPV8F7J22zkbgNj+IFwOXmlnaOff9/JWcc5uATQAbN250fX19s6z28fr7+6lkeWFSW+pTvbcle+/PWbtqMX195x2/cNe9cD+ctfEVHNjdVLYdCx7sp3dhL31951enshVS7/tkJtSW+lPLdgQJ4/uB08xsLbAHuBy4Mn8F59za3Pdmdgtwx/QgFpHqGi11zNi/+pY3TH2obFnd7S2MaJhapGbKhrFzLm1m1+PNko4BNzvnHjWz6/zlOk4sEjLnHPFkmp6iYexdl9qbwFU+jHvbmxkZT1WugiJSUpCeMc65zcDmaa8VDGHn3NUnXi0RmYlEMoNzpS6FmQvjRYHK62lvZt/QeIVqJyLl6ApcIhEwWva61IPQPg+aS1+XOqenrUU9Y5EaUhiLREDu+G7Ri34kDgY6rSmnp72Z4TEdMxapFYWxSARM9oyLXQ5zdGZh3NvRwlgqQyqTrUT1RKQMhbFIBJS/feLBwMeL4VgPOz6u3rFILSiMRSIgF5pFh6lHB2c4TN0CwIjCWKQmFMYiETDsh2avH6JTZDMwdjjQTSJycqE+rElcIjWhMBaJgOExLzTndRQI49GD4LLQvSxweQpjkdpSGItEwJAfxt2Fhqnj+73HGYRxr4apRWpKYSwSAcPjKXramok1FbhRi8JYpO4pjEUiYGgsRW+hIWo4FsY9Mx+m1oU/RGpDYSwSAcNj6fJh3LU0cHndk2GsnrFILSiMRSJgeDxFb7HTmkb2Q1svtHYGLq8l1kRHS2xyYpiIVJfCWCQChsdShWdSg9cznsHx4pzejmb1jEVqRGEsEgHDJY8ZH5hVGPe0tzAyoZ6xSC0ojEUiYHg8XfiCHwDx56A7+PHinJ529YxFakVhLNLg0pks8Yl0iWHqA9Bz0ozL7Wlvmbyyl4hUl8JYpMHleq+9HQUmcE3EIRmffc9YE7hEakJhLNLgclffKjhMPYsLfuT0trfocpgiNaIwFmlwucAsOEwdP+A9ziKM53e2cDSRwjl3ItUTkQAUxiINbngsN0xd2Z7xgs4W0lk3ea9kEakehbFIgxsqdcemEwrjVgCOJjRULVJtCmORBpcbpi44gSu+HywGnYtmXG4ujA+PJk+ofiJSnsJYpMENl5rANbLfm0ndNPMf9QVdXnlHEgpjkWpTGIs0uKGxFM1NRmdr7PiFI3uhd8Wsyp2vYWqRmlEYizS4I4kkC7paMStwL+Ph2YfxQg1Ti9SMwlikwR0eTU4G53GG9kDvylmV29vRghkc1TC1SNUpjEUa3JHR1OTx3SnGhyE5MuuecazJmNfRwhENU4tUncJYpMEdGp1gYVeBnvHwXu9xlj1j8IaqD6tnLFJ1gcLYzC42syfMbLuZfbDA8reZ2cP+171mdl7lqyoihRxJpCZPQ5pieI/3eAJh7F2FS2EsUm1lw9jMYsAXgEuAs4ArzOysaas9A7zSOfcC4BPApkpXVESOl8k6jiaSLCrZM57dMDV45xofGdUwtUi1BekZXwBsd87tcM4lgduAy/JXcM7d65w74j+9D1hV2WqKSCHDYymyDhYUDGO/Z9yzfNblz+9s1XnGIjUQJIxXArvzng/4rxXzTuDOE6mUiARzyD/tqPAx4z3QtRSai8y0DmBhV4vCWKQGClw/7zgFTl6k4G1czOxVeGF8YZHl1wLXAixbtoz+/v5gtQwgHo9XtLwwqS31qR7b8uSRDAC7t2+j/+hTU5adu3MrLU29/GFanWfSjiPPJRlPZbnr53fTFiv0qyBc9bhPZkttqT81bYdzruQX8FLgrrznHwI+VGC9FwBPA6eXK9M5x4YNG1wl3X333RUtL0xqS32qx7b8eOs+d8oNd7hHBo4ev/DGFzt365XHvTyTdnzzvl3ulBvucHuOJE6gltVTj/tkttSW+lONdgAPuAKZGGSY+n7gNDNba2atwOXA7fkrmNnJwH8Ab3fOPVmhvxNEpIzBkQkAFne3TV3gHBzdBfNPPqHyF3d7Q9wH4xMnVI6IlFZ2mNo5lzaz64G7gBhws3PuUTO7zl9+E/ARYBHwRf+SfGnn3MbqVVtEwAtjM1jUPe248OhBSCVgwZoTKn9JT9vkdkSkeoIcM8Y5txnYPO21m/K+fxfwrspWTUTKOTAywcLOVlpi0wa5juz0HuefckLlL+1tn9yOiFSPrsAl0sAGRyYme69THN3lPZ5gzzg3TK2esUh1KYxFGthgvEgYH3nGezzBY8ZtzTHmd7YojEWqTGEs0sAGh8eLhPEu6F4GrZ0nvI0l3W0cGBk/4XJEpDiFsUiDcs4xGJ9gaU/78QuP7DzhIeqcpb1t6hmLVJnCWKRBHU2kSGVc8WPGFQpjr2esMBapJoWxSIPa7w8dL+udFsaZFAwNnPBM6pxlve0cGJ4gmy144T0RqQCFsUiD2nNkDIAV8zumLjj8DLgsLDq1IttZMb+DZCY7eR1sEak8hbFIg9p71AvjVdPDePBx73HJ8yuynZV++Xv87YlI5SmMRRrUwNExWmNNx18K8+AT3uPi0yuynVzPe6/CWKRqFMYiDWrv0XGWz2+nqWna3ZQGn4B5J0NrV0W2s3KB3zM+ojAWqRaFsUiD2nMkwYp5HccvGHyiYkPUAL3tzXS3NWuYWqSKFMYiDWrP0bHJXuukbBYOPlXRMDYzVs7vUBiLVJHCWKQBjacyHBiZYNX0MB56FtJjFQ1jgNULO9h9OFHRMkXkGIWxSAN65uAozsG6Jd1TFxzY5j0uOaOi21u7uItnDo7qXGORKlEYizSgHYOjAKxbPG2S1t4tYE2w7JyKbm/dkm4m0lkNVYtUicJYpAE9czAOeD3WKfb+ERY/vyI3iMiXC/1nDo5WtFwR8SiMRRrQjsFRTuptp6uteeqCfVtgxfqKb2/tki5/u/GKly0iCmORhrR9MM66JdN6xUMDEN8Pyysfxku62+hpb2a7wlikKhTGIg0mlcny+HMjnL2id+qCXb/1Hk95acW3aWactbyXR/cOV7xsEVEYizSc7QfiJNNZzlk5b+qCZ++Ftt6KT97KOWflPLbtGyadyValfJG5TGEs0mC27hkC4OwV08J4529g9YuhKVaV7Z6zspfxVJYdmsQlUnEKY5EG8/DAEJ2tsamnNR3Z5d0gYt0rq7bdc/2e+EO7j1ZtGyJzlcJYpMHct+MQG9csnHqDiKd+4j2efknVtrtucTcLu1q5b8fhqm1DZK5SGIs0kAMj4zx1IM7LTl00dcHjP4KF62Dx86q27aYm46XrFnHv0wdxTlfiEqkkhbFIA7l3+yGAqWE8vBee+SWc81+qvv2XnrqIfUPjPD2o48YilaQwFmkgdzy8j5N626dO3tryLXBZOO+Kqm//NWcuxQw2P7Kv6tsSmUsUxiIN4mgiyS+fPMAbz1tOLHe8ODUOv98E6/pg0alVr8PyeR1csGYh39+yR0PVIhWkMBZpEF//7S5SGcebz1917MX7v+RddevC99WsHn+6YRU7Bkf55ZODNdumSNQpjEUawOHRJF/+zTO85oylnLncv/LW4R3Q/2k47XWw9hU1q8tl61eyYl47//unT5LRLRVFKiJQGJvZxWb2hJltN7MPFlhuZvZ5f/nDZnZ+5asqMjdlso4P/PvDjE6k+Z8XP997cWQ/3HoFxFrg0n8Es9KFVFBrcxM3XHIGDw0M8X9+8VTNtisSZWXD2MxiwBeAS4CzgCvM7Kxpq10CnOZ/XQv8c4XrKTIn7R8e57pvPMjPtu3nw5eeyRlLOmHrd2FTHxx9Fv7s67DglJrX60/OW8GbX7iSz/3sKf5u8zZGxlM1r4NIlDSXX4ULgO3OuR0AZnYbcBnwWN46lwFfc96MjvvMbL6ZLXfO1WTK5eEDexjafi9/TA7M+L2lJqEYDii+vOR7nSvxTig192Vo4FkeTDxd7J0lquTy/i1UJyhaq7KfQ6lVir93eO8eHhh+rOB7T6RcSrzX8t5XqAQr9eG74u8dfm4f9x/6Y/GalZzQ5IquZtNKS6azDI2n2Dc0zq5Do5xCitvPiPGC3d+CX98LY4dhyZlwxa1VuV1iEGbGZ/70BbS3xtj0qx18875dvPx5izl9WQ+Lu1vpbm+hJWbEmozmJsPMqEbffev+NMlHn6tCyV4ba2nrgTSpx/bXdJvVEpW2bD2Q5uWZLC2x6h/RDRLGK4Hdec8HgBcHWGclMCWMzexavJ4zy5Yto7+/f4bVLWxo1xYuG/iMt9UI2ACwN+xaVFB1fleGo9ZzlvzLTGd3NZNoX8rQvPUcWvdiDi5+ETx5FJ7sn3GR8Xi8Yj97r1sAp7+0nV/sTvPQzgP8bNt+an4Y+Y8P1niDVfSHB8KuQeVEpC3PX/BLOluq/4dZkDAuVIvpP25B1sE5twnYBLBx40bX19cXYPPlxYfP44ebuzjzzOmj53kVLPVXrpX6q8dKHo4rtcyZUbw/ULzcrY9s5Zxzzy36TleqlzFZaPGKFa2zv6DwzvSXFX9zwfc99NBDnLd+fYntlqtvbrvF61tOqe0WbE+Rcv/w4B84f8OG0qvlvXj8J2KFVjtuWUvMmNfRcqxusVaaOhbQaUYnsLxg7YLr7++nUj97OVf7j9ms43AiSWIiQzqbJZN1pLOuahO9HnzwATZs2FiVsmvtgQceYONGtaWePPDAA7z+NX3HTiWsoiBhPACsznu+iuP7bUHWqZru3gX0nHQazzvv5bXaZFXt3B9nzZkbyq/YALYPHGTVqWeHXY2K6Hx6gJNWV/9c3kbW1GQs7m6D7tps7+BTseNvJdmg1Jb6c/CpWE2CGILNpr4fOM3M1ppZK3A5cPu0dW4HrvJnVb8EGKrV8WIREZFGV7Zn7JxLm9n1wF14R7Buds49ambX+ctvAjYDlwLbgQRwTfWqLCIiEi1Bhqlxzm3GC9z8127K+94Bf1HZqomIiMwNugKXiIhIyBTGIiIiIVMYi4iIhExhLCIiEjKFsYiISMgUxiIiIiFTGIuIiITMSt15qKobNhsEdlWwyMXAwQqWFya1pT5FpS1RaQeoLfUqKm2pRjtOcc4tmf5iaGFcaWb2gHOu8a9MjtpSr6LSlqi0A9SWehWVttSyHRqmFhERCZnCWEREJGRRCuNNYVeggtSW+hSVtkSlHaC21KuotKVm7YjMMWMREZFGFaWesYiISENqqDA2s7eY2aNmljWzjdOWfcjMtpvZE2b2+iLvX2hmPzWzp/zHBbWpeWlm9m0z2+J/7TSzLUXW22lmj/jrPVDregZhZh8zsz157bm0yHoX+/tqu5l9sNb1DMLM/sHMHjezh83se2Y2v8h6dblfyn3G5vm8v/xhMzs/jHqWY2arzexuM9vm//z/ZYF1+sxsKO//3UfCqGsQ5f6/NMJ+MbPn533WW8xs2MzeO22dut0nZnazmR0ws615rwXKh6r97nLONcwXcCbwfKAf2Jj3+lnAQ0AbsBZ4GogVeP/fAx/0v/8g8Jmw21Sgjv8EfKTIsp3A4rDrWKb+HwPeX2admL+P1gGt/r47K+y6F6jn64Bm//vPFPv/Uo/7JchnDFwK3AkY8BLgd2HXu0hblgPn+9/3AE8WaEsfcEfYdQ3YnpL/Xxplv+TVNwY8h3f+bEPsE+AVwPnA1rzXyuZDNX93NVTP2Dm3zTn3RIFFlwG3OecmnHPPANuBC4qs91X/+68C/7k6NZ0dMzPgz4Bbw65LlV0AbHfO7XDOJYHb8PZNXXHO/cQ5l/af3gesCrM+MxTkM74M+Jrz3AfMN7Plta5oOc65fc65P/jfjwDbgJXh1qqqGmK/5HkN8LRzrpIXcaoq59yvgMPTXg6SD1X73dVQYVzCSmB33vMBCv+wLnPO7QPvBxxYWoO6zcRFwH7n3FNFljvgJ2b2oJldW8N6zdT1/vDazUWGeoLur3ry53i9lULqcb8E+Ywbbj+Y2RrghcDvCix+qZk9ZGZ3mtnZNa3YzJT7/9Jo++VyincgGmWfQLB8qNq+aa5EIZVkZj8DTiqw6K+dcz8o9rYCr9XVNPGA7bqC0r3ilzvn9prZUuCnZva4/xdeTZVqC/DPwCfwPv9P4A27//n0Igq8N5T9FWS/mNlfA2ngm0WKqYv9Mk2Qz7hu9kMQZtYNfBd4r3NueNriP+ANk8b9eQrfB06rdR0DKvf/pWH2i5m1An8CfKjA4kbaJ0FVbd/UXRg75147i7cNAKvznq8C9hZYb7+ZLXfO7fOHfQ7Mpo6zUa5dZtYMvBnYUKKMvf7jATP7Ht6QSc1/6QfdR2b2JeCOAouC7q+qC7Bf3gG8AXiN8w8aFSijLvbLNEE+47rZD+WYWQteEH/TOfcf05fnh7NzbrOZfdHMFjvn6u76yAH+vzTMfgEuAf7gnNs/fUEj7RNfkHyo2r6JyjD17cDlZtZmZmvx/vr6fZH13uF//w6gWE87DK8FHnfODRRaaGZdZtaT+x5vctHWQuuGadqxrTdRuI73A6eZ2Vr/L+vL8fZNXTGzi4EbgD9xziWKrFOv+yXIZ3w7cJU/e/clwFBumK6e+HMpvgxsc859tsg6J/nrYWYX4P1uO1S7WgYT8P9LQ+wXX9HRvEbZJ3mC5EP1fneFOaNtpl94v9wHgAlgP3BX3rK/xpvl9gRwSd7r/4o/8xpYBPwceMp/XBh2m/LqeQtw3bTXVgCb/e/X4c3cewh4FG8YNfR6F2jH14FHgIf9/6TLp7fFf34p3qzYp+u4Ldvxjg9t8b9uaqT9UugzBq7L/T/DG3L7gr/8EfLOUKinL+BCvKHAh/P2xaXT2nK9//k/hDfZ7mVh17tIWwr+f2nQ/dKJF67z8l5riH2C9wfEPiDlZ8o7i+VDrX536QpcIiIiIYvKMLWIiEjDUhiLiIiETGEsIiISMoWxiIhIyBTGIiIiIVMYi4iIhExhLCIiEjKFscgcYWYv8m/g0e5fCepRMzsn7HqJCLroh8hcYmafBNqBDmDAOfe/Qq6SiKAwFplT/Ovp3g+M412eMBNylUQEDVOLzDULgW6gB6+HLCJ1QD1jkTnEzG4HbgPW4t3E4/qQqyQi1OH9jEWkOszsKiDtnPuWmcWAe83s1c65X4RdN5G5Tj1jERGRkOmYsYiISMgUxiIiIiFTGIuIiIRMYSwiIhIyhbGIiEjIFMYiIiIhUxiLiIiETGEsIiISsv8LPejPF2GBUokAAAAASUVORK5CYII=\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "plot_normal_pdf_cdf(mu=0, sigma=0.3)\n", "# notice how the pdf can be larger than 1 " ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "## Bayes Rule\n", "---\n", "Suppose that the events $B_1, ..., B_k$ are mutually exclusive and form a partiotion of the sample space (i.e. one of them must occur), then for any event $Pr(A)$, **Bayes Rule**:\n", "$$Pr(B_i|A) = \\frac{Pr(A,B_i)}{Pr(A)} = \\frac{Pr(A|B_i)Pr(B_i)}{Pr(A)} = \\frac{Pr(A|B_i)Pr(B_i)}{\\sum_{j=1}^k Pr(A|B_j)Pr(B_j)}$$" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "subslide" } }, "source": [ "#### REMEMBER THESE!\n", "* **Posterior Distribution** - $Pr(B_i|A)$\n", "* **Liklihood Distribution** - $Pr(A|B_i)$\n", "* **Prior Distribution** - $Pr(B_i)$\n", "* **Evidence** - $Pr(A)$" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "subslide" } }, "source": [ "\n", "\n", "* Image Source" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "subslide" } }, "source": [ "### Example \n", "---\n", "Given a dataset where each sample is a male or a female with their height, what is the probability that given a certain height, that person is a female, that is, calculate: $Pr(\\textit{Gender} = \\textit{Female} | \\textit{Height} = X cm)$?" ] }, { "cell_type": "code", "execution_count": 6, "metadata": { "slideshow": { "slide_type": "subslide" } }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Number of rows in the dataset: 10000\n" ] }, { "data": { "text/html": [ "
\n", "\n", "\n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", "
GenderHeight
0Male187.571423
1Male174.706036
2Male188.239668
3Male182.196685
4Male177.499761
5Male170.822660
6Male174.714106
7Male173.605229
8Male170.228132
9Male161.179495
\n", "
" ], "text/plain": [ " Gender Height\n", "0 Male 187.571423\n", "1 Male 174.706036\n", "2 Male 188.239668\n", "3 Male 182.196685\n", "4 Male 177.499761\n", "5 Male 170.822660\n", "6 Male 174.714106\n", "7 Male 173.605229\n", "8 Male 170.228132\n", "9 Male 161.179495" ] }, "execution_count": 6, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# load the data\n", "dataset = pd.read_csv('./datasets/heights_dataset.csv')\n", "# use only the heights\n", "dataset = dataset.drop('Weight', axis=1)\n", "# inch -> cm\n", "dataset['Height'] = dataset['Height'] * 2.54\n", "## print the number of rows in the data set\n", "number_of_rows = len(dataset)\n", "print('Number of rows in the dataset: {}'.format(number_of_rows))\n", "## show the first 10 rows\n", "dataset.head(10)" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "### Histogram\n", "---\n", "* A histogram is an accurate representation of the distribution of numerical data. \n", "* It is an estimate of the probability distribution (PDF) of a continuous variable. \n", "* To construct a histogram, the first step is to \"bin\" (or \"bucket\") the range of values—that is, divide the entire range of values into a series of intervals—and then count how many values fall into each interval. \n", " * The bins are usually specified as consecutive, non-overlapping intervals of a variable. \n", " * The bins (intervals) must be adjacent, and are often (but are not required to be) of equal size." ] }, { "cell_type": "code", "execution_count": 7, "metadata": { "slideshow": { "slide_type": "subslide" } }, "outputs": [ { "data": { "text/plain": [ "Text(0.5, 0, 'Height(cm)')" ] }, "execution_count": 7, "metadata": {}, "output_type": "execute_result" }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYsAAAEGCAYAAACUzrmNAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8QZhcZAAAgAElEQVR4nO3df5iVdZ3/8ecLVAYFETRnEdwd8iJKg1AHta9tDGph1jK2leHaBmsuuy2aFlnYbol24eW2amtYFkVJ/ppMSzFrV2Q92XZJ4A9kVERoZW2S0FhFhgAdfH//OPfgYTgz9znD3DPnMK/Hdc019/25P/f9eX9mzjnvc//63IoIzMzMujKgrwMwM7PK52RhZmapnCzMzCyVk4WZmaVysjAzs1QH9HUAWTniiCOirq4u0za2bdvGIYcckmkbWav2PlR7/FD9faj2+KH6+9CT8T/66KN/jIi3dCzfb5NFXV0djzzySKZt5HI5GhoaMm0ja9Xeh2qPH6q/D9UeP1R/H3oyfkn/W6zch6HMzCyVk4WZmaVysjAzs1T77TkLM+vfXn/9dVpaWtixY0dq3WHDhrFmzZpeiCob3Ym/pqaG0aNHc+CBB5ZU38nCzPZLLS0tDB06lLq6OiR1WXfr1q0MHTq0lyLreeXGHxFs3ryZlpYWxowZU9I6PgxlZvulHTt2cPjhh6cmiv5IEocffnhJe13tnCzMbL/lRNG5cv82ThZmZpbK5yzMrF+om3tfj25vw9UfTK0jiU984hPcfPPNALS1tTFy5EhOPvlkfvazn3W6Xi6X45prrumyTm9zsjDbT/X0h2MxpXxg9meHHHIITz75JNu3b2fw4MEsXbqUUaNG9XVY3eJkYZVj3rDy1xl3BcxrLHH7W8rfvtk++sAHPsB9993HRz/6UW6//XbOPfdcfvWrXwGwYsUKLrnkkt3J5Ac/+AHjxo3bY/1t27Zx0UUX0dzcTFtbG/PmzaOxscTXfA/yOQszswxNnz6dpqYmduzYwerVqzn55JN3L3v729/OQw89xOOPP86VV17Jl770pb3Wnz9/PqeddhorV67kwQcf5NJLL2Xbtm292QXAexZmZpmaMGECGzZs4Pbbb+ess87aY9mWLVuYMWMG69atQxKvv/76Xuvff//9LFmyhGuuuQbIXxL8/PPP8453vKNX4m/nZGFmlrFp06bx+c9/nlwux+bNm3eXf/nLX2bKlCn89Kc/ZcOGDUVHjo0I7rrrrr0OT/U2H4YyM8vY+eefz1e+8hXGjx+/R/mWLVt2n/C+6aabiq47depUFixYQEQA8Pjjj2caa2e8Z2HWyzpepTRnfBsze+HKpf6uqyu3sh7uY/To0Vx88cV7lX/hC19gxowZXHfddZx22mlF1/3yl7/MJZdcwoQJE4gI6urq+uSSWicLM7OMtLa27lXW0NCw+3DTu9/9bp599tndy7761a/uVWfw4MF85zvfyTzWND4MZWZmqZwszMwslQ9DWf/RnZv+St62b/iz/Zv3LMzMLJWThZmZpXKyMDOzVD5nYWb9QxfnrLp1h0UJ56kGDhy4x414d999N3V1dd1pLdVNN93EI488wg033JDJ9p0szMwyMnjwYFatWtXXYfQIH4YyM+tFu3bt4tJLL2XSpElMmDBh9w13uVyOyZMnc8455/C2t72NuXPncuutt3LSSScxfvx4fvvb3wJw7733cvLJJ3P88cdzxhlnsGnTpr3aeOmll/jIRz7CpEmTmDRpEr/+9a/3Oe7MkoWk70t6UdKTBWX/JukZSasl/VTSYQXLLpO0XtJaSVMLyk+U1Jws+4b8UF0zqxLbt29n4sSJTJw4kQ9/+MMALFq0iGHDhrFy5UpWrlzJd7/7XZ577jkAnnjiCa6//nqam5u5+eabefbZZ1mxYgUXXHABCxYsAOA973kPy5cv5/HHH2f69Ol87Wtf26vdiy++mM9+9rOsXLmSu+66iwsuuGCf+5LlYaibgBuAHxaULQUui4g2Sf8KXAZ8UdKxwHTgOOAo4AFJb4uIXcCNwCxgOfBz4EzgFxnGbWbWI4odhrr//vtZvXo1d955J5AfTHDdunUcdNBBTJo0iZEjRwJwzDHH8P73vx+A8ePH8+CDDwLQ0tLCxz/+cTZu3Mhrr73GmDFj9mr3gQce4Omnn949/+qrr+7z+FeZJYuIeEhSXYey+wtmlwMfTaYbgaaI2Ak8J2k9cJKkDcChEfEwgKQfAmfjZGFmVSoiWLBgAVOnTt2jPJfLMWjQoN3zAwYM2D0/YMAA2traALjooov43Oc+x7Rp08jlcsybN2+vNt544w0efvhhBg8e3GNx9+UJ7vOBHyXTo8gnj3YtSdnryXTH8qIkzSK/F0JtbS25XK4Hw91ba2tr5m1kraL6MO6KsldpHXQUuW6s1+PK+BvOGd+2x3zt4L3LqkUul6us11CBYcOGsXXr1t3zPT2mbOG2y6k3efJkFixYwKRJkzjwwANZt24dRx11FH/6059oa2vbXX/Xrl1s27aNrVu37rHs5Zdf5rDDDmPr1q1873vfY9euXezatYsdO3bw2muvsXXrVqZMmcK11167e6Tb1atXM2HChL1i27FjR8n/uz5JFpL+GWgDbm0vKlItuigvKiIWAgsB6uvro9iDRHpSLpcr+rCSalJRfSj1WdoFcuOuoGHt5RkEU6ZzSx/uo+Nw5HPGt3Ftc3VemLjhvIbKeg0VWLNmzZ6HXbq41LU7h2hKrd1xuxdeeCF/+MMfmDx5MhHBW97yFu6++24OPvhgDjjggN31Bw4cyCGHHMLQoUP3WHbllVcyc+ZMRo0axSmnnEJLSwsDBw6kpqaGgw46iKFDh3LjjTcye/ZsTj31VNra2njve9/Lt7/97b1iq6mp4fjjjy+pH73+CpU0A/gQcHq0P80jv8dwdEG10cALSfnoIuVmZhWv2BDlAwYM4KqrruKqq67ao7xwWHJgj2/8hcsaGxtpbNzzi9XWrVuZOXMmM2fOBOCII47gRz/6ET2pVy+dlXQm8EVgWkT8qWDREmC6pEGSxgBjgRURsRHYKumU5CqoTwL39GbMZmaW4Z6FpNuBBuAISS3A5eSvfhoELE2ugF0eEf8YEU9JugN4mvzhqdnJlVAAnyZ/ZdVg8ie2fXLbzKyXZXk11LlFihd1UX8+ML9I+SPAO3swNDPrJyIC35pV3JtnAUrjO7jNbL9UU1PD5s2by/5Q7A8igs2bN1NTU1PyOtV5CYaZWYrRo0fT0tLCSy+9lFp3x44dZX1wVpruxF9TU8Po0aPTKyacLMxsv3TggQcWvbu5mFwuV/IlpJWoN+L3YSgzM0vlZGFmZqmcLMzMLJWThZmZpXKyMDOzVE4WZmaWysnCzMxSOVmYmVkqJwszM0vlZGFmZqmcLMzMLJWThZmZpXKyMDOzVE4WZmaWysnCzMxSOVmYmVkqJwszM0vlZGFmZqmcLMzMLFVmyULS9yW9KOnJgrIRkpZKWpf8Hl6w7DJJ6yWtlTS1oPxESc3Jsm9IUlYxm5lZcVnuWdwEnNmhbC6wLCLGAsuSeSQdC0wHjkvW+Zakgck6NwKzgLHJT8dtmplZxjJLFhHxEPB/HYobgcXJ9GLg7ILypojYGRHPAeuBkySNBA6NiIcjIoAfFqxjZma9RPnP4Iw2LtUBP4uIdybzr0TEYQXLX46I4ZJuAJZHxC1J+SLgF8AG4OqIOCMp/0vgixHxoU7am0V+L4Ta2toTm5qasuoaAK2trQwZMiTTNrJWUX3YuKrsVVoHHcWQnS9kEEyZRk4suWrz77fsMV87GDZt7+mAesf4UcMq6zXUTdXeh56Mf8qUKY9GRH3H8gN6ZOv7rth5iOiivKiIWAgsBKivr4+GhoYeCa4zuVyOrNvIWkX1YV5j2avkxl1Bw9rLMwimTOduSa+TmDn3vj3m54xv49rmSnkrlmfDeQ2V9RrqpmrvQ2/E39tXQ21KDi2R/H4xKW8Bji6oNxp4ISkfXaTczMx6UW8niyXAjGR6BnBPQfl0SYMkjSF/IntFRGwEtko6JbkK6pMF65iZWS/JbN9X0u1AA3CEpBbgcuBq4A5JnwKeBz4GEBFPSboDeBpoA2ZHxK5kU58mf2XVYPLnMX6RVcxmZlZcZskiIs7tZNHpndSfD8wvUv4I8M4eDM3MzMrkO7jNzCxVdV6CYWYVoW7ufcwZ37bXFV49acPVH8xs21Y671mYmVkqJwszM0vlZGFmZqmcLMzMLJWThZmZpXKyMDOzVE4WZmaWysnCzMxSOVmYmVkqJwszM0vlZGFmZqmcLMzMLJWThZmZpXKyMDOzVE4WZmaWqqTnWUh6Z0Q8mXUwVgXmDevrCMysD5S6Z/FtSSsk/ZOkwzKNyMzMKk5JySIi3gOcBxwNPCLpNknvyzQyMzOrGCWfs4iIdcC/AF8EJgPfkPSMpL/OKjgzM6sMJSULSRMkfR1YA5wG/FVEvCOZ/nqG8ZmZWQUodc/iBuAx4F0RMTsiHgOIiBfI722URdJnJT0l6UlJt0uqkTRC0lJJ65LfwwvqXyZpvaS1kqaW256Zme2bUpPFWcBtEbEdQNIASQcDRMTN5TQoaRTwGaA+It4JDASmA3OBZRExFliWzCPp2GT5ccCZwLckDSynTTMz2zelJosHgMEF8wcnZd11ADBY0gHJtl4AGoHFyfLFwNnJdCPQFBE7I+I5YD1w0j60bWZmZVJEpFeSVkXExLSykhuVLgbmA9uB+yPiPEmvRMRhBXVejojhkm4AlkfELUn5IuAXEXFnke3OAmYB1NbWntjU1NSd8ErW2trKkCFDMm0ja2X3YeOq7ILphtZBRzFk5wt9HQaMLP2t0Pz7LXvM1w6GTdt7OqDek3X840dlf29Ptb+XezL+KVOmPBoR9R3LS7opD9gm6YT2cxWSTiT/QV+25FxEIzAGeAX4saRPdLVKkbKiGS4iFgILAerr66OhoaE7IZYsl8uRdRtZK7sP8xozi6U7cuOuoGHt5X0dBpy7Jb1OYubc+/aYnzO+jWubS30rVp6s499wXkNm225X7e/l3oi/1P/wJeQ/1Nu/wo0EPt7NNs8AnouIlwAk/QT4f8AmSSMjYqOkkcCLSf0W8vd3tBtN/rCVmZn1kpKSRUSslPR2YBz5b/rPRMTr3WzzeeCU5AT5duB04BFgGzADuDr5fU9Sfwlwm6TrgKOAscCKbrZtlo0yhkHZULPnfG7AFWyo6XrvqG7Hbd2JyqzHlLPvOAmoS9Y5XhIR8cNyG4yI30i6k/yluG3A4+QPHQ0B7pD0KfIJ5WNJ/ack3QE8ndSfHRG7ym3XzMy6r9SBBG8GjgFWAe0f1AGUnSwAIuJyoONXqZ3k9zKK1Z9P/oS4mZn1gVL3LOqBY6OUS6fMzGy/U+p9Fk8Cf5ZlIGZmVrlK3bM4Anha0gryh4sAiIhpmURlZmYVpdRkMS/LIMzMrLKVeunsLyX9BTA2Ih5ILnv1+ExmZv1EqUOU/z1wJ/CdpGgUcHdWQZmZWWUp9QT3bOBU4FXY/SCkI7MKyszMKkupyWJnRLzWPpOMFuvLaM3M+olSk8UvJX2J/LDi7wN+DNybXVhmZlZJSk0Wc4GXgGbgH4Cf040n5JmZWXUq9WqoN4DvJj9mZtbPlDo21HMUOUcREW/t8YjMzKzilDM2VLsa8iPCjuj5cMzMrBKVdM4iIjYX/Pw+Iv4dOC3j2MzMrEKUehjqhILZAeT3NIZmEpGZmVWcUg9DXVsw3QZsAM7p8WjMzKwilXo11JSsAzEzs8pV6mGoz3W1PCKu65lwzMysEpVzNdQkYEky/1fAQ8DvsgjKzMwqSzkPPzohIrYCSJoH/DgiLsgqMDMzqxylDvfx58BrBfOvAXU9Ho2ZmVWkUvcsbgZWSPop+Tu5Pwz8MLOozMysopR6U9584O+Al4FXgL+LiKu626ikwyTdKekZSWskvVvSCElLJa1Lfg8vqH+ZpPWS1kqa2t12zcyse0o9DAVwMPBqRFwPtEgasw/tXg/8R0S8HXgXsIb8yLbLImIssCyZR9KxwHTgOOBM4FuS/EhXM7NeVOpjVS8HvghclhQdCNzSnQYlHQq8F1gEEBGvRcQrQCOwOKm2GDg7mW4EmiJiZ0Q8B6wHTupO22Zm1j2l7ll8GJgGbAOIiBfo/nAfbyX/bIwfSHpc0vckHQLURsTGZPsbefOxraPY8xLdlqTMzMx6iSLSn44qaUVEnCTpsYg4IflwfzgiJpTdoFQPLAdOjYjfSLqe/LO9L4qIwwrqvRwRwyV9M2nrlqR8EfDziLiryLZnAbMAamtrT2xqaio3vLK0trYyZMiQTNvIWtl92Lgqu2C6oXXQUQzZ+UJfh7FPSulD8xv7ctQ3W7WDYdP27LY/ftSw7DaeqPb3ck/GP2XKlEcjor5jealXQ90h6TvAYZL+Hjif7j8IqQVoiYjfJPN3kj8/sUnSyIjYKGkk8GJB/aML1h8NFH1nRcRCYCFAfX19NDQ0dDPE0uRyObJuI2tl92FeY2axdEdu3BU0rL28r8PYJ6X0YeaO23opmvLNGd/Gtc2lfpSUb8N5DZltu121v5d7I/5Sr4a6hvyH+l3AOOArEbGgOw1GxB+A30kalxSdDjxN/u7wGUnZDOCeZHoJMF3SoOSk+lhgRXfaNjOz7kn9OpBcefSfEXEGsLSH2r0IuFXSQcD/kL8sdwD5PZhPAc+Tf8ASEfGUpDvIJ5Q2YHZE7OqhOMzMrASpySIidkn6k6RhEbGlJxqNiFXs+fS9dqd3Un8+ML8n2jYzs/KVeqBxB9AsaSnJFVEAEfGZTKIyM7OKUmqyuC/5MTOzfqjLZCHpzyPi+YhY3FU9MzPbv6VdDXV3+4Skve5rMDOz/iEtWahg+q1ZBmJmZpUrLVlEJ9NmZtaPpJ3gfpekV8nvYQxOpknmIyIOzTQ6MzOrCF0mi4jwUOBmZlbW8yzMzKyfcrIwM7NUThZmZpbKycLMzFI5WZiZWSonCzMzS+VkYWZmqZwszMwslZOFmZmlcrIwM7NUThZmZpbKycLMzFI5WZiZWSonCzMzS9VnyULSQEmPS/pZMj9C0lJJ65LfwwvqXiZpvaS1kqb2VcxmZv1VX+5ZXAysKZifCyyLiLHAsmQeSccC04HjgDOBb0nyczbMzHpRnyQLSaOBDwLfKyhuBBYn04uBswvKmyJiZ0Q8B6wHTuqtWM3MrO/2LP4d+ALwRkFZbURsBEh+H5mUjwJ+V1CvJSkzM7Neoojo3QalDwFnRcQ/SWoAPh8RH5L0SkQcVlDv5YgYLumbwMMRcUtSvgj4eUTcVWTbs4BZALW1tSc2NTVl2pfW1laGDBmSaRtZK7sPG1dlF0w3tA46iiE7X+jrMPZJKX1ofmNML0VTvtrBsGl7dtsfP2pYdhtPVPt7uSfjnzJlyqMRUd+xvMtncGfkVGCapLOAGuBQSbcAmySNjIiNkkYCLyb1W4CjC9YfDRR9Z0XEQmAhQH19fTQ0NGTUhbxcLkfWbWSt7D7Ma8wslu7IjbuChrWX93UY+6SUPszccVsvRVO+OePbuLY5u4+SDec1ZLbtdtX+Xu6N+Hv9MFREXBYRoyOijvyJ6/+KiE8AS4AZSbUZwD3J9BJguqRBksYAY4EVvRy2mVm/1hd7Fp25GrhD0qeA54GPAUTEU5LuAJ4G2oDZEbGr78I0M+t/+jRZREQOyCXTm4HTO6k3H5jfa4GZmdkefAe3mZmlcrIwM7NUThZmZpaqkk5wm5ntpW7ufZm3MWd8Gw2Zt1LdvGdhZmapnCzMzCyVk4WZmaVysjAzs1Q+wb2/mVfmoGvjrqi48Z7MrPJ4z8LMzFI5WZiZWSonCzMzS+VzFmZVYEPN32S27boKflaGVQ7vWZiZWSonCzMzS+VkYWZmqZwszMwslZOFmZmlcrIwM7NUThZmZpbKycLMzFI5WZiZWSonCzMzS9XryULS0ZIelLRG0lOSLk7KR0haKmld8nt4wTqXSVovaa2kqb0ds5lZf9cXexZtwJyIeAdwCjBb0rHAXGBZRIwFliXzJMumA8cBZwLfkjSwD+I2M+u3ej1ZRMTGiHgsmd4KrAFGAY3A4qTaYuDsZLoRaIqInRHxHLAeOKl3ozYz698UEX3XuFQHPAS8E3g+Ig4rWPZyRAyXdAOwPCJuScoXAb+IiDuLbG8WMAugtrb2xKampkzjb21tZciQIZm2UbaNq8qq3jroKIbsfCGjYLJX7fFD3/eh+Y0x+7R+7WDYtL2HgukjtYPhyBFlPmWygvTkZ9GUKVMejYj6juV9NkS5pCHAXcAlEfGqpE6rFikrmuEiYiGwEKC+vj4aGhp6INLO5XI5sm6jbGU+IjU37goa1l6eUTDZq/b4oe/7MHMfhyifM76Na5ur+2kHc8a3cU6lvZfL0BufRX1yNZSkA8knilsj4idJ8SZJI5PlI4EXk/IW4OiC1UcD1f1V0sysyvTF1VACFgFrIuK6gkVLgBnJ9AzgnoLy6ZIGSRoDjAVW9Fa8ZmbWN4ehTgX+FmiW1H6A/UvA1cAdkj4FPA98DCAinpJ0B/A0+SupZkfErt4P28ys/+r1ZBER/03x8xAAp3eyznxgfmZBmZlZl3wHt5mZpXKyMDOzVE4WZmaWysnCzMxSOVmYmVkqJwszM0vlZGFmZqmcLMzMLJWThZmZpXKyMDOzVE4WZmaWysnCzMxSOVmYmVmq6n68VTWaV72PbjTbn9XNvS/zNjZc/cHM28iK9yzMzCyV9yzM+rkNNX+zT+vnBlzBhprOnyFet4/P+LbK4D0LMzNL5WRhZmapnCzMzCyVk4WZmaVysjAzs1ROFmZmlqpqLp2VdCZwPTAQ+F5EXJ1ZY6XeODfuCpjXmFkYZmaVoiqShaSBwDeB9wEtwEpJSyLi6b6NzMysdFndJT5nfBszk21ndZd4VSQL4CRgfUT8D4CkJqARcLIwq3D7etNfV3zDX+9RRPR1DKkkfRQ4MyIuSOb/Fjg5Ii7sUG8WMCuZHQeszTi0I4A/ZtxG1qq9D9UeP1R/H6o9fqj+PvRk/H8REW/pWFgtexYqUrZXlouIhcDC7MPJk/RIRNT3VntZqPY+VHv8UP19qPb4ofr70BvxV8vVUC3A0QXzo4EX+igWM7N+p1qSxUpgrKQxkg4CpgNL+jgmM7N+oyoOQ0VEm6QLgf8kf+ns9yPiqT4OC3rxkFeGqr0P1R4/VH8fqj1+qP4+ZB5/VZzgNjOzvlUth6HMzKwPOVmYmVkqJ4suSPq+pBclPVlk2eclhaQjCsouk7Re0lpJU3s32r0Vi1/SPEm/l7Qq+TmrYFlFxQ+d/w8kXZTE+ZSkrxWUV1QfOvkf/Kjg779B0qqCZRUVP3Tah4mSlid9eETSSQXLKqoPncT/LkkPS2qWdK+kQwuWVVr8R0t6UNKa5PV+cVI+QtJSSeuS38ML1un5PkSEfzr5Ad4LnAA82aH8aPIn2/8XOCIpOxZ4AhgEjAF+CwystPiBecDni9StuPi76MMU4AFgUDJ/ZKX2obPXUMHya4GvVGr8XfwP7gc+kEyfBeQqtQ+dxL8SmJxMnw98tYLjHwmckEwPBZ5N4vwaMDcpnwv8a5Z98J5FFyLiIeD/iiz6OvAF9rwxsBFoioidEfEcsJ78MCV9pov4i6m4+KHTPnwauDoidiZ1XkzKK64PXf0PJAk4B7g9Kaq4+KHTPgTQ/m18GG/e91Rxfegk/nHAQ8n0UuAjyXQlxr8xIh5LprcCa4BR5GNdnFRbDJydTGfSByeLMkmaBvw+Ip7osGgU8LuC+ZakrBJdKGl1snvevutaTfG/DfhLSb+R9EtJk5LyauoDwF8CmyJiXTJfTfFfAvybpN8B1wCXJeXV0ocngWnJ9Md486bfio5fUh1wPPAboDYiNkI+oQBHJtUy6YOTRRkkHQz8M/CVYouLlFXidck3AscAE4GN5A+DQPXED/n7g4YDpwCXAnck39KrqQ8A5/LmXgVUV/yfBj4bEUcDnwUWJeXV0ofzgdmSHiV/aOe1pLxi45c0BLgLuCQiXu2qapGyfe6Dk0V5jiF/DPAJSRvIDzvymKQ/o0qGJImITRGxKyLeAL7Lm7unVRF/ogX4SeStAN4gP5Ba1fRB0gHAXwM/KiiumviBGcBPkukfU2Wvo4h4JiLeHxEnkk/Yv00WVWT8kg4knyhujYj2v/smSSOT5SOB9sOxmfTByaIMEdEcEUdGRF1E1JH/p5wQEX8gP/zIdEmDJI0BxgIr+jDcotpfXIkPk98dhyqJP3E3cBqApLcBB5EfcbOa+nAG8ExEtBSUVVP8LwCTk+nTgPZDaVXRB0lHJr8HAP8CfDtZVHHxJ3vNi4A1EXFdwaIl5JM2ye97Csp7vg99eZa/0n/If+PYCLxOPjF8qsPyDSRXQyXz/0z+G8pakitFKi1+4GagGVidvKhGVmr8XfThIOAW8onuMeC0Su1DZ68h4CbgH4vUr6j4u/gfvAd4lPxVN78BTqzUPnQS/8Xkryp6FriaZDSLCo3/PeQPI60GViU/ZwGHA8vIJ+plwIgs++DhPszMLJUPQ5mZWSonCzMzS+VkYWZmqZwszMwslZOFmZmlcrIwAyS1dpifKemGlHWmSZqbUqdB0s86WXZJMipA+7wk/VfhCKjdJemBwlFIzfaVk4VZN0XEkoi4eh82cQlwcMH8WcAT0fVQDqW6GfinHtiOGeBkYZZK0lsk3SVpZfJzalK+e+9D0jHJ8x1WSrqyw57KEEl3SnpG0q3JHsRngKOAByU9mNQ7jzfvwkXSJ5MBH5+QdHNSdpOkG5PnG/yPpMnJgJBrJN1U0OYS8mNPmfWIA/o6ALMKMVgFDyECRpD/wAW4Hvh6RPy3pD8n/yyTd3RY/3rg+oi4XdI/dlh2PHAc+SEyfg2cGhHfkPQ5YEpE/DGpdyrwDwCSjiN/F+6pEfFHSSMKtjec/BAb04B7k/UuAFZKmhgRqyLi5WS4h8MjYnM3/yZmuzlZmOVtj4iJ7TOSZgL1yewZwLH5IXoAOFTS0A7rv5s3nydwG/lhu9utiGQMqCQh1QH/XSSGEZF/XgHkk8Gd7YkkIgqfx3BvRISkZvJDnDcn234q2XZ70nuR/N6Lk4XtM5KKZP0AAAEdSURBVCcLs3QDgHdHxPbCwoLkkWZnwfQuOn/ftUkaEPkRgUXnw0q3b++NDtt+o8O2a4A9YjbrLp+zMEt3P3Bh+4ykiUXqLOfNp61NL3G7W8k/S6HdWuCtyfQy4BxJhydtjqAMyUilf0Z+sEuzfeZkYZbuM0B9crL5aaDjOQnIX9n0OUkryD8zeUsJ210I/KLgBPd9QANARDwFzAd+KekJ4LqiW+jcicDyiGgrcz2zojzqrFkPSO6X2J6cS5gOnBsRjWVuYyTww4h4Xw/Ecz2wJCKW7eu2zMDnLMx6yonADcnhn1fIP7azLBGxUdJ3JR3aA/daPOlEYT3JexZmZpbK5yzMzCyVk4WZmaVysjAzs1ROFmZmlsrJwszMUv1/8SX6wmcPIAQAAAAASUVORK5CYII=\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "# let's plot the histogram\n", "figure = plt.figure()\n", "ax = figure.add_subplot(1,1,1)\n", "male_ds = dataset[:5000].rename(index=str, columns={\"Height\": \"Male\"}).plot.hist(ax=ax)\n", "female_ds = dataset[5000:].rename(index=str, columns={\"Height\": \"Female\"}).plot.hist(ax=ax)\n", "ax.grid()\n", "ax.set_xlabel('Height(cm)')" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "subslide" } }, "source": [ "* Assume that the height is a *discrete* variable (we use quantization on the dataset, only integers).\n", "* $Pr(Male) = Pr(Female) = 0.5$\n", "* $Pr(170cm|Female) \\approx \\frac{800}{5000} = 0.16$ (in the presentation $0.1$)\n", "* $Pr(170cm|Male) \\approx \\frac{1300}{5000} = 0.26$ (in the presentation $0.3$)\n", "* Using **Bayes rule**: $$ Pr(Female|170cm) = \\frac{Pr(170cm|Female)Pr(Female)}{Pr(170cm)} = \\frac{Pr(170cm|Female)Pr(Female)}{Pr(170cm|Male)Pr(Male) + Pr(170cm|Female)Pr(Female)} = \\frac{0.16 * 0.5}{0.16 * 0.5 + 0.26 * 0.5} = 0.38$$ ($0.25$ in the presentation)" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "## Mean & Variance\n", "---\n", "### Mean (Expectation) - $\\mu$\n", "---\n", "The mean is the proability weighted average of all possible values.\n", "* Discrete Variables:\n", " * $\\mathbb{E}[X] = \\sum_{x \\in X} xp(x)$\n", " * $\\mathbb{E}[f(X)] = \\sum_{x \\in X} f(x)p(x)$\n", "* Continuous Variables:\n", " * $\\mathbb{E}[X] = \\int_x xp(x)$\n", " * $\\mathbb{E}[f(X)] = \\int_{x \\in X} f(x)p(x)$\n", "* Example: the mean of a fair six-sided dice: $$\\mathbb{E}[X] = 1*\\frac{1}{6} + 2*\\frac{1}{6} + 3*\\frac{1}{6} + 4*\\frac{1}{6} + 5*\\frac{1}{6} + 6*\\frac{1}{6} = 3.5$$" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "subslide" } }, "source": [ "* **The Law of Total Expectation (Smoothing Theorem)**: $$\\mathbb{E}[X] = \\mathbb{E}\\big[\\mathbb{E}[X|Y] \\big] $$\n", " * Proof: $$\\mathbb{E}\\big[\\mathbb{E}[X|Y] \\big] = \\mathbb{E}\\big[ \\sum_x x \\cdot P(X=x|Y) \\big] = \\sum_y \\big[ \\sum_x x \\cdot P(X=x|Y) \\big] \\cdot P(Y=y)$$ $$ = \\sum_y \\big[ \\sum_x x \\cdot P(X=x|Y) \\cdot P(Y=y)\\big] = \\sum_y \\big[ \\sum_x x \\cdot P(X=x, Y=y)\\big] $$ $$ = \\sum_x x\\cdot \\big[ \\sum_y \\cdot P(X=x, Y=y)\\big] = \\sum_x x\\cdot \\big[P(X=x)\\big] = \\mathbb{E}[X]$$\n", " * Example: Suppose that two factories supply light bulbs to the market. Factory X's bulbs work for an average of 5000 hours, whereas factory Y's bulbs work for an average of 4000 hours. It is known that factory X supplies 60% of the total bulbs available. What is the expected length of time that a purchased bulb will work for? $$\\mathbb{E}[L] = \\mathbb{E}\\big[\\mathbb{E}[L|factory] \\big] = \\mathbb{E}[L|X] \\cdot P(X) + \\mathbb{E}[L|Y] \\cdot P(Y) = 5000 \\cdot 0.6 + 4000 \\cdot 0.4 = 4600$$" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "subslide" } }, "source": [ "### Variance - $\\sigma^2$\n", "---\n", "The variance is a measure of the \"spread\" of the distribution (can also be considered as confidence).\n", "\n", "* $var[X] = \\mathbb{E}[(X - \\mu)^2] = \\sum (x-\\mu)^2 p(x) = \\sum x^2 p(x) + \\mu^2 \\sum p(x) -2 \\mu \\sum x p(x) = \\mathbb{E}[X^2] - \\mu^2$\n", "* **The Standard Deviation** - $std[X] = \\sqrt{var[X]}$\n", "* Example: the variance of a fair six-sided dice: $$var[X] = \\sum_{i=1}^6 \\frac{1}{6} (i-3.5)^2$$ $$ E[X^2] = 1^2 *\\frac{1}{6} + 2^2 *\\frac{1}{6} + 3^2 *\\frac{1}{6} + 4^2*\\frac{1}{6} + 5^2*\\frac{1}{6} + 6^2*\\frac{1}{6} = \\frac{91}{6} $$ $$ var[X] = E[X^2] -\\mu ^2 = \\frac{91}{6} - 3.5^2 \\approx 2.92 $$" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "subslide" } }, "source": [ "### Example cont.\n", "---\n", "What is the mean and variance of the heights of males? females? combined together?" ] }, { "cell_type": "code", "execution_count": 8, "metadata": { "slideshow": { "slide_type": "subslide" } }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "the mean of males' height is: 175.327 cm\n", "the variance of males' height is: 52.896 cm^2\n", "the std of males' height is: 7.273 cm\n" ] } ], "source": [ "# easy with pandas\n", "print(\"the mean of males' height is: {:.3f} cm\".format(dataset[:5000].Height.mean()))\n", "print(\"the variance of males' height is: {:.3f} cm^2\".format(dataset[:5000].Height.var()))\n", "print(\"the std of males' height is: {:.3f} cm\".format(dataset[:5000].Height.std()))" ] }, { "cell_type": "code", "execution_count": 9, "metadata": { "slideshow": { "slide_type": "subslide" } }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "the mean of females' height is: 161.820 cm\n", "the variance of females' height is: 46.903 cm^2\n", "the std of females' height is: 6.849 cm\n" ] } ], "source": [ "print(\"the mean of females' height is: {:.3f} cm\".format(dataset[5000:].Height.mean()))\n", "print(\"the variance of females' height is: {:.3f} cm^2\".format(dataset[5000:].Height.var()))\n", "print(\"the std of females' height is: {:.3f} cm\".format(dataset[5000:].Height.std()))" ] }, { "cell_type": "code", "execution_count": 10, "metadata": { "slideshow": { "slide_type": "subslide" } }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "the mean of total height is: 168.574 cm\n", "the variance of total height is: 95.506 cm^2\n", "the std of total' height is: 9.773 cm\n" ] } ], "source": [ "print(\"the mean of total height is: {:.3f} cm\".format(dataset.Height.mean()))\n", "print(\"the variance of total height is: {:.3f} cm^2\".format(dataset.Height.var()))\n", "print(\"the std of total' height is: {:.3f} cm\".format(dataset.Height.std()))" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "## Correlation\n", "---\n", "Correlation is a measure of linear dependency between two variables.\n", "\n", "We define **correlation** between two Random Variables (RV) as $\\sigma_{xy}$:\n", "* $\\sigma_{xy} = Cov(X,Y) = \\mathbb{E}[(X - \\mu_x)(Y - \\mu_y)] = \\mathbb{E}[XY] - \\mu_x \\mu_y$\n", "* $X, Y$ are **uncorrelated** if $\\sigma_{xy} = 0 \\leftrightarrow \\mu_{xy} = \\mu_x \\mu_y$" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "subslide" } }, "source": [ " **REMEMBER**: \n", "\n", "Independence $\\rightarrow$ Uncorrelated **BUT** Uncorrelated $\\nrightarrow$ Independence" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "subslide" } }, "source": [ "### Pearson's Correlation Coefficient (Pearson's r)\n", "---\n", "It is a measure of linear correlation between two variables X and Y, denoted $\\rho$ or $r_{xy}$:\n", "* $\\rho = r_{xy} = \\frac{\\sigma_{xy}}{\\sigma_x \\sigma_y}$\n", "* From **Cauchy-Schwarz** inequality: $$-1 \\leq \\rho \\leq 1$$\n", " * $\\rho = 0$ - no correlation\n", " * $\\rho = 1$ - positive linear correlation\n", " * $\\rho = -1$ - negative linear correlation\n", " * Reminder: Cauchy-Schwarz inequality:$|\\langle x,y \\rangle|^2 \\leq \\langle x,x \\rangle \\cdot \\langle y,y \\rangle$, equality iff $x,y$ are linealy dependent.\n", " * $\\mathbb{E}(X,Y) = \\langle X, Y \\rangle \\rightarrow |\\mathbb{E}(X,Y)|^2 \\leq \\mathbb{E}(X^2) \\mathbb{E}(Y^2) $\n", " * $$|Cov(X,Y)|^2 = |\\mathbb{E}[(X-\\mu_x) (Y-\\mu_y)]|^2 = |\\langle X-\\mu_x, Y-\\mu_y \\rangle|^2 \\leq \\langle X-\\mu_x, X-\\mu_x \\rangle \\cdot \\langle Y-\\mu_y, Y-\\mu_y \\rangle$$ $$ = \\mathbb{E}[(X - \\mu_x)^2] \\mathbb{E}[(Y - \\mu_y)^2] = Var(X) Var(Y) $$" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "subslide" } }, "source": [ "Illustration:\n", "\n", "By Kiatdd - Own work, CC BY-SA 3.0, Link" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "subslide" } }, "source": [ "#### Example \n", "---\n", "* Let $X \\in \\{-1,0,1\\} \\sim U(\\frac{1}{3})$, $Y=X^2$\n", "* $X,Y$ are clearly dependent **BUT**:\n", " * $\\mu_x = 0$, $\\mu_y = \\frac{2}{3}$\n", " * $\\mu_{xy} = \\mathbb{E}[X^3] = 0$\n", " * $Cov(X,Y) = \\mu_{xy} - \\mu_x \\mu_y = 0 - 0 * \\frac{2}{3} = 0$" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "subslide" } }, "source": [ "### Correlation DOES NOT Imply Causation\n", "Below are examples (from the presentation) that show correlated variables, but they are not neccessariy caused by one another:\n", "\n", "" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "subslide" } }, "source": [ "" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "### Outliers\n", "---\n", "* In statistics, an outlier is an observation point that is distant from other observations. \n", "* An outlier may be due to variability in the measurement or it may indicate experimental error. The latter are sometimes excluded from the data set. \n", "* An outlier can cause serious problems in statistical analyses." ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "subslide" } }, "source": [ "#### Example\n", "---\n", "If one is calculating the average temperature of 10 objects in a room, and nine of them are between 20 and 25 degrees Celsius, but an oven is at 175 °C, the median of the data will be between 20 and 25 °C but the mean temperature will be between 35.5 and 40 °C. In this case, the median better reflects the temperature of a randomly sampled object (but not the temperature in the room) than the mean. Naively interpreting the mean as \"a typical sample\", equivalent to the median, is incorrect. As illustrated in this case, outliers may indicate data points that belong to a different population than the rest of the sample set." ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "subslide" } }, "source": [ "### Correlation - Sensitive to Ouliers\n", "All the examples (from the presentation) below share the same Pearson's r:\n", "" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "## Vectors of Random Variables\n", "---\n", "* Let $\\overline{X}$ be a d-dimensional **random** vector\n", " * $\\overline{X} = [x_1, x_2, ..., x_d]$\n", "* The d-dimensional **mean vector $\\overline{\\mu}$** is:\n", " * $\\overline{\\mu} = \\mathbb{E}[\\overline{X}] = [\\mathbb{E}[x_1], \\mathbb{E}[x_2], ..., \\mathbb{E}[x_d]] = [\\mu_1, \\mu_2, ..., \\mu_d]$" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "subslide" } }, "source": [ "* The **covariance matrix $\\Sigma$** is defined as the (**square**) matrix, where each component $\\sigma_{ij}$ is the covariance of $x_i, x_j$:\n", " * $\\sigma_{ij} = \\mathbb{E}[(x_i - \\mu_i)(x_j - \\mu_j)]$\n", " * $$\\Sigma = \n", " \\begin{pmatrix}\n", " \\sigma_1 ^2 & \\sigma_{1,2} & \\cdots & \\sigma_{1,d} \\\\\n", " \\sigma_{2,1} & \\sigma_2^2 & \\cdots & \\sigma_{2,d} \\\\\n", " \\vdots & \\vdots & \\ddots & \\vdots \\\\\n", " \\sigma_{d,1} & \\sigma_{d,2} & \\cdots & \\sigma_d^2 \n", " \\end{pmatrix}$$\n", " " ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "## Multivariate Normal Distribution\n", "---\n", "* $ x \\sim N_d (\\mu, \\Sigma)$\n", "* $$f(x) = \\frac{1}{(2\\pi)^{\\frac{d}{2}} |\\Sigma|^{\\frac{1}{2}}} e^{- \\frac{1}{2}(x - \\mu)^{T} \\Sigma^{-1} (x - \\mu)}$$" ] }, { "cell_type": "code", "execution_count": 11, "metadata": { "slideshow": { "slide_type": "subslide" } }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "mu:\n", "[0. 0. 0. 0. 0.]\n", "Sigma:\n", "[[0.04800934 0.0189633 0.13069076 0.03377338 0.08095299]\n", " [0.0189633 0.00760599 0.05172954 0.01336806 0.03204252]\n", " [0.13069076 0.05172954 0.35660823 0.0921296 0.22082975]\n", " [0.03377338 0.01336806 0.0921296 0.02390833 0.05706729]\n", " [0.08095299 0.03204252 0.22082975 0.05706729 0.13688724]]\n", "draw a sample from each variable:\n", "[-0.24298271 -0.09838847 -0.72393916 -0.176606 -0.41587365]\n" ] } ], "source": [ "num_samples = 1000\n", "num_variables = 5\n", "mu = np.random.random(size=(1, num_variables))\n", "sigma = mu * mu.T + np.eye(num_variables) * 1e-4\n", "# generate multivariate distribution\n", "mult_var = np.random.multivariate_normal(np.zeros(num_variables), sigma)\n", "print(\"mu:\")\n", "print(np.zeros(num_variables))\n", "print(\"Sigma:\")\n", "print(sigma)\n", "print(\"draw a sample from each variable:\")\n", "print(mult_var)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Mahalanobis Distance\n", "---\n", "* $d = (x - \\mu)^{T} \\Sigma^{-1} (x - \\mu)$\n", "* Measures the distance from $x$ to $\\mu$ in terms of $\\Sigma$, that is, the distance between a point $p$ and a distribution $D$.\n", " * If $p$ is the mean of $D$, the distance is 0.\n", " * The distance grows as $p$ moves away from the mean along each principal component axis.\n", " * Note: if $\\Sigma = I$ , that is the **Euclidean** distance.\n", "* It normalizes the expression for difference in variances and correlations." ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "subslide" } }, "source": [ "#### Example - Bivariate\n", "---\n", "* $d=2$\n", "* $\\Sigma = \n", " \\begin{pmatrix}\n", " \\sigma_1 ^2 & \\rho \\sigma_1 \\sigma_2 \\\\\n", " \\rho \\sigma_1 \\sigma_2 & \\sigma_2^2\n", " \\end{pmatrix}$\n", "* $\\rho$ is the Pearson coefficient" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "subslide" } }, "source": [ "* Reminder: 2x2 matrix inversion $$ \\begin{pmatrix}\n", " a & b \\\\\n", " c & d\n", " \\end{pmatrix}^{-1} = \\frac{1}{ad-bc} \\begin{pmatrix}\n", " d & -b \\\\\n", " -c & a\n", " \\end{pmatrix} $$\n", "
\n", "* $$\\Sigma^{-1} = \\frac{1}{\\sigma_1^2 \\sigma_2^2 (1 - \\rho^2)} \\begin{pmatrix} \\sigma_2 ^2 & - \\rho \\sigma_1 \\sigma_2 \\\\ - \\rho \\sigma_1 \\sigma_2 & \\sigma_1^2 \\end{pmatrix}$$\n", "
\n", "* $$f(x_1, x_2) = \\frac{1}{2\\pi \\sigma_1 \\sigma_2 \\sqrt{1 - \\rho^2}} e^{- \\frac{1}{2(1-\\rho^2)}(z_1^2 -2\\rho z_1 z_2 + z_2^2)}$$\n", " * $z_i = \\frac{x_i - \\mu_i}{\\sigma_i}$ (also called **standardization**)" ] }, { "cell_type": "code", "execution_count": 12, "metadata": { "slideshow": { "slide_type": "skip" } }, "outputs": [], "source": [ "from matplotlib import cm\n", "from mpl_toolkits.mplot3d import Axes3D\n", "from scipy.stats import multivariate_normal\n", "\n", "def plot_3d_normal_dist():\n", " # Our 2-dimensional distribution will be over variables X and Y\n", " N = 60\n", " X = np.linspace(-3, 3, N)\n", " Y = np.linspace(-3, 4, N)\n", " X, Y = np.meshgrid(X, Y)\n", "\n", " # Mean vector and covariance matrix\n", " mu = np.array([0., 1.])\n", " Sigma = np.array([[ 1. , -0.5], [-0.5, 1.5]])\n", "\n", " # Pack X and Y into a single 3-dimensional array\n", " pos = np.empty(X.shape + (2,))\n", " pos[:, :, 0] = X\n", " pos[:, :, 1] = Y\n", "\n", " F = multivariate_normal(mu, Sigma)\n", " Z = F.pdf(pos)\n", "\n", " # Create a surface plot and projected filled contour plot under it.\n", " fig = plt.figure(figsize=(8,5))\n", " ax = fig.gca(projection='3d')\n", " ax.plot_surface(X, Y, Z, rstride=3, cstride=3, linewidth=1, antialiased=True,\n", " cmap=cm.viridis)\n", "\n", " # cset = ax.contourf(X, Y, Z, zdir='z', offset=-0.15, cmap=cm.viridis)\n", "\n", " # Adjust the limits, ticks and view angle\n", " ax.set_zlim(-0.15,0.2)\n", " ax.set_zticks(np.linspace(0,0.2,5))\n", " ax.view_init(27, -21)\n", "\n", " plt.show()" ] }, { "cell_type": "code", "execution_count": 13, "metadata": { "slideshow": { "slide_type": "subslide" } }, "outputs": [ { "data": { "application/javascript": [ "/* Put everything inside the global mpl namespace */\n", "window.mpl = {};\n", "\n", "\n", "mpl.get_websocket_type = function() {\n", " if (typeof(WebSocket) !== 'undefined') {\n", " return WebSocket;\n", " } else if (typeof(MozWebSocket) !== 'undefined') {\n", " return MozWebSocket;\n", " } else {\n", " alert('Your browser does not have WebSocket support. ' +\n", " 'Please try Chrome, Safari or Firefox ≥ 6. ' +\n", " 'Firefox 4 and 5 are also supported but you ' +\n", " 'have to enable WebSockets in about:config.');\n", " };\n", "}\n", "\n", "mpl.figure = function(figure_id, websocket, ondownload, parent_element) {\n", " this.id = figure_id;\n", "\n", " this.ws = websocket;\n", "\n", " this.supports_binary = (this.ws.binaryType != undefined);\n", "\n", " if (!this.supports_binary) {\n", " var warnings = document.getElementById(\"mpl-warnings\");\n", " if (warnings) {\n", " warnings.style.display = 'block';\n", " warnings.textContent = (\n", " \"This browser does not support binary websocket messages. \" +\n", " \"Performance may be slow.\");\n", " }\n", " }\n", "\n", " this.imageObj = new Image();\n", "\n", " this.context = undefined;\n", " this.message = undefined;\n", " this.canvas = undefined;\n", " this.rubberband_canvas = undefined;\n", " this.rubberband_context = undefined;\n", " this.format_dropdown = undefined;\n", "\n", " this.image_mode = 'full';\n", "\n", " this.root = $('
');\n", " this._root_extra_style(this.root)\n", " this.root.attr('style', 'display: inline-block');\n", "\n", " $(parent_element).append(this.root);\n", "\n", " this._init_header(this);\n", " this._init_canvas(this);\n", " this._init_toolbar(this);\n", "\n", " var fig = this;\n", "\n", " this.waiting = false;\n", "\n", " this.ws.onopen = function () {\n", " fig.send_message(\"supports_binary\", {value: fig.supports_binary});\n", " fig.send_message(\"send_image_mode\", {});\n", " if (mpl.ratio != 1) {\n", " fig.send_message(\"set_dpi_ratio\", {'dpi_ratio': mpl.ratio});\n", " }\n", " fig.send_message(\"refresh\", {});\n", " }\n", "\n", " this.imageObj.onload = function() {\n", " if (fig.image_mode == 'full') {\n", " // Full images could contain transparency (where diff images\n", " // almost always do), so we need to clear the canvas so that\n", " // there is no ghosting.\n", " fig.context.clearRect(0, 0, fig.canvas.width, fig.canvas.height);\n", " }\n", " fig.context.drawImage(fig.imageObj, 0, 0);\n", " };\n", "\n", " this.imageObj.onunload = function() {\n", " fig.ws.close();\n", " }\n", "\n", " this.ws.onmessage = this._make_on_message_function(this);\n", "\n", " this.ondownload = ondownload;\n", "}\n", "\n", "mpl.figure.prototype._init_header = function() {\n", " var titlebar = $(\n", " '
');\n", " var titletext = $(\n", " '
');\n", " titlebar.append(titletext)\n", " this.root.append(titlebar);\n", " this.header = titletext[0];\n", "}\n", "\n", "\n", "\n", "mpl.figure.prototype._canvas_extra_style = function(canvas_div) {\n", "\n", "}\n", "\n", "\n", "mpl.figure.prototype._root_extra_style = function(canvas_div) {\n", "\n", "}\n", "\n", "mpl.figure.prototype._init_canvas = function() {\n", " var fig = this;\n", "\n", " var canvas_div = $('
');\n", "\n", " canvas_div.attr('style', 'position: relative; clear: both; outline: 0');\n", "\n", " function canvas_keyboard_event(event) {\n", " return fig.key_event(event, event['data']);\n", " }\n", "\n", " canvas_div.keydown('key_press', canvas_keyboard_event);\n", " canvas_div.keyup('key_release', canvas_keyboard_event);\n", " this.canvas_div = canvas_div\n", " this._canvas_extra_style(canvas_div)\n", " this.root.append(canvas_div);\n", "\n", " var canvas = $('');\n", " canvas.addClass('mpl-canvas');\n", " canvas.attr('style', \"left: 0; top: 0; z-index: 0; outline: 0\")\n", "\n", " this.canvas = canvas[0];\n", " this.context = canvas[0].getContext(\"2d\");\n", "\n", " var backingStore = this.context.backingStorePixelRatio ||\n", "\tthis.context.webkitBackingStorePixelRatio ||\n", "\tthis.context.mozBackingStorePixelRatio ||\n", "\tthis.context.msBackingStorePixelRatio ||\n", "\tthis.context.oBackingStorePixelRatio ||\n", "\tthis.context.backingStorePixelRatio || 1;\n", "\n", " mpl.ratio = (window.devicePixelRatio || 1) / backingStore;\n", "\n", " var rubberband = $('');\n", " rubberband.attr('style', \"position: absolute; left: 0; top: 0; z-index: 1;\")\n", "\n", " var pass_mouse_events = true;\n", "\n", " canvas_div.resizable({\n", " start: function(event, ui) {\n", " pass_mouse_events = false;\n", " },\n", " resize: function(event, ui) {\n", " fig.request_resize(ui.size.width, ui.size.height);\n", " },\n", " stop: function(event, ui) {\n", " pass_mouse_events = true;\n", " fig.request_resize(ui.size.width, ui.size.height);\n", " },\n", " });\n", "\n", " function mouse_event_fn(event) {\n", " if (pass_mouse_events)\n", " return fig.mouse_event(event, event['data']);\n", " }\n", "\n", " rubberband.mousedown('button_press', mouse_event_fn);\n", " rubberband.mouseup('button_release', mouse_event_fn);\n", " // Throttle sequential mouse events to 1 every 20ms.\n", " rubberband.mousemove('motion_notify', mouse_event_fn);\n", "\n", " rubberband.mouseenter('figure_enter', mouse_event_fn);\n", " rubberband.mouseleave('figure_leave', mouse_event_fn);\n", "\n", " canvas_div.on(\"wheel\", function (event) {\n", " event = event.originalEvent;\n", " event['data'] = 'scroll'\n", " if (event.deltaY < 0) {\n", " event.step = 1;\n", " } else {\n", " event.step = -1;\n", " }\n", " mouse_event_fn(event);\n", " });\n", "\n", " canvas_div.append(canvas);\n", " canvas_div.append(rubberband);\n", "\n", " this.rubberband = rubberband;\n", " this.rubberband_canvas = rubberband[0];\n", " this.rubberband_context = rubberband[0].getContext(\"2d\");\n", " this.rubberband_context.strokeStyle = \"#000000\";\n", "\n", " this._resize_canvas = function(width, height) {\n", " // Keep the size of the canvas, canvas container, and rubber band\n", " // canvas in synch.\n", " canvas_div.css('width', width)\n", " canvas_div.css('height', height)\n", "\n", " canvas.attr('width', width * mpl.ratio);\n", " canvas.attr('height', height * mpl.ratio);\n", " canvas.attr('style', 'width: ' + width + 'px; height: ' + height + 'px;');\n", "\n", " rubberband.attr('width', width);\n", " rubberband.attr('height', height);\n", " }\n", "\n", " // Set the figure to an initial 600x600px, this will subsequently be updated\n", " // upon first draw.\n", " this._resize_canvas(600, 600);\n", "\n", " // Disable right mouse context menu.\n", " $(this.rubberband_canvas).bind(\"contextmenu\",function(e){\n", " return false;\n", " });\n", "\n", " function set_focus () {\n", " canvas.focus();\n", " canvas_div.focus();\n", " }\n", "\n", " window.setTimeout(set_focus, 100);\n", "}\n", "\n", "mpl.figure.prototype._init_toolbar = function() {\n", " var fig = this;\n", "\n", " var nav_element = $('
');\n", " nav_element.attr('style', 'width: 100%');\n", " this.root.append(nav_element);\n", "\n", " // Define a callback function for later on.\n", " function toolbar_event(event) {\n", " return fig.toolbar_button_onclick(event['data']);\n", " }\n", " function toolbar_mouse_event(event) {\n", " return fig.toolbar_button_onmouseover(event['data']);\n", " }\n", "\n", " for(var toolbar_ind in mpl.toolbar_items) {\n", " var name = mpl.toolbar_items[toolbar_ind][0];\n", " var tooltip = mpl.toolbar_items[toolbar_ind][1];\n", " var image = mpl.toolbar_items[toolbar_ind][2];\n", " var method_name = mpl.toolbar_items[toolbar_ind][3];\n", "\n", " if (!name) {\n", " // put a spacer in here.\n", " continue;\n", " }\n", " var button = $('');\n", " button.click(method_name, toolbar_event);\n", " button.mouseover(tooltip, toolbar_mouse_event);\n", " nav_element.append(button);\n", " }\n", "\n", " // Add the status bar.\n", " var status_bar = $('');\n", " nav_element.append(status_bar);\n", " this.message = status_bar[0];\n", "\n", " // Add the close button to the window.\n", " var buttongrp = $('
');\n", " var button = $('');\n", " button.click(function (evt) { fig.handle_close(fig, {}); } );\n", " button.mouseover('Stop Interaction', toolbar_mouse_event);\n", " buttongrp.append(button);\n", " var titlebar = this.root.find($('.ui-dialog-titlebar'));\n", " titlebar.prepend(buttongrp);\n", "}\n", "\n", "mpl.figure.prototype._root_extra_style = function(el){\n", " var fig = this\n", " el.on(\"remove\", function(){\n", "\tfig.close_ws(fig, {});\n", " });\n", "}\n", "\n", "mpl.figure.prototype._canvas_extra_style = function(el){\n", " // this is important to make the div 'focusable\n", " el.attr('tabindex', 0)\n", " // reach out to IPython and tell the keyboard manager to turn it's self\n", " // off when our div gets focus\n", "\n", " // location in version 3\n", " if (IPython.notebook.keyboard_manager) {\n", " IPython.notebook.keyboard_manager.register_events(el);\n", " }\n", " else {\n", " // location in version 2\n", " IPython.keyboard_manager.register_events(el);\n", " }\n", "\n", "}\n", "\n", "mpl.figure.prototype._key_event_extra = function(event, name) {\n", " var manager = IPython.notebook.keyboard_manager;\n", " if (!manager)\n", " manager = IPython.keyboard_manager;\n", "\n", " // Check for shift+enter\n", " if (event.shiftKey && event.which == 13) {\n", " this.canvas_div.blur();\n", " event.shiftKey = false;\n", " // Send a \"J\" for go to next cell\n", " event.which = 74;\n", " event.keyCode = 74;\n", " manager.command_mode();\n", " manager.handle_keydown(event);\n", " }\n", "}\n", "\n", "mpl.figure.prototype.handle_save = function(fig, msg) {\n", " fig.ondownload(fig, null);\n", "}\n", "\n", "\n", "mpl.find_output_cell = function(html_output) {\n", " // Return the cell and output element which can be found *uniquely* in the notebook.\n", " // Note - this is a bit hacky, but it is done because the \"notebook_saving.Notebook\"\n", " // IPython event is triggered only after the cells have been serialised, which for\n", " // our purposes (turning an active figure into a static one), is too late.\n", " var cells = IPython.notebook.get_cells();\n", " var ncells = cells.length;\n", " for (var i=0; i= 3 moved mimebundle to data attribute of output\n", " data = data.data;\n", " }\n", " if (data['text/html'] == html_output) {\n", " return [cell, data, j];\n", " }\n", " }\n", " }\n", " }\n", "}\n", "\n", "// Register the function which deals with the matplotlib target/channel.\n", "// The kernel may be null if the page has been refreshed.\n", "if (IPython.notebook.kernel != null) {\n", " IPython.notebook.kernel.comm_manager.register_target('matplotlib', mpl.mpl_figure_comm);\n", "}\n" ], "text/plain": [ "" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/html": [ "" ], "text/plain": [ "" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "# let's see how the MLE performs\n", "plot_normal_mle()" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "### Exercise 2.5 - MLE for m-Dimensional Gaussian\n", "---\n", "Given $\\{x_i\\}_{i=1}^n$ i.i.d samples of $X \\sim N(\\mu, \\Sigma)$, what is the MLE?" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "subslide" } }, "source": [ "### Solution 2.5\n", "---\n", "The final results are pretty much the same, but with vectors and matrices, though the math is a little more complicated.\n", "$$ \\hat{\\overline{\\mu}}_{MLE} = \\frac{1}{n} \\sum_{i=1}^n \\overline{x_i} $$\n", "$$ \\hat{\\Sigma}_{MLE} = \\frac{1}{n} \\sum_{i=1}^n (\\overline{x_i} - \\hat{\\overline{\\mu}}_{MLE}) (\\overline{x_i} - \\hat{\\overline{\\mu}}_{MLE})^{T}$$" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "subslide" } }, "source": [ "#### Vector & Matrix Deriviatives\n", "---\n", "* $\\nabla_x Ax = A^{T}$\n", "* $\\nabla_x x^{T} A x = (A + A^{T}) x$ \n", "* $\\frac{\\partial}{\\partial A} \\ln |A| = A^{-T}$\n", "* $\\frac{\\partial}{\\partial A} Tr[AB] = B^{T}$" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "subslide" } }, "source": [ "Using the above, we will use the following:\n", "1. $\\nabla_{\\mu} {\\mu}^{T} \\Sigma^{-1} x_i = \\Sigma^{-1} x_i$ \n", "2. $\\nabla_{\\mu} {\\mu}^{T} \\Sigma^{-1} \\mu = (\\Sigma^{-1} + {\\Sigma}^{-T}) \\mu$\n", "3. $\\frac{\\partial}{\\partial \\Sigma^{-1}} \\ln |\\Sigma^{-1}| = \\Sigma^{T} = \\Sigma$\n", "4. $\\frac{\\partial}{\\partial \\Sigma^{-1}} Tr[\\Sigma^{-1} \\sum_{i=1}^n (\\overline{x_i} - \\overline{\\mu}) (\\overline{x_i} - \\overline{\\mu})^{T}] = \\sum_{i=1}^n (\\overline{x_i} - \\overline{\\mu}) (\\overline{x_i} - \\overline{\\mu})^{T}$" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "subslide" } }, "source": [ "#### Solve for the d-dimensional case\n", "* $p(x|\\mu, \\Sigma) = \\frac{1}{(2\\pi)^{\\frac{nd}{2}} |\\Sigma|^{\\frac{n}{2}}} e^{- \\frac{1}{2}\\sum_{i=1}^n (x_i - \\mu)^{T} \\Sigma^{-1} (x_i - \\mu)}$\n", "\n", "
\n", "\n", "* $\\ln p(x|\\mu, \\Sigma) \\propto -\\frac{n}{2} \\ln |\\Sigma^{-1}| -\\frac{1}{2} \\sum_{i=1}^n (\\overline{x_i} - \\overline{\\mu})^{T} \\Sigma^{-1} (\\overline{x_i} - \\overline{\\mu}) $\n", "\n", "
\n", "\n", "* $\\nabla_{\\mu} \\sum_{i=1}^n (\\overline{x_i} - \\overline{\\mu})^{T} \\Sigma^{-1} (\\overline{x_i} - \\overline{\\mu}) = \\sum_{i=1}^{n} (-2\\Sigma^{-1} x_i + (\\Sigma^{-1} + {\\Sigma}^{-T}) \\mu) = 0 \\rightarrow \\hat{\\overline{\\mu}}_{MLE} = \\frac{1}{n} \\sum_{i=1}^n \\overline{x_i} $" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "subslide" } }, "source": [ "\n", "* **The Trace Trick** - $\\sum_{i=1}^n (\\overline{x_i} - \\overline{\\mu})^{T} \\Sigma^{-1} (\\overline{x_i} - \\overline{\\mu}) = \\sum_{i=1}^n \\textit{Trace}\\big((\\overline{x_i} - \\overline{\\mu})^{T} \\Sigma^{-1} (\\overline{x_i} - \\overline{\\mu})\\big) = \\textit{Trace}\\big(\\Sigma^{-1} \\sum_{i=1}^n (\\overline{x_i} - \\overline{\\mu}) (\\overline{x_i} - \\overline{\\mu})^{T} \\big)$\n", "\n", "
\n", "\n", "* $\\frac{\\partial}{\\partial \\Sigma^{-1}}\\big( \\frac{n}{2} \\ln |\\Sigma^{-1}| -\\frac{1}{2} \\sum_{i=1}^n (\\overline{x_i} - \\overline{\\mu})^{T} \\Sigma^{-1} (\\overline{x_i} - \\overline{\\mu}) \\big) = \\frac{\\partial}{\\partial \\Sigma^{-1}}\\big( \\frac{n}{2} \\ln |\\Sigma^{-1}| -\\frac{1}{2} \\Sigma^{-1} \\sum_{i=1}^n (\\overline{x_i} - \\overline{\\mu}) (\\overline{x_i} - \\overline{\\mu})^{T} \\big) = $ $$ \\frac{n}{2} \\Sigma - \\frac{1}{2} \\sum_{i=1}^n (\\overline{x_i} - \\overline{\\mu}) (\\overline{x_i} - \\overline{\\mu})^{T} = 0 \\rightarrow \\hat{\\Sigma}_{MLE} = \\frac{1}{n} \\sum_{i=1}^n (\\overline{x_i} - \\hat{\\overline{\\mu}}_{MLE}) (\\overline{x_i} - \\hat{\\overline{\\mu}}_{MLE})^{T} $$" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "### Exercise 3 - MLE for Geometric Distribution\n", "---\n", "Given $\\{x_i\\}_{i=1}^n$ i.i.d samples of $X \\sim \\textit{Geom}(\\theta)$, what is the MLE?\n", "\n", "Assume:\n", "* $f(x;\\theta) = Pr(X=x) = \\theta(1-\\theta)^{x-1}$\n", "* $ 0 < \\theta < 1$\n", "* $\\mathbb{E}(X) = \\frac{1}{\\theta}$\n", "* $Var(X) = \\frac{1 - \\theta}{\\theta^2}$" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "subslide" } }, "source": [ "### Solution 3\n", "---\n", "* $L(x_1, x_2, ..., x_n; \\theta) = \\prod_{i=1}^n f(x_i;\\theta) = \\theta^n (1 -\\theta)^{\\sum_{i=1}^n (x_i - 1)}$\n", "* $l(\\theta) = \\ln L(x_1, x_2, ..., x_n; \\theta) = \\sum_{i=1}^n \\ln f(x_i;\\theta) = n\\ln(\\theta) +\\ln (1-\\theta) \\sum_{i=1}^n (x_i - 1) $\n", "* $\\theta_{MLE} = \\underset{0 < \\theta < 1}{\\mathrm{argmax}} l(\\theta)$\n", "* First derivative: $$\\frac{\\partial l(\\theta)}{\\partial \\theta} = \\frac{n}{\\theta} - \\frac{1}{1 - \\theta} \\sum_{i=1}^n (x_i - 1)$$" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "subslide" } }, "source": [ "* Second derivative: $$\\frac{\\partial^2 l(\\theta)}{\\partial \\theta^2} = -\\frac{n}{\\theta^2} - (\\frac{1}{1 - \\theta})^2 \\sum_{i=1}^n (x_i - 1)$$\n", "* $$ \\frac{n}{\\theta} - \\frac{1}{1 - \\theta} \\sum_{i=1}^n (x_i - 1) = 0 \\rightarrow \\theta_{MLE} = \\frac{1}{n}\\sum_{i=1}^n x_i$$" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "subslide" } }, "source": [ "* Plug in $\\theta_{MLE}$ in the second deriviative and keep in mind that $0 < \\theta_{MLE} < 1$:\n", " * $\\sum_{i=1}^n (x_i - 1) = n(\\frac{1}{\\theta_{MLE}} - 1) = n \\frac{1 - \\theta_{MLE}}{\\theta_{MLE}}$\n", " * $$ -\\frac{n}{\\theta_{MLE}^2} - (\\frac{1}{1 - \\theta_{MLE}})^2 \\sum_{i=1}^n (x_i - 1) = -\\frac{n}{\\theta_{MLE}^2} - (\\frac{1}{1 - \\theta_{MLE}})^2 n \\frac{1 - \\theta_{MLE}}{\\theta_{MLE}} = ... = - \\frac{n}{\\theta_{MLE}^2(1-\\theta_{MLE})} < 0$$\n", " * Since $0 < \\theta_{MLE} < 1$, we have a maximum." ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "subslide" } }, "source": [ "### Recommended Videos\n", "---\n", "#### Warning!\n", "* These videos do not replace the lectures and tutorials.\n", "* Please use these to get a better understanding of the material, and not as an alternative to the written material.\n", "\n", "#### Video By Subject\n", "* Basic Probability - Math Antics - Basic Probability\n", "* Probability for ML - Machine Learning 1/5: Probability\n", "* Maximum Likelihood Estimation (MLE)\n", " * Simple Version (6 min) - StatQuest\n", " * Complete Lecture (50 min) - Cornell CS4780" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "skip" } }, "source": [ "## Credits\n", "---\n", "* Icons from Icon8.com - https://icons8.com\n", "* Datasets from Kaggle - https://www.kaggle.com/" ] } ], "metadata": { "celltoolbar": "Slideshow", "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.6.9" } }, "nbformat": 4, "nbformat_minor": 2 }