1ex-Exercicio1-com-solucao.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {
    "run_control": {
     "frozen": false,
     "read_only": false
    }
   },
   "source": [
    "# Exercício 1"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "run_control": {
     "frozen": false,
     "read_only": false
    }
   },
   "source": [
    "Uma boa parte do trabalho das disciplinas de Física Experimental é o de calcular ajustes de reta para um dado conjunto de dados e uma teoria.\n",
    "Considere um conjunto de dados $\\{(x_1,y_1),(x_2,y_2),...,(x_n,y_n)\\}$ para o qual a teoria diz que a relação entre $x$ e $y$ é dada por:\n",
    "\n",
    "$$\n",
    "y = a + bx\n",
    "$$\n",
    "\n",
    "Sabendo que as medidas estão sujeitas a várias formas de imprecisão e que a teoria não conta com todas tais contribuições, consideremos o modelo pára os dados\n",
    "\n",
    "$$\n",
    "y_i = a + bx_i + \\xi_i\n",
    "$$\n",
    "\n",
    "onde $\\xi_i$ são valores de ruído que assumimos ser Gausiano.\n",
    "\n",
    "Os valores de $a$ e $b$ que melhor descrevem essa teoria devem ser aqueles que maximizam a verossimilhança de obsevar um certo valor de ruído, ou seja, minimizam $\\sum_i\\xi_i^2 = \\sum_i(y_i - a - bx_i)^2$.\n",
    "\n",
    "A solução desse problema é dada por:\n",
    "$$\n",
    "\\hat{a} = \\bar{y} - \\hat{b}\\bar{x}\\\\\n",
    "\\hat{b} = \\frac{\\sum_i(x_i-\\bar{x})(y_i-\\bar{y})}{\\sum_i(x_i-\\bar{x})^2}\n",
    "$$\n",
    "\n",
    "onde\n",
    "$$\n",
    "\\bar{x} = \\frac{1}{n}\\sum_i x_i\\\\\n",
    "\\bar{y} = \\frac{1}{n}\\sum_i y_i\n",
    "$$\n",
    "são as médias amostrais de $x$ e $y$\n",
    "\n",
    "Com seu conhecimento de Python, implemente funções para calcular a os coeficientes de ajuste da reta e aplique-os para o seguinte conjunto de dados:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2017-07-19T14:22:17.529203Z",
     "start_time": "2017-07-19T14:22:17.518896Z"
    },
    "run_control": {
     "frozen": false,
     "read_only": false
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\n",
      "+--------+--------+\n",
      "|    x   |    y   |\n",
      "+--------+--------+\n",
      "|  -19.8 |  -68.1 |\n",
      "|  -19.1 |  -64.9 |\n",
      "|  -17.7 |  -61.6 |\n",
      "|  -16.2 |  -59.2 |\n",
      "|  -16.1 |  -55.0 |\n",
      "|  -15.1 |  -51.2 |\n",
      "|  -13.2 |  -46.2 |\n",
      "|  -12.6 |  -43.4 |\n",
      "|  -12.2 |  -42.0 |\n",
      "|  -10.7 |  -36.4 |\n",
      "|  -10.2 |  -33.6 |\n",
      "|  -9.23 |  -30.4 |\n",
      "|  -7.88 |  -25.6 |\n",
      "|  -7.96 |  -21.7 |\n",
      "|  -6.86 |  -18.3 |\n",
      "|  -5.28 |  -16.6 |\n",
      "|  -4.51 |  -12.6 |\n",
      "|  -2.84 |  -8.42 |\n",
      "|  -2.45 |  -4.27 |\n",
      "|  -1.71 |  -2.23 |\n",
      "|  0.733 |   1.56 |\n",
      "|  0.887 |   4.14 |\n",
      "|   2.03 |   7.55 |\n",
      "|   2.29 |   13.1 |\n",
      "|   3.73 |   17.1 |\n",
      "|   5.06 |   19.2 |\n",
      "|   5.42 |   23.8 |\n",
      "|   7.19 |   26.6 |\n",
      "|    7.7 |   29.1 |\n",
      "|   8.85 |   33.6 |\n",
      "|    9.7 |   38.3 |\n",
      "|   11.9 |   40.2 |\n",
      "|   12.0 |   45.3 |\n",
      "|   12.5 |   44.6 |\n",
      "|   14.4 |   51.6 |\n",
      "|   14.4 |   54.3 |\n",
      "|   16.1 |   57.5 |\n",
      "|   16.0 |   61.3 |\n",
      "|   17.3 |   62.8 |\n",
      "|   19.1 |   68.0 |\n",
      "|   20.4 |   72.1 |\n",
      "+--------+--------+\n"
     ]
    }
   ],
   "source": [
    "from numpy.random import randn, seed\n",
    "seed(42)\n",
    "a_real, b_real = 1.75, 3.5\n",
    "amostra = {x: y for (x, y) in zip([x_i + 0.5*randn() for x_i in range(-20,21)],\n",
    "                                  [a_real + b_real*x_i + randn() for x_i in range(-20,21)])}\n",
    "\n",
    "print(\"\"\"\n",
    "+--------+--------+\n",
    "|    x   |    y   |\n",
    "+--------+--------+\n",
    "\"\"\" + \"\\n\".join([f\"| {float(x):6.3} | {y:6.3} |\" for (x,y) in amostra.items()]) +\\\n",
    "\"\"\"\n",
    "+--------+--------+\"\"\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "run_control": {
     "frozen": false,
     "read_only": false
    }
   },
   "source": [
    "No dicionário chamado `amostra`, as chaves são os medidas de $x$ e os valores são as de $y$. \n",
    "\n",
    "A linha \n",
    "\n",
    "```python\n",
    "from numpy.random import randn, seed\n",
    "```\n",
    "\n",
    "está importando uma função de uma biblioteca para gerar números aleatórios com distribuição normal e uma função para definir o estado inicial do gerador de números (pseudo)aleatórios, para que todos tenham os mesmos dados.\n",
    "\n",
    "Note que, como os dados são fictícios, sabemos exatamente os valores dos coeficientes linear e angular definidos nas variáveis `a_real` e `b_real`.\n",
    "\n",
    "Para entender a função `zip`, veja o exemplo a seguir:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2017-07-19T14:20:46.146502Z",
     "start_time": "2017-07-19T14:20:46.133071Z"
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\n",
      "tupla1 = (1, 2, 3, 4, 5)\n",
      "tupla2 = (10, 20, 30, 40, 50)\n",
      "zip12 = [(1, 10), (2, 20), (3, 30), (4, 40), (5, 50)]\n",
      "unzip1 = (1, 2, 3, 4, 5) == tupla1\n",
      "unzip2 = (10, 20, 30, 40, 50) == tupla2\n",
      "\n"
     ]
    }
   ],
   "source": [
    "tupla1 = (1,2,3,4,5)\n",
    "tupla2 = (10,20,30,40,50)\n",
    "zip12 = list(zip(tupla1, tupla2))\n",
    "unzip1, unzip2 = zip(*zip12)\n",
    "print(f\"\"\"\n",
    "tupla1 = {tupla1}\n",
    "tupla2 = {tupla2}\n",
    "zip12 = {zip12}\n",
    "unzip1 = {unzip1} == tupla1\n",
    "unzip2 = {unzip2} == tupla2\n",
    "\"\"\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "run_control": {
     "frozen": false,
     "read_only": false
    }
   },
   "source": [
    "## Solução do Felippe"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "run_control": {
     "frozen": false,
     "read_only": false
    }
   },
   "source": [
    "A forma como vou resolver aqui, embora envolva funções das quais não falei em aula, serve mais para ilustrar como é possível resolver tarefas complexas com, relativamente, poucas linhas de código em python. \n",
    "(Lembrem-se, parte da vantagem de usar o jupyter notebook é a facilidade de descobrir coisas novas na linguagem)\n",
    "\n",
    "Vou colocar a solução aqui e em seguida ilustro algumas coisas"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2017-07-19T14:22:20.997779Z",
     "start_time": "2017-07-19T14:22:20.974270Z"
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\n",
      "+---+--------+--------+\n",
      "|   | Real   | Ajuste |\n",
      "+---+--------+--------+\n",
      "| a |   1.75 |   2.05 |\n",
      "| b |    3.5 |   3.53 |\n",
      "+---+--------+--------+\n",
      "\n"
     ]
    }
   ],
   "source": [
    "def coeficientes_da_reta(dados):\n",
    "    x_media, y_media = (sum(a)/len(a) for a in zip(*dados.items()))\n",
    "    b_hat = sum((x_media - x_i)*(y_media - y_i) \n",
    "                for (x_i, y_i) in amostra.items())/sum((x_media - x_i)**2\n",
    "                                                       for x_i in amostra.keys())\n",
    "    a_hat = y_media - b_hat*x_media\n",
    "    return a_hat, b_hat\n",
    "\n",
    "a_chapeu, b_chapeu = coeficientes_da_reta(amostra)\n",
    "print(f\"\"\"\n",
    "+---+--------+--------+\n",
    "|   | Real   | Ajuste |\n",
    "+---+--------+--------+\n",
    "| a | {a_real:6.3} | {a_chapeu:6.3} |\n",
    "| b | {b_real:6.3} | {b_chapeu:6.3} |\n",
    "+---+--------+--------+\n",
    "\"\"\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Para mostrar os dados e a reta, vou usar uma biblioteca para plotar gráficos que veremos logo no curso:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2017-07-19T14:20:51.314418Z",
     "start_time": "2017-07-19T14:20:50.914783Z"
    }
   },
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "/home/felippe/miniconda3/lib/python3.6/site-packages/matplotlib/font_manager.py:278: UserWarning: Matplotlib is building the font cache using fc-list. This may take a moment.\n",
      "  'Matplotlib is building the font cache using fc-list. '\n"
     ]
    },
    {
     "data": {
      "text/plain": [
       "<matplotlib.legend.Legend at 0x7f913a1bffd0>"
      ]
     },
     "execution_count": 4,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAcMAAAFpCAYAAADgAUNWAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMi4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvhp/UCwAAIABJREFUeJzt3Xl4VOXd//H3dxLCqhA2QfawGhQtQY0VFcVdBATZhYAiVqV1+fm4drPV6qNVW/tgraLsyKYsItoiLgU1KlHZt4CEHQKECAJJJnP//phBA4Q1k5zJzOd1XbnInHNmzrenAx/vc+7FnHOIiIjEMp/XBYiIiHhNYSgiIjFPYSgiIjFPYSgiIjFPYSgiIjFPYSgiIjFPYSgiIjFPYSgiIjEvLGFoZg+Y2TIzW2pmb5lZJTNrZmZfmtkaM5tsZgnhOJeIiEi4lTgMzawB8Bugg3PuXCAO6Av8L/CSc64lkAPcUdJziYiIlIb4MH5OZTMrAKoAW4GrgP6h/WOAPwL/PN6H1K5d2zVt2jRMJYmIxIaMjIydzrk6XtdRnpU4DJ1zm83sr8AG4ADwHyAD2OOc84cO2wQ0ONFnNW3alIULF5a0JBGRmGJmWV7XUN6F4zZpItANaAacDVQFbijm0GJnBDezYWa20MwWZmdnl7QcERGRUxaODjRXA98757KdcwXAO8AvgRpmdqjl2RDYUtybnXOvOec6OOc61KmjVr6IiJS9cIThBiDVzKqYmQGdgeXAx8CtoWPSgJlhOJeIiEjYlTgMnXNfAtOAb4Aloc98DXgEeNDMMoFawBslPZeIiEhpCEtvUufcH4A/HLF5HXBROD5fRESkNGkGGhERiXkKQxERiXkKQxERiXkKQxERiXkKQxERiXkKQxERiXkKQxERr62ZC/t3e11FTFMYioh4afFUmNgH5v3J60pimsJQRMQr306Ad+6EJr+Ea5/yupqYpjAUEfHCwlEw8x5I6sQ3l73GiM+3kZGV43VVMStci/uKiMjJSn8VPngEWl7HN6l/p//o78j3B0iI9zFhaCopTRK9rjDmqGUoIlKWPvt7MAjbdIE+4/liw4/k+wMEHBT4A6Sv2+V1hTFJYSgiUlY+fQ7m/h7a9oBeoyE+gdSkWiTE+4gzqBDvIzWpltdVxiTdJhURKW3OwUdPwfy/wvn9oNsI8MUBkNIkkQlDU0lft4vUpFq6ReoRhaGISGlyDub+Dj7/B7QfBF3+Dr7Db8qlNElUCHpMYSgiUlqcg/cfga/+BRfeCTc8d1QQSmRQGIqIlIZAAN57ADJGwyXDg+MIzbyuSo5BYSgiEm6BQpg5HBZNhI4PQuffKwgjnNrrIiLhVOiHd4bBools+cUDjPD1J2PDHq+rkhNQy1BEJFz8+fD2HbBiFptSHuHqr35Bvn81CR9najB9hFPLUEQkHPx5MGUQrJgF1z3DzGq9NZi+HFHLUESkpAoOwKQBsHYe3PQCXDiU1KwcEuJ9FPgDGkxfDigMRURKIv/H4BJM6xdA1/+D9gMBDaYvbxSGIiKn6+APMLE3bPwSbvkXnN/nsN0aTF9+KAxFRE7HgT0wvids/Q5ufRPa3uJ1RVICCkMRkVO1fzeM6w7bl0PvsdDmJq8rkhJSGIqInIp92TC2G+zKhL4TodW1XlckYaAwFBE5WXu3wZiusGcD9J8Mza/0uiIJE4WhiMjJyN0MY24OBuJt06BpR68rkjDSoHsRkRPJySLv9evIy93GymvGKgijkMJQROR4dq0lf+T1HNy7mz4HHqX7uwVkZOV4XZWEmcJQRORYslfDqBspzPuR/gVP8F2guaZWi1IKQxGJaRlZOYz4OPPo1t725TD6RnAB1nWZwtq4JOIMTa0WpcLSgcbMagAjgXMBB9wOrAImA02B9UBv55zuLYhIxMjIymHAyHTy/QES4n0/ryyxdRGM7Q7xFSHtXdrWbsmExBxNrRbFwtUy/DvwgXOuDXA+sAJ4FJjnnGsJzAu9FhGJGOnrdh29ssSmjGCv0YSqMGQO1G4JBKdWu/fKFgrCKFXiMDSzM4HLgTcAnHP5zrk9QDdgTOiwMUD3kp5LRCScUpNqkRDv++n2Z+eq3wcH1FdODAZhzSSvS5QyEo7bpElANjDKzM4HMoD7gLOcc1sBnHNbzaxucW82s2HAMIDGjRuHoRwRkZNTdGWJqyuvpvXcNDizPgyaBdUbeF2elKFw3CaNB9oD/3TO/QL4kVO4Jeqce80518E516FOnTphKEdE5OSlNEnk3sYbaP3h7VCjEQx+T0EYg8LRMtwEbHLOfRl6PY1gGG43s/qhVmF9YEcYziUiEl4vt4fdayGxWTAIq9b2uiLxQIlbhs65bcBGM2sd2tQZWA7MAtJC29KAmSU9l4hIWK2YHQxCgKHzFIQxLFxzk/4amGBmCcA6YAjBoJ1iZncAG4BeYTqXiEjJLX0b3r6TvGoNmXDeKM7f6SOlqtdFiVfCEobOue+ADsXs6hyOzxcRCatFk2DG3eyt24FOW+4m55NsEhbs+nmcocQczUAjIrHlm7Ew/VfQ9DImtnyRHH/Fw8cZSkxSGIpI7PjqdZj1a2jRGfpPpkPLhoeNM9Q0a7FL6xmKSGz4YgT8+3FofSP0Gg3xFUlpUvmncYaaZi22KQxFJPrNfxHmPQnJ3aHnSIir8NOulCaJCkHRbVIRiR5HrUDhHHzybDAIz+sNPd84LAhFDlHLUESiwlErUNxxMSmZ/4AFL8IFt0HXl8EX53WZEqEUhiISFQ5fgaKQuLm/hc0ToMPtcOML4NONMDk2haGIRIVDK1D4/X6erDCWCzb/By6+G65/Bsy8Lk8inMJQRKJCSpNEJg5pT+tJv6Rq/k649D64+kkFoZwUhaGIRIdCP+3f7wr5O6HFNQpCOSW6iS4i5V9hAbwzFLJXsrNlb0Y0eJaMDXu8rkrKEbUMRaR88+fBtNth5Ww2XvgE16SfR/7SVcEepZprVE6SWoYiUn4VHITJt8HK2XDD88yq0qNIj1LNNSonTy1DESmf8vfDpP6w7hPo8jfoMITUrBwS4n0U+AOaa1ROicJQRMqfvH0wsQ9s+By6vwIX9AeCPUo116icDoWhiJQvB3NhQi/YtBB6vA7n3XrYbs01KqdDYSgi5ceBHBjXA7YtCa48kdzV64okSigMRaR8+HEXjOsG2augz3hofb3XFUkUURiKSOTbtwPGdoPd66DfpODivCJhpDAUkcj2wxYY0xV+2AwDpkKzy72uSKKQwlBEIteejTDmZvhxJ9z2DjS5xOuKJEopDEUkMuWsDwbhgVwYNAMadvC6IoliCkMRiTw7M2FsVyjYD2mz4OwLvK5IopzCUEQ8l5GV8/NA+crbg0EYKITB78FZbb0uT2KAwlBEPJWRlcOAkenk+wOcG7+Rt6s+R4X4CjBkDtRp7XV5EiM0UbeIeCp93S7y/QGSWccY35856OJZeu1bjFgaR0ZWjtflSYxQy1BEPJWaVIsL49fyuu8ZfqAqGamjeWTqdvL9W7UMk5QZtQxFxFMprGBixWdxlWuyu/cMNlFXyzBJmVPLUES8s+4TeKsfcdUbUD1tFu3OPJuCalqGScqewlBEvLHmQ5g8AGomwaCZUK0uoGWYxBsKQxEpeyvnwNS0YG/RgTOh6uGtPy3DJGVNzwxFpGwtmwFTBsJZ50Lau0cFoYgXwhaGZhZnZt+a2ezQ62Zm9qWZrTGzyWaWEK5ziUg5tXgqTBsCDVKCU6xVVutPIkM4W4b3ASuKvP5f4CXnXEsgB7gjjOcSkfLm2/Hwzp3Q5NLgpNuVqntdkchPwhKGZtYQuAkYGXptwFXAtNAhY4Du4TiXiJRDC9+EmfdCUifoPwUqVjtsd0ZWDiM+ztQge/FMuDrQ/A14GDgj9LoWsMc55w+93gQ0CNO5RKQ8SX8VPngEWl4LvcdBhUqH7S46HZsG2YtXStwyNLMuwA7nXEbRzcUc6o7x/mFmttDMFmZnZ5e0HBGJJJ/9PRiEbbpAnwlHBSH8PB2bBtmLl8Jxm/RSoKuZrQcmEbw9+jeghpkdank2BLYU92bn3GvOuQ7OuQ516tQJQzkiEhE+fQ7m/h7a9oBeoyG++D50qUm1SIj3EWdokL14psS3SZ1zjwGPAZhZJ+Ah59wAM5sK3EowINOAmSU9l4hEnsOWX2qSCM7BR0/B/L9Cu77Q/RXwxR3z/RpkL5GgNAfdPwJMMrOngG+BN0rxXCLigaOe991xMSmrX4TP/wHtB0GXvx03CA/RIHvxWljD0Dn3CfBJ6Pd1wEXh/HwRiSyHP+8rpMLcR2HzZLjwTrjhOfBpXg8pHzQdm4ictkPP+/x+P09XGEW7zfPgkuFw7VNgxfWjE4lMCkMROW0pTRKZcPuFVPngPs7ZPg86Pgidf68glHJHYSgip6/QT0rGI7B9NnR6HK54WEEo5ZLCUEROjz8f3r4DVsyCq/8IHR/wuiKR06YwFJFT58+DKWmw+n247hm45B6vKxIpEYWhiJyaggMwaQCsnQc3vQAXDvW6IpESUxiKyMnL/xHe6gvfz4eu/wftB3pdkUhYKAxF5IQysnLIWJ1Fv8yHOGNHBtzyLzi/j9dliYSNRsSKyHFlZOXwq5HzuHD+HVTa9g1rr/iHglCijsJQRI7r21VredP+TFv7nuH++/jApXpdkkjYKQxF5Nj2ZdN/xb20ss3c7X+QT30XaVUJiUp6Zigixdu7DcZ0pcreDay+fhTtD57DPVpVQqKUwlBEjpa7GcbcHAzE26bRqmlHWnldk0gpUhiKyOFysoJBeCAHBk6Hxhd7XZFIqVMYisjPdq2FMV0hfy8MmgENUryuSKRMKAxFJCh7dbBFWJgPabOhfjuvKxIpM+pNKiIs+y6d/a9dR0GhHwa/pyCUmKMwFIlxy7+Zz9nTb+WHfEfXfY+TcbC+1yWJlDmFoUgs25xB0px+7KcivfN/z2p/PdLX7fK6KpEyp2eGIrFqQzqMvxWrVIOBeQ+zmVpUiPdpUL3EJIWhSCz6fj5M7ANn1qfioFk8v6cK6et2kapB9RKjFIYisWbtRwQm9mNPxbPZePVEzq/egJTqKAQlpumZoUgsWf1vAhP7sNpfl+ty/oc+E9eRkZXjdVUinlMYisSKFe/CpAHsrNyc/vlPkO2qU+APqMOMCApDkdiw9G2YkgZnX8CWrpPYH1+dOEMdZkRC9MxQJNotmgQz7oZGqTBgChdUPIMJQ6urw4xIEQpDkWj2zViY9Rtodhn0mwQJVYFgZxmFoMjPdJtUpBzLyMphxMeZxXeC+ep1mPVraNEZ+k/5KQhF5GhqGYqUUxlZOQwYmU6+P0BCvI8JQ1N/bu19MQL+/Ti0vhF6jYb4ip7WKhLp1DIUKafS1+0i3x8g4Di8V+j8F4NBmNwNeo9VEIqcBLUMRcqp1KRaJMT7KPAHgr1Cm9WEj5+BT5+F83pB91chTn/FRU6G/qaIlFMpTRKZMDQ12Cu0WU1SMl+GBS/BBbdB15fBF+d1iSLlhsJQpBxLaZJISuMa8O8nIH0EdLgdbnwBfHoCInIqSvw3xswamdnHZrbCzJaZ2X2h7TXNbK6ZrQn9qX7cIuEWCMCch4JBePHdcNOLCkKR0xCOvzV+4P85584BUoF7zSwZeBSY55xrCcwLvRaRcAkUwuz74OuRcOl9cP0zYOZ1VSLlUonD0Dm31Tn3Tej3vcAKoAHQDRgTOmwM0L2k5xKRkEI/zLgnOKj+8ofh6icVhCIlENZnhmbWFPgF8CVwlnNuKwQD08zqhvNcIjGrsADeuROWTYerfguX/4/XFYmUe2F7uGBm1YC3gfudcz+cwvuGmdlCM1uYnZ0drnJEopM/D6YOhmXT+SzpPjKaDPW6IpGoEJYwNLMKBINwgnPundDm7WZWP7S/PrCjuPc6515zznVwznWoU6dOOMoRiU4FB2HybbByNn8ODGbgiosZMDJd6xGKhEE4epMa8Aawwjn3YpFds4C00O9pwMySnkskZuXvh7f6wJq5fNzyCUYVXHv0zDMictrC8czwUmAgsMTMvgttexx4FphiZncAG4BeYTiXSOzJ2wcT+8CGz6H7K5yZeAMJK9N/nnlG6xGKlFiJw9A5twA4Vje2ziX9fJFYkpGVc/g6gwdzYUIv2LQQerwO591KCvw884zWIxQJC81AIxIhjlyFYtLANlzwyRDYtgR6jQpOvB2i9QhFwkthKBIhiq5CUc2/hwaz+sCB76HPeGh9g9fliUQ1haFIhDi0CsWZ/hzGJzxNrQPZ0O8taHG116WJRD2FoUiESGmSyJR+TWn07iOcWbALX/+pkHSF12WJxASFoUik2LORdnP7QeFuGDgdmlzidUUiMUNhKBIJdn8PY7oGe48OmgENO3hdkUhMURiKeG1nJoy5GfwHIG0WnH2B1xWJxByFoYiXdqyEsV2DyzGlzYZ653pdkUhM0iqgIqUsIyuHER9nHj2H6LYlMPpGwGDIHAWhiIfUMhQpRUcOpJ8wNDU4WH7LtzC2OyRUhbR3oVZzr0sViWlqGYqUoqID6X+aVHvj1zCmG1Q6M9giVBCKeE5hKFKKDg2kjzOoEO/j6iprYVx3qFITBs+BxKZelygi6DapSKlKaZL406Ta11RaSasPb4fqDWHQLDizvtfliUiIwlCklKU0SSQlPwMm3wE1k2DQTKhW1+uyRKQI3SYVKW0r58CkflC7ZXD4hIJQJOIoDEVK0/KZMGUgnHVusNdoVS3EKxKJFIYipWXxVJg6BBqkBKdYq6z1B0UilcJQpDR8OwHeuROa/BJuewcqVfe6IhE5DoWhSLgtHAUz74GkTtB/ClSs5nVFInICCkORcEp/FWbfDy2v45tL/8mIz7YcPQ2biEQcDa0QCZfP/g5zfw9tuvDNRS/Sf9Q3R0/DJiIRSS1DkXD49PlgELbtAb1G80XW3qOnYRORiKWWoUhJOAcfPQXz/wrt+kL3V8AX99M0bAX+ABXifaQmaUiFSCRTGIqcLudg7u/g839A+0HQ5W/giwMOn4YtNamWbpGKRDiFocjpcA7efwS++hdceCfc8Bz4Dn/qkNIkUSEoUk4oDEVOVSAA7z0AGaPhkuFw7VNg5nVVIlICCkORUxEohJnDYdFEMhoNgdb/jxQFoUi5p96kIier0A/vDINFE/lH4FZ6ZV7NgDe+1DhCkSigMBQ5Gf58mDYElk7ji2bDeamgBwFnGjYhEiUUhiIn4s+DKYNgxSy47hkSOj102Or1GjYhUv7pmaHI8RQcgEkDYO08uOkFuHAoKaBhEyJRRmEociz5P8LEPrB+AXQNjSUM0bAJkeiiMBQpzsEfYGJv2Pgl3PIvOL/PT7sysnLUKhSJMqUehmZ2PfB3IA4Y6Zx7trTPKVIiB/bA+J6w5Vvo+Qac2+OnXRlZOQwYma4JuEWiTKl2oDGzOGAEcAOQDPQzs+TSPKdIiezfDWO7wtZF0HvsYUEIkL5ulybgFolCpd2b9CIg0zm3zjmXD0wCupXyOUVOz75sGHMz7FgJfSfCOV2OOuTQBNzqSSoSXUr7NmkDYGOR15uAi0v5nCKnbu82GNMV9myA/pOh+ZXFHqYJuEWiU2mHYXHzVLnDDjAbBgwDaNy4cSmXI1KM3M3BFuHebXDbNGja8biHqyepSPQp7dukm4BGRV43BLYUPcA595pzroNzrkOdOnVKuRyRI+Rkwagb4MdsGDj9hEEoItGptMPwa6ClmTUzswSgLzCrlM8pcpSMrBxGfJx5+Dyiu9bCqBvh4B4YNAMa6w6+SKwq1dukzjm/mQ0H/k1waMWbzrllpXlOkSMVOxyiSnaw16g/D9JmQ/12XpcpIh4q9XGGzrk5wJzSPo/IsRw5HGL1ki9JWTEcMBj8Hpyl0T4isU4TdUvUKzocol18Fr2W/Ap88TBkjoJQRABNxyYx4NBwiHXf/Zdblj1DfKXqkDYLaiZ5XZqIRAiFoUSd4uYOTbHVpCwfDlVrQtq7UEPDeETkZwpDiSrFdpZxy2BCbzijXjAIqzfwukwRiTB6ZihR5cjOMpsWvgfjb4UajYLPCBWEIlIMhaFElaKdZTrHL+Lm5Q9ArRbBXqNn1PO6PBGJULpNKlHlUGeZ7V9N4/oVL+Kr2zY4s0yVml6XJiIRTGEoUSdl70ew/BFo0B4GTIPKNbwuSUQinG6TSnRZNAneHgqNLg62CBWEInISFIYSPb4ZC9N/FZxs+7ZpUPEMrysSkXJCYShRYcMHL8OsX5Pb4HLoPwUSqnpdkoiUIwpDKfc2znmexum/48PC9nTcMJSMLQe9LklEyhmFoZRv81+g0VdP8X7hRdxdcD/7/XGkr9vldVUiUs6oN6mUT87BJ8/Cp8+yK6kbD63pTcCMCvE+UpNqeV2diJQzCkMpf5yDeU/CgpfgggHU6voPxm784aj5SEVETpbCUMoX5+Dfj0P6K5AyBG56EXw+UpokKgRF5LQpDKX8CARgzkOw8A24+Fdw/bNg5nVVIhIFFIZSPgQK4d374Ntx8MvfwDV/UhCKSNgoDCXyFfph5j2weDJc/jBc+biCUETCSmEoEeewxXkbVgtOr7Z8Blz5W7jif7wuT0SikMJQIkrRxXmrxhcyv9lYamyYC9f8GS79jdfliUiUUhhKRDm0OG8Fl88/eIkaGxbBDc/BxXd5XZqIRDHNQCMRJTWpFmfGF/BGwl+53LeYrF/+RUEoIqVOLUOJKCn14pnf4BWqbltOVsfnaXb1nV6XJCIxQGEokeNgLoy/lTO2Z0DP12l23q1eVyQiMUJhKJFh/24Y3wO2LYFeoyC5m9cViUgMURiK937cBeO6QfYq6DMeWt/gdUUiEmMUhuKZjKwcFq1cTb8Vw6m8bwP0ewtaXO11WSISg9SbVDyRkZXDAyPncMVnabicLFZ3flNBKCKeURiKJ5YuX8pY35OcZTkMKXiEuQdbe12SiMQwhaGUvd3f03fZXdRkL4MKHmdRXLIW5BURT+mZoYTNYXOKHmttwZ2ZMOZmKhYeYF3XyXT+4Wye0IK8IuIxhaGERdE5RRPifUwYmnp0wO1YCWO7BpdjSpvNOfXO5RxvyhUROUyJbpOa2fNmttLMFpvZdDOrUWTfY2aWaWarzOy6kpcqkezQnKIBBwX+AOnrdh1+wLYlMPrG4O+D34N655Z9kSIix1DSZ4ZzgXOdc+2A1cBjAGaWDPQF2gLXA6+YWVwJzyURLDWpFgnxPuIMKsT7Dn8GuOVbGN0F4ivBkPehbhvvChURKUaJbpM65/5T5GU6cGj+rG7AJOdcHvC9mWUCFwFflOR8ErlSmiQyYWjq0c8MN34N43tC5eqQ9i4kNvW0ThGR4oTzmeHtwOTQ7w0IhuMhm0LbJIqlNEk8/Dlh1ucwoRdUrRMMwhqNvCtOROQ4ThiGZvYhUK+YXU8452aGjnkC8AMTDr2tmOPdMT5/GDAMoHHjxidRspQL6z6Bt/pB9YYwaBacWd/rikREjumEYeicO+60IGaWBnQBOjvnDgXeJqBoM6AhsOUYn/8a8BpAhw4dig1MiQwnNXQCYM2HMHkA1EyCQTOhWt2yK1JE5DSU6DapmV0PPAJc4ZzbX2TXLGCimb0InA20BL4qybnEWyc1dAJg5RyYmgZ1WsPAmVBVg+lFJPKVtDfp/wFnAHPN7DszexXAObcMmAIsBz4A7nXOFZbwXOKhEw6dAFg+E6YMhLPODT4jVBCKSDlR0t6kLY6z72ng6ZJ8vkSOQ0MnCvyBo4dOACyeCtPvgoYdYMA0qHSmN4WKiJwGzUAjJ+WYQycAvp0AM++Fph2h3ySoWM27QkVEToPCUE7aUUMnABaOgtn3Q9KV0HciJFTxpjgRkRLQqhVy+r78VzAIW14XbBEqCEWknFIYyun57GV4/2Fo0wX6jIcKlbyuSETktOk2qZy0Q+MMu/8wkQbfvgBte0CP1yCugteliYiUiMJQTkpwnOEX/JY3aBD3Ibua96BWz5Hg0/zrIlL+KQzlpKSv3cnKuH4AfFx4AcsbPMa9CkIRiRJ6Zign5hy9sv/vp5f38DCpzet4WJCISHipZSjHFwjAew9Qd8Votre9g2m17mZ889rHn5tURKScURjKsQUKYeZwWDQROj7IWZ1/z71W3IIkIiLlm8JQilfoD06vtnQadHocrngYFIQiEqUUhnI0fz68fQesmAWd/wCXPeh1RSIipUphKIfz58GUNFj9Plz3DFxyj9cViYiUOoWh/KzgAEwaAGvnwU0vwIVDva5IRKRMKAyFjKwcFq7ZRL/M/+HMbenQ9R/QfpDXZYmIlBmNM4xxGVk5DBv5Ee3/ewdVt6bz/eUvKghFJOYoDGPcN6vXM9L+wi9sDff7f80cu9zrkkREypxuk8ay/bvpv3I4Fex7hvvv4xPfxQw+cgV7EZEYoDCMVfuyYWw3qu7JZM01IzmvoB13HrmCvYhIjFAYxqK922BMV9izAfpPpmXzK2npdU0iIh5SGMaa3M0w5uZgIN42DZp29LoiERHPKQxjREZWDsuWL6HPsnuoWJALA6dD44u9LktEJCKoN2kMyMjK4bGRM+icPoQDe3ez4trxCkIRkSIUhjFg5dKFjPU9SSXyuK3gCT76oYHXJYmIRBSFYbTbvpzeS35FHAFuK/gtmXFJpGr4hIjIYfTMMJptXQRju1OhQkW23ziZLrtqkKrhEyIiR1EYRqvNGTDuFqh4JqTN4tyaSZzrdU0iIhFKt0mj0YZ0Ckd3JZdqLLn2LaiZ5HVFIiIRTWEYbb6fT+HYW9iQfwY35D5Kr0mbyMjK8boqEZGIpjCMJms/ggm9yE2oR9/837LF1aLAHyB93S6vKxMRiWgKw2ix+t8wsS/UasHGblPJja9FnEGFeJ96j4qInIA60ESDFe/C1CFwVlsYOJ3zq9RkwtBapK/bpd6jIiInISwtQzN7yMycmdUOvTYze9nMMs1ssZm1D8d5pBhL34YpaXD2BTBoJlSpCUBKk0TuvbKFglBE5CSUOAzNrBFwDbChyOYbgJahn2HAP0t6HinGoknw9lBodHFwrtHKNbyuSESkXApHy/Al4GHAFdnWDRjrgtKBGmZWPwznkkMyxsD0XwVXnbhtGlQ8w+uKRETKrRKFoZl1BTY75xYdsatjVuQvAAATtklEQVQBsLHI602hbRIOX70O7/4GWnSG/lMgoarXFYmIlGsn7EBjZh8C9YrZ9QTwOHBtcW8rZpsrZhtmNozgrVQaN258onLkixHw78eh9Y3QazTEV/S6IhGRcu+EYeicu7q47WZ2HtAMWGRmAA2Bb8zsIoItwUZFDm8IbDnG578GvAbQoUOHYgNTQua/APP+BMndoOcbEFfB64pERKLCad8mdc4tcc7Vdc41dc41JRiA7Z1z24BZwKBQr9JUINc5tzU8Jccg5+DjZ4JBeF4v6PmmglBEJIxKa5zhHOBGIBPYDwwppfNEP+dg3pOw4CVW1OvK/vbPkhKn4aEiIuEUthloQi3EnaHfnXPuXudcc+fcec65heE6T0xxLvh8cMFLTApcTZes3gx482vNNSoiEmaaji1SBQIw5yFIf4UPqt3C4wVDKHQ+zTUqIlIKFIaRKFAIs++Dr0fyeuBm7t55KwFn+DTXqIhIqdDDp0hT6IeZ98LiSXzdeCjPrLkSh+EDLm1Rm/uvbqUp1kREwkwtw0hSWADv3AmLJ8FVv8XX+bckxMcRZ5BQwacgFBEpJWoZRgp/Hky7HVbOhmv+DJf+hhRgwtBUrT4hIlLKFIaRoOAgTBkIa/4DNzwHF9/1066UJokKQRGRUqbbpF7L3w9v9YE1c8n65V8Ysb+zhk6IiJQxtQy9lLcPJvaBrM9Y3/E5rv+0Efn+VSTE+5gwNFUtQhGRMqKWoVcO5sL4HrDhC+jxOu/FXUW+P0DAobGEIiJlTC1DLxzIgXE9YNti6DUKkruRmpVDQnxwUL3GEoqIlC2FYVn7cReM6wbZq6DPeGh9AxDsKKOeoyIi3lAYlqV9O2BMV8j5Hvq+BS0PXx1LPUdFRLyhMCwrP2wJBuEPm4Or0ydd4XVFIiISojAsC3s2wpib4cdsuO1taPJLrysSEZEiFIalbff3wRbhwVwYOAMaXeh1RSIicgSFYWnamQlju0LBfkibCWf/wuuKRESkGArD0rJjZTAIA4WQ9i7UO8/rikRE5Bg06L4ULP/2c/a/fj0FhQEY/J6CUEQkwikMw2xFxn85e0YvcvONm/c9RsaBs7wuSURETkBhGE4bvyZpTj/2UYne+b9jjb+eplUTESkHFIbhkvU5jOuOq1KLgYE/soWzNK2aiEg5oQ404bDuU3irL1RvSKVBs/hrTiVNqyYiUo4oDEtqzYcweQDUTIJBM6FaXVLORCEoIlKO6DZpSax6Hyb1g9otIW02VKvrdUUiInIaFIana/lMmHwbnHVucBxhVT0bFBEprxSGp2HdR6MJTBnCvtoXBG+NVtYtURGR8kxheIrWf/gaTT+9n68Crbh863Aythd6XZKIiJSQwvBULHyTpgv+h89dWwbnP0yuP0HjCEVEooDC8GSlvwqzHyC34ZXc4x6mwCpqHKGISJTQ0IqT8dnfYe7voU0Xqt86ilGbf9Q4QhGRKKIwPJFPn4ePn4K2t0CP1yGuAilNEhSCIiJRRGF4LM7Bx0/Df5+Hdn2h2wiI0+USEYlG+te9OM4Fb4t+/jL8YiDc/HfwxXldlYiIlJISd6Axs1+b2SozW2ZmzxXZ/piZZYb2XVfS85QZ5+CDR4NBeOFQuPnlw4IwIyuHER9nkpGV42GRIiISTiVqGZrZlUA3oJ1zLs/M6oa2JwN9gbbA2cCHZtbKORfZg/ICAXjvQcgYBan3wnVPg9lPuzOychgwMp18f4CEeB8Thqbq2aGISBQoacvwbuBZ51wegHNuR2h7N2CScy7POfc9kAlcVMJzla5AIcwaHgzCjg8eFYQA6et2ke8PEHBQ4A9ojKGISJQo6TPDVsBlZvY0cBB4yDn3NdAASC9y3KbQtshU6IcZv4IlU6HTY3DFI0cFIUBqUi0S4n0U+AMaYygi5UJGRkbd+Pj4kcC5xO7Y8gCw1O/3D01JSdlR3AEnDEMz+xCoV8yuJ0LvTwRSgQuBKWaWBBydJOCO8fnDgGEAjRs3PlE54efPh3eGBife7vwHuOzBYx6a0iSRCUNTNcZQRMqN+Pj4kfXq1TunTp06OT6fr9h/h6NdIBCw7Ozs5G3bto0EuhZ3zAnD0Dl39bH2mdndwDvOOQd8ZWYBoDbBlmCjIoc2BLYc4/NfA14D6NChQ9n+H+XPgylpsPp9FjR/kMqNh5BygrekNElUCIpIeXJuLAchgM/nc3Xq1Mndtm3bucc8poTnmAFcBWBmrYAEYCcwC+hrZhXNrBnQEviqhOcKr4IDMKk/rH6fJwO3M2h5BwaMTFcvURGJNr5YDsJDQtfgmJlX0jB8E0gys6XAJCDNBS0DpgDLgQ+AeyOqJ2n+jzCxN2TO46NWv2NMwdXqFCMiUkri4uJS2rRpk9yiRYu2rVu3Tv7jH/94VmHhqUVCz549m44aNarUbsuVqAONcy4fuO0Y+54Gni7J55eKgz8Eg3Djl3DLv6he41oSVqSrU4yISCmpWLFiYOXKlcsBNm/eHN+rV6+k3NzcuJdeeqnYx2deiK2eRQf2wLhbYONX0PMNOL/PT51iHry2tcYNiogACzJ3Vv3fD1bWW5C5s2q4P7tBgwb+kSNHrh81alTdQCDAqlWrElJSUlonJyefk5ycfM7cuXOrAgQCAQYNGtS4efPmbTt16tRi586dPzXeZs6cecY555yT3KpVq+RevXo1PXDggAHcc889DZo3b962VatWycOGDWt4KnXFznRs+3fDuO6wfTn0HgvndPlplzrFiIgELcjcWfX20V+38hcGfG8s+D7w5uALV3dsUfvHcJ4jOTk5PxAIsHnz5vizzz7bP3/+/NVVqlRxS5YsqdivX7+kpUuXrhg3blyNzMzMiqtWrVq2adOmCuedd17bwYMH79q/f7/dddddzf7zn/+sateuXd4tt9zS9Pnnn69z11137ZozZ07iunXrlvp8Pnbu3HlKc2jGRstwXzaMuRl2rIS+Ew8LQhER+dlnmTvP8BcGfAEH/sKA77PMnWeUxnmCgxAgPz/f+vfv3zTUymu+du3aSgCffvrpGb17994dHx9P06ZNCy655JK9AIsWLarUsGHDvHbt2uUBDB48eNeCBQvOqFmzZmHFihUDffv2bTJmzJga1apVC5xKPVEdhhlZOYz64AsOvH497FoL/SdBq2u9LktEJGJd2qL23vg4X8BnEB/nC1zaovbecJ9j+fLlCXFxcTRo0MD/9NNPn1W3bt2CFStWLF+yZMnygoKCn3LJipn85FCIHqlChQp89913K3r27LlnxowZNTp16tTyVGqK2jDMyMrhwZHv0enzwQT2bGTVNaOg+VVelyUiEtE6tqj945uDL1x91xXNN5fGLdItW7bE33nnnU2GDBmyw+fzkZubG1e/fv2CuLg4XnnllVqHepleccUVe6dOnVrT7/eTlZVVIT09/QyACy644ODmzZsTli5dWhFg7NixtS677LK9ubm5vt27d8f16dMn99VXX924YsWKKqdSV9Q+M1y2fAnjfH+kBvsYXPAYnfa3pLXXRYmIlAMdW9T+MZwhmJeX52vTpk2y3++3uLg416dPn11/+MMftgPcf//9O3r27Nl8xowZiR07dtxbuXLlAMDAgQP3zJs378zWrVu3bdas2cGLLrpoL0CVKlXcq6++ur5Xr17NCwsLOf/88/c/9NBD2Tt27Ijv0qVLi7y8PAN46qmnNp5KjXasJqcXOnTo4BYuXFjyD9q1lvw3u3BgXy5pBY+yMq6leoqKSNQyswznXIfi9i1atGj9+eefv7Osa4pEixYtqn3++ec3LW5f9LUMs1fD2K4kBPLI7DaFa3Lr8TvNIyoiIscRXWG4fTmM7QoYDH6P5LOSSfa6JhERiXjR04Fm6yIYfRP44mHwe3CWYlBERE5OVIThioWfcPCNm8j3VQoGYZ1WXpckIiLlSLkPw4ysHH4/cwlr82tyww+PkrGvptcliYhIOVPuwzB93S4y/El0yX+a9f7aWnVCREROWbkPw9SkWiTE+/CZT6tOiIhEqLFjx9Yws5Rvv/22EsD69esrXH/99Umn81kvv/xyrfXr11cIZ33lPgy16oSISOSbNGlSzfbt2+8bN25cTYCmTZsWfPDBB+tO57PGjx9fe8OGDQrDI6U0SeTeK1soCEVEIlBubq5v4cKF1UaNGrV++vTpiQCrVq1KaNmyZVsItvQGDRrU+NDxV155ZYvZs2ef4ff76dmzZ9OWLVu2bdWqVfKTTz5Zd9SoUYlLly6tMmjQoKQ2bdok79u3z+bPn1/lwgsvbN22bdtzOnbs2DIrK+uUgzK6xhmKiMixzbi3ETuWn9KcnSdUN3k/3Uccd+qzCRMm1OjUqVNuu3bt8mrUqFG4YMGCKnXq1PGf6KO/+OKLKlu3bq2wZs2aZQA7d+6Mq127duE///nPun/96183Xn755fvz8vLsN7/5TeP33nsv8+yzz/a//vrriQ899FCDqVOnrj+V/xlR0TIUEZHINWXKlJr9+vXLAejZs+fuQ7dKT6RNmzZ5GzdurJiWltZo2rRpZyYmJhYeeczixYsrrlmzpvJVV13Vqk2bNsnPP/98/S1btqhlKCIix3CCFlxp2LZtW1x6evqZq1evrjx8+HAKCwvNzNwDDzyw49Ax8fHxLhD4efnBvLw8H0CdOnUKly5dunz69OlnvvLKK3UnT55c88gWn3POWrRoceC7775bWZI61TIUEZFSM27cuMQePXrs2rJly5LNmzcv2bZt2+KGDRvmr1+/PuHQMc2bN89ftmxZlcLCQjIzMyssXry4KsDWrVvjCwsLGTx48J6nnnpq85IlS6oAVKtWrTA3NzcOoF27dgd3794d/+GHH1YFyMvLs4ULF1Y61ToVhiIiUmqmTp1aq0ePHjlFt3Xr1i3nL3/5S30zcwDXXHPNvkaNGuW1bt267X333dcoOTl5PwSHX3Ts2LF1mzZtkm+//fZmf/rTnzYBDBo0aOevf/3rJqFloZg0adLaRx99tGHr1q2T27Ztm/zpp59WO9U6o3MJJxGRGFIel3CaP39+lQcffLDR119/vaqsznm8JZzUMhQRkTL13//+t8rAgQOThg8fvt3rWg5RBxoRESlTl19++f7169cv9bqOotQyFBGRmKcwFBGJboFAIGBeF+G10DUIHGu/wlBEJLotzc7Orh7LgRgIBCw7O7s6cMxbsxHVm9TMsoGs03hrbSDSektFYk2guk6V6jo1quvUhKuuJs65OsXtyMjIqBsfHz8SOJfYbQAFgKV+v39oSkrKjuIOiKgwPF1mtvBY3Yq9Eok1geo6Varr1KiuUxOpdcWiWP2vBBERkZ8oDEVEJOZFSxi+5nUBxYjEmkB1nSrVdWpU16mJ1LpiTlQ8MxQRESmJaGkZioiInLZyG4Zm9ryZrTSzxWY23cxqFNn3mJllmtkqM7uujOvqZWbLzCxgZh2KbG9qZgfM7LvQz6uRUFdon2fX64g6/mhmm4tcoxu9qiVUz/Wha5JpZo96WUtRZrbezJaErpFnM9ub2ZtmtsPMlhbZVtPM5prZmtCfiRFSl6ffLTNrZGYfm9mK0N/D+0LbPb9eElRuwxCYC5zrnGsHrAYeAzCzZKAv0Ba4HnjFzOLKsK6lQA/gv8XsW+ucuyD086syrOmYdUXA9TrSS0Wu0RyvighdgxHADUAy0C90rSLFlaFr5GW3/NEEvzNFPQrMc861BOaFXpe10RxdF3j73fID/885dw6QCtwb+j5FwvUSynEYOuf+45zzh16mAw1Dv3cDJjnn8pxz3wOZwEVlWNcK51yZLUlyso5Tl6fXK4JdBGQ659Y55/KBSQSvlYQ45/4L7D5iczdgTOj3MUD3Mi2KY9blKefcVufcN6Hf9wIrgAZEwPWSoHIbhke4HXg/9HsDYGORfZtC2yJBMzP71sw+NbPLvC4mJNKu1/DQre83Pb5lFGnXpSgH/MfMMsxsmNfFHOEs59xWCAYAUNfjeoqKiO+WmTUFfgF8SWRfr5gS0Us4mdmHQL1idj3hnJsZOuYJgrcgJhx6WzHHh7XL7MnUVYytQGPn3C4zSwFmmFlb59wPHtdV6tfrsJMdp0bgn8CfQ+f/M/ACwf/Q8UKZXpdTdKlzbouZ1QXmmtnKUGtIji0ivltmVg14G7jfOfeDWcxOFxpxIjoMnXNXH2+/maUBXYDO7ucxIpuARkUOawhsKcu6jvGePCAv9HuGma0FWgFh6wBxOnVRBterqJOt0cxeB2aXVh0noUyvy6lwzm0J/bnDzKYTvKUbKWG43czqO+e2mll9oNh5IMuac+6nRWS9+m6ZWQWCQTjBOfdOaHNEXq9YVG5vk5rZ9cAjQFfn3P4iu2YBfc2sopk1A1oCX3lRY1FmVudQxxQzSyJY1zpvqwIi6HqF/jE45BaOM8N8GfgaaGlmzcwsgWAno1ke1gOAmVU1szMO/Q5ci7fX6UizgLTQ72nAse5IlCmvv1sWbAK+Aaxwzr1YZFdEXq+Y5Jwrlz8EO3psBL4L/bxaZN8TwFpgFXBDGdd1C8FWRR6wHfh3aHtPYBmwCPgGuDkS6vL6eh1R4zhgCbCY4D8S9T3+jt1IsKfyWoK3mj2rpUhNSaHv0KLQ98mzuoC3CN7+Lwh9t+4AahHsFbkm9GfNCKnL0+8W0JHgLdrFRf7NujESrpd+gj+agUZERGJeub1NKiIiEi4KQxERiXkKQxERiXkKQxERiXkKQxERiXkKQxERiXkKQxERiXkKQxERiXn/H3feCiY0SU3CAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 432x432 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "%matplotlib inline\n",
    "import matplotlib.pyplot as plt\n",
    "\n",
    "fig, ax = plt.subplots(figsize=(6,6))\n",
    "x, y = zip(*amostra.items())\n",
    "ax.plot(x, y, '.',label='Dados')\n",
    "ax.plot(x, [a_chapeu + b_chapeu*x_i for x_i in x], label='Ajuste')\n",
    "ax.legend(loc=(1.02,0))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "run_control": {
     "frozen": false,
     "read_only": false
    }
   },
   "source": [
    "### Recuperando os valores de $x$ e $y$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Sabendo que os valores de $x$ são as chaves e $y$ os valores do dicionário `amostra`, existem alguns modos de definir essas variáveis em uma linha de código:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2017-07-19T14:20:54.365319Z",
     "start_time": "2017-07-19T14:20:54.353314Z"
    },
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# modo 1: usando zip e o metodo amostra.items\n",
    "x_zip, y_zip = zip(*amostra.items())\n",
    "\n",
    "# modo 2: usando tuple e os metodos amostra.keys e amostra.values\n",
    "x_method, y_method = tuple(amostra.keys()), tuple(amostra.values())\n",
    "\n",
    "# modo 3: usando compreensão de lista\n",
    "x_compr, y_compr = tuple(x for x in amostra), tuple(amostra[x] for x in amostra)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Cada uma tem suas vantagens e desvantagens, mas o ponto é que existem muitas formas de se fazer algo facilmente em python."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Obtendo as médias amostrais"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Novamente, existem alguns modos de se fazer isso em uma linha:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2017-07-19T14:20:55.953695Z",
     "start_time": "2017-07-19T14:20:55.942468Z"
    }
   },
   "outputs": [],
   "source": [
    "# modo 1: usando zip, amostra.items, sum, len e compreensão de tuplas\n",
    "x_tup, y_tup = (sum(a)/len(a) for a in zip(*amostra.items()))\n",
    "\n",
    "# modo 2: usando amostra.keys, amostra.values, map e uma função lambda\n",
    "x_map, y_map = map(lambda a: sum(a)/len(a), (amostra.keys(), amostra.values()))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Modo honesto (com as ferramentas que mostrei em aula)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Sabendo atribuir variáveis, aritmética e os `for` loops, é possível calcular os coeficientes lembrando como iterar sobre um dicionário:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2017-07-19T14:20:59.320572Z",
     "start_time": "2017-07-19T14:20:59.265885Z"
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\n",
      "+---+--------+--------+\n",
      "|   | Real   | Ajuste |\n",
      "+---+--------+--------+\n",
      "| a |   1.75 |   1.66 |\n",
      "| b |    3.5 |   3.44 |\n",
      "+---+--------+--------+\n",
      "\n"
     ]
    }
   ],
   "source": [
    "x, y = [], []\n",
    "x_media, y_media = 0, 0\n",
    "n = 0\n",
    "for (x_i, y_i) in amostra.items():\n",
    "    x.append(x_i)    \n",
    "    y.append(y_i)\n",
    "    x_media += x_i\n",
    "    y_media += y_i\n",
    "    n += 1\n",
    "x_media /= n\n",
    "y_media /= n\n",
    "num, den = 0, 0\n",
    "for i in range(n):\n",
    "    num += (x_media - x[i])*(y_media - y[i])\n",
    "    den += (x_media - x[i])**2\n",
    "b_chapeu = num/den\n",
    "a_chapeu = y_media - b_chapeu*x_media\n",
    "\n",
    "print(f\"\"\"\n",
    "+---+--------+--------+\n",
    "|   | Real   | Ajuste |\n",
    "+---+--------+--------+\n",
    "| a | {a_real:6.3} | {a_chapeu:6.3} |\n",
    "| b | {b_real:6.3} | {b_chapeu:6.3} |\n",
    "+---+--------+--------+\n",
    "\"\"\")"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.6.4"
  },
  "toc": {
   "colors": {
    "hover_highlight": "#DAA520",
    "running_highlight": "#FF0000",
    "selected_highlight": "#FFD700"
   },
   "moveMenuLeft": true,
   "nav_menu": {
    "height": "48px",
    "width": "252px"
   },
   "navigate_menu": true,
   "number_sections": true,
   "sideBar": true,
   "threshold": 4,
   "toc_cell": false,
   "toc_section_display": "block",
   "toc_window_display": false
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}