{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Analyse du risque de défaillance des joints toriques de la navette Challenger"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Le 27 Janvier 1986, veille du décollage de la navette *Challenger*, eu\n",
    "lieu une télé-conférence de trois heures entre les ingénieurs de la\n",
    "Morton Thiokol (constructeur d'un des moteurs) et de la NASA. La\n",
    "discussion portait principalement sur les conséquences de la\n",
    "température prévue au moment du décollage de 31°F (juste en dessous de\n",
    "0°C) sur le succès du vol et en particulier sur la performance des\n",
    "joints toriques utilisés dans les moteurs. En effet, aucun test\n",
    "n'avait été effectué à cette température.\n",
    "\n",
    "L'étude qui suit reprend donc une partie des analyses effectuées cette\n",
    "nuit là et dont l'objectif était d'évaluer l'influence potentielle de\n",
    "la température et de la pression à laquelle sont soumis les joints\n",
    "toriques sur leur probabilité de dysfonctionnement. Pour cela, nous\n",
    "disposons des résultats des expériences réalisées par les ingénieurs\n",
    "de la NASA durant les 6 années précédant le lancement de la navette\n",
    "Challenger.\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Chargement des données\n",
    "Nous commençons donc par charger ces données:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>Date</th>\n",
       "      <th>Count</th>\n",
       "      <th>Temperature</th>\n",
       "      <th>Pressure</th>\n",
       "      <th>Malfunction</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>4/12/81</td>\n",
       "      <td>6</td>\n",
       "      <td>66</td>\n",
       "      <td>50</td>\n",
       "      <td>0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>11/12/81</td>\n",
       "      <td>6</td>\n",
       "      <td>70</td>\n",
       "      <td>50</td>\n",
       "      <td>1</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>3/22/82</td>\n",
       "      <td>6</td>\n",
       "      <td>69</td>\n",
       "      <td>50</td>\n",
       "      <td>0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td>11/11/82</td>\n",
       "      <td>6</td>\n",
       "      <td>68</td>\n",
       "      <td>50</td>\n",
       "      <td>0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4</th>\n",
       "      <td>4/04/83</td>\n",
       "      <td>6</td>\n",
       "      <td>67</td>\n",
       "      <td>50</td>\n",
       "      <td>0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>5</th>\n",
       "      <td>6/18/82</td>\n",
       "      <td>6</td>\n",
       "      <td>72</td>\n",
       "      <td>50</td>\n",
       "      <td>0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>6</th>\n",
       "      <td>8/30/83</td>\n",
       "      <td>6</td>\n",
       "      <td>73</td>\n",
       "      <td>100</td>\n",
       "      <td>0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>7</th>\n",
       "      <td>11/28/83</td>\n",
       "      <td>6</td>\n",
       "      <td>70</td>\n",
       "      <td>100</td>\n",
       "      <td>0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>8</th>\n",
       "      <td>2/03/84</td>\n",
       "      <td>6</td>\n",
       "      <td>57</td>\n",
       "      <td>200</td>\n",
       "      <td>1</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>9</th>\n",
       "      <td>4/06/84</td>\n",
       "      <td>6</td>\n",
       "      <td>63</td>\n",
       "      <td>200</td>\n",
       "      <td>1</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>10</th>\n",
       "      <td>8/30/84</td>\n",
       "      <td>6</td>\n",
       "      <td>70</td>\n",
       "      <td>200</td>\n",
       "      <td>1</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>11</th>\n",
       "      <td>10/05/84</td>\n",
       "      <td>6</td>\n",
       "      <td>78</td>\n",
       "      <td>200</td>\n",
       "      <td>0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>12</th>\n",
       "      <td>11/08/84</td>\n",
       "      <td>6</td>\n",
       "      <td>67</td>\n",
       "      <td>200</td>\n",
       "      <td>0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>13</th>\n",
       "      <td>1/24/85</td>\n",
       "      <td>6</td>\n",
       "      <td>53</td>\n",
       "      <td>200</td>\n",
       "      <td>2</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>14</th>\n",
       "      <td>4/12/85</td>\n",
       "      <td>6</td>\n",
       "      <td>67</td>\n",
       "      <td>200</td>\n",
       "      <td>0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>15</th>\n",
       "      <td>4/29/85</td>\n",
       "      <td>6</td>\n",
       "      <td>75</td>\n",
       "      <td>200</td>\n",
       "      <td>0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>16</th>\n",
       "      <td>6/17/85</td>\n",
       "      <td>6</td>\n",
       "      <td>70</td>\n",
       "      <td>200</td>\n",
       "      <td>0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>17</th>\n",
       "      <td>7/29/85</td>\n",
       "      <td>6</td>\n",
       "      <td>81</td>\n",
       "      <td>200</td>\n",
       "      <td>0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>18</th>\n",
       "      <td>8/27/85</td>\n",
       "      <td>6</td>\n",
       "      <td>76</td>\n",
       "      <td>200</td>\n",
       "      <td>0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>19</th>\n",
       "      <td>10/03/85</td>\n",
       "      <td>6</td>\n",
       "      <td>79</td>\n",
       "      <td>200</td>\n",
       "      <td>0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>20</th>\n",
       "      <td>10/30/85</td>\n",
       "      <td>6</td>\n",
       "      <td>75</td>\n",
       "      <td>200</td>\n",
       "      <td>2</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>21</th>\n",
       "      <td>11/26/85</td>\n",
       "      <td>6</td>\n",
       "      <td>76</td>\n",
       "      <td>200</td>\n",
       "      <td>0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>22</th>\n",
       "      <td>1/12/86</td>\n",
       "      <td>6</td>\n",
       "      <td>58</td>\n",
       "      <td>200</td>\n",
       "      <td>1</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "        Date  Count  Temperature  Pressure  Malfunction\n",
       "0    4/12/81      6           66        50            0\n",
       "1   11/12/81      6           70        50            1\n",
       "2    3/22/82      6           69        50            0\n",
       "3   11/11/82      6           68        50            0\n",
       "4    4/04/83      6           67        50            0\n",
       "5    6/18/82      6           72        50            0\n",
       "6    8/30/83      6           73       100            0\n",
       "7   11/28/83      6           70       100            0\n",
       "8    2/03/84      6           57       200            1\n",
       "9    4/06/84      6           63       200            1\n",
       "10   8/30/84      6           70       200            1\n",
       "11  10/05/84      6           78       200            0\n",
       "12  11/08/84      6           67       200            0\n",
       "13   1/24/85      6           53       200            2\n",
       "14   4/12/85      6           67       200            0\n",
       "15   4/29/85      6           75       200            0\n",
       "16   6/17/85      6           70       200            0\n",
       "17   7/29/85      6           81       200            0\n",
       "18   8/27/85      6           76       200            0\n",
       "19  10/03/85      6           79       200            0\n",
       "20  10/30/85      6           75       200            2\n",
       "21  11/26/85      6           76       200            0\n",
       "22   1/12/86      6           58       200            1"
      ]
     },
     "execution_count": 1,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "import numpy as np\n",
    "import pandas as pd\n",
    "data = pd.read_csv(\"shuttle.csv\")\n",
    "data"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Le jeu de données nous indique la date de l'essai, le nombre de joints\n",
    "toriques mesurés (il y en a 6 sur le lançeur principal), la\n",
    "température (en Farenheit) et la pression (en psi), et enfin le\n",
    "nombre de dysfonctionnements relevés. "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Inspection graphique des données\n",
    "Les vols où aucun incident n'est relevé n'apportant aucun information\n",
    "sur l'influence de la température ou de la pression sur les\n",
    "dysfonctionnements, nous nous concentrons sur les expériences où au\n",
    "moins un joint a été défectueux.\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "#data = data[data.Malfunction > 0]\n",
    "#data"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Très bien, nous avons une variabilité de température importante mais\n",
    "la pression est quasiment toujours égale à 200, ce qui devrait\n",
    "simplifier l'analyse.\n",
    "\n",
    "Comment la fréquence d'échecs varie-t-elle avec la température ?\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "/opt/conda/lib/python3.6/site-packages/pandas/plotting/_core.py:203: UserWarning: 'color' and 'colormap' cannot be used simultaneously. Using 'color'\n",
      "  warnings.warn(\"'color' and 'colormap' cannot be used \"\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYwAAAEKCAYAAAAB0GKPAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMi4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvIxREBQAAHNBJREFUeJzt3X2QXXWd5/H3pzvPDw5sIq2TkEk0KTUrIYY2wIBOZxUW3BoybsYBV8XSYSKrGQvWXaFcl2UcZ0pYH3aoQUNksYQtJqNiNLtGA6gtg4okODEBhLE3QdJEEtIGSUMeutPf/eOe1pvO7e7f7b7nPvXnVdXV5/zO73f6+72nb3/7PNxzFBGYmZmNpqXWAZiZWWNwwTAzsyQuGGZmlsQFw8zMkrhgmJlZEhcMMzNL4oJhZmZJXDDMzCyJC4aZmSWZVOsAKmnu3LmxcOHCWodxkhdffJGZM2fWOoyKc16Np1lzc17j88gjjxyMiJen9G2qgrFw4UK2b99e6zBO0tnZSUdHR63DqDjn1XiaNTfnNT6Sfpna14ekzMwsSa4FQ9Ilkp6U1CXp+hLLV0vaKWmHpO2SLixa9pSkXYPL8ozTzMxGl9shKUmtwK3ARUA3sE3S5oh4vKjbd4HNERGSlgFfAV5btHxVRBzMK0YzM0uX5x7GSqArInZHxHFgI7C6uENE9Mbv7q8+E/C91s3M6lSeBWMesLdovjtrO4mkt0t6AvgW8P6iRQHcK+kRSWtzjNPMzBIorwcoSXoH8G8j4qps/j3Ayoj4y2H6vxm4ISLems3/fkTsk3QGcB/wlxHxQIlxa4G1AG1tbeds3Lgxl3zGqre3l1mzZtU6jIpzXo2nWXNzXuOzatWqRyKiPaVvnpfVdgNnFs3PB/YN1zkiHpD0aklzI+JgROzL2g9I2kThENcpBSMiNgAbANrb26PeLq/zJX+NpVnzgubNzXlVT56HpLYBSyQtkjQFuALYXNxB0mJJyqZXAFOAHkkzJc3O2mcCFwOP5hirmZmNIrc9jIjol7QO2Aq0AndExGOSrs6WrwfWAFdK6gOOAJdnV0y1AZuyWjIJuDsivpNXrGZmNrpcP+kdEVuALUPa1hdN3wTcVGLcbuDsPGMzM7Py+JPeZmaWxAUD6Ok9xs/2Pk9P77Fah2JmZejpPcaRvhN+71bJhC8Y39zxDBfc9D3efftPuOCm77F5xzO1DsnMEgy+d/c896Lfu1UyoQtGT+8xrrtnJ0f7Bjh8rJ+jfQN89J6d/m/FrM4Vv3dPRPi9WyUTumB0HzrC5JaTX4LJLS10HzpSo4jMLIXfu7UxoQvG/NOn0zcwcFJb38AA80+fXqOIzCyF37u1MaELxpxZU7l5zTKmTW5h9tRJTJvcws1rljFn1tRah2ZmIyh+77ZKfu9WSVM9cW8sLls+jwsWz6X70BHmnz7dv3BmDWLwvfvwjx/kh5dd6PduFUz4ggGF/1b8y2bWeObMmsr0ya1+/1bJhD4kZWZm6VwwzMwsiQuGmZklccEwM7MkLhhmZpbEBcPMzJK4YJiZWRIXDDMzS+KCYWZmSVwwzMwsSa4FQ9Ilkp6U1CXp+hLLV0vaKWmHpO2SLkwda2Zm1ZVbwZDUCtwKXAosBd4paemQbt8Fzo6I5cD7gdvLGGtmZlWU5x7GSqArInZHxHFgI7C6uENE9EZEZLMzgUgda2Zm1ZVnwZgH7C2a787aTiLp7ZKeAL5FYS8jeayZmVVPnrc3V4m2OKUhYhOwSdKbgb8G3po6FkDSWmAtQFtbG52dnWONNxe9vb11F1MlOK/G06y5Oa/qybNgdANnFs3PB/YN1zkiHpD0aklzyxkbERuADQDt7e3R0dExzrArq7Ozk3qLqRKcV+Np1tycV/XkeUhqG7BE0iJJU4ArgM3FHSQtlqRsegUwBehJGWtmZtWV2x5GRPRLWgdsBVqBOyLiMUlXZ8vXA2uAKyX1AUeAy7OT4CXH5hWrmZmNLtdHtEbEFmDLkLb1RdM3ATeljjUzs9rxJ73NzCyJC4aZmSVxwTAzsyQuGGZmlsQFw8zMkrhgmJlZEhcMMzNL4oJhZmZJXDDMzCyJC4aZmSVxwTAzsyQuGGZmlsQFw8zMkrhgmJlZEhcMMzNL4oJhZmZJXDDMzCyJC4aZmSVxwTAzsyS5FgxJl0h6UlKXpOtLLH+XpJ3Z148knV207ClJuyTtkLQ9zzjNzGx0k/JasaRW4FbgIqAb2CZpc0Q8XtRtD/BHEXFI0qXABuDcouWrIuJgXjGamVm6PPcwVgJdEbE7Io4DG4HVxR0i4kcRcSibfQiYn2M8ZmY2DnkWjHnA3qL57qxtOH8OfLtoPoB7JT0iaW0O8ZmZWRlyOyQFqERblOworaJQMC4sar4gIvZJOgO4T9ITEfFAibFrgbUAbW1tdHZ2jjvwSurt7a27mCrBeTWeZs3NeVVPngWjGzizaH4+sG9oJ0nLgNuBSyOiZ7A9IvZl3w9I2kThENcpBSMiNlA490F7e3t0dHRUMIXx6+zspN5iqgTn1XiaNTfnVT15HpLaBiyRtEjSFOAKYHNxB0kLgK8D74mIfylqnylp9uA0cDHwaI6xmpnZKHLbw4iIfknrgK1AK3BHRDwm6eps+XrgBmAO8HlJAP0R0Q60AZuytknA3RHxnbxiNTOz0eV5SIqI2AJsGdK2vmj6KuCqEuN2A2cPbTczs9rxJ73NzCyJC4aZmSVxwTAzsyQuGGZmlsQFw8zMkrhgmJlZEhcMMzNL4oJhZmZJXDDMzCyJC4aZmSVxwTAzsyQuGGZmlsQFw8zMkrhgmJlZEhcMMzNL4oJhZmZJXDDMzCxJUsGQ9Pq8AzEzs/qWuoexXtLDkj4o6bRcIzIzs7qUVDAi4kLgXcCZwHZJd0u6aLRxki6R9KSkLknXl1j+Lkk7s68fSTo7dayZmVVX8jmMiPgF8HHgOuCPgFskPSHp35fqL6kVuBW4FFgKvFPS0iHd9gB/FBHLgL8GNpQx1szMqij1HMYySZ8Dfg78G+CPI+J12fTnhhm2EuiKiN0RcRzYCKwu7hARP4qIQ9nsQ8D81LFmZlZdkxL7/T3wReBjEXFksDEi9kn6+DBj5gF7i+a7gXNH+Bl/Dny73LGS1gJrAdra2ujs7BzhR1Rfb29v3cVUCc6r8TRrbs6relILxtuAIxFxAkBSCzAtIl6KiLuGGaMSbVGyo7SKQsG4sNyxEbGB7FBWe3t7dHR0DJdDTXR2dlJvMVWC82o8zZqb86qe1HMY9wPTi+ZnZG0j6aZwknzQfGDf0E6SlgG3A6sjoqecsWZmVj2pBWNaRPQOzmTTM0YZsw1YImmRpCnAFcDm4g6SFgBfB94TEf9SzlgzM6uu1ENSL0paERE/BZB0DnBkpAER0S9pHbAVaAXuiIjHJF2dLV8P3ADMAT4vCaA/ItqHGzuG/MzMrEJSC8Y1wFclDR4WeiVw+WiDImILsGVI2/qi6auAq1LHmplZ7SQVjIjYJum1wGsonJB+IiL6co3MzMzqSuoeBsAbgYXZmDdIIiLuzCUqMzOrO0kFQ9JdwKuBHcCJrDkAFwwzswkidQ+jHVgaESU/C2FmZs0v9bLaR4FX5BmImZnVt9Q9jLnA45IeBo4NNkbEZblEZWZmdSe1YNyYZxBmZlb/Ui+r/YGkPwCWRMT9kmZQ+ECdmZlNEKm3N/8L4GvAbVnTPOAbeQVlZmb1J/Wk94eAC4AX4LcPUzojr6DMzKz+pBaMY9mDjACQNIlhbjdu+erpPcbP9j5PT++xivZtNM2cWx669h/m0Et9dO0/XOtQrIGlnvT+gaSPAdOzZ3l/EPg/+YVlpXxzxzNcd89OJre00DcwwM1rlnHZ8nnj7ttomjm3PNzwjV3c+dDTfOSsfq793ANcef4CPrH6rFqHZQ0odQ/jeuA5YBfwAQo3BRzuSXuWg57eY1x3z06O9g1w+Fg/R/sG+Og9O0v+h11O30bTzLnloWv/Ye586OmT2u788dPe07AxSSoYETEQEV+MiHdExJ9m0z4kVUXdh44wueXkzTW5pYXuQ6feZb6cvo2mmXPLw469z5fVbjaS1HtJ7aHEOYuIeFXFI7KS5p8+nb6BgZPa+gYGmH/69HH1bTTNnFselp95WlntZiNJPSTVTuFutW8E3gTcAvzvvIKyU82ZNZWb1yxj2uQWZk+dxLTJLdy8ZhlzZk0dV99G08y55WFx22yuPH/BSW1Xnr+AxW2zaxSRNbLUD+71DGn6n5IepPDEPKuSy5bP44LFc+k+dIT5p08f8Y9kOX0bTTPnlodPrD6LK89byK5HHuL+a89zsbAxSz0ktaJotoXCHod/62pgzqypyX8gy+nbaJo5tzwsbptN94zJLhY2LqmX1X6maLofeAr4s4pHY2ZmdSv1kNSqsaxc0iXA31G479TtEfGpIctfC3wJWAH814j4dNGyp4DDFB7Y1B8R7WOJwczMKiP1kNR/Gml5RHy2xJhW4FbgIqAb2CZpc0Q8XtTt18CHgT8ZZtWrIuJgSoxmZpavcq6S+o8Ubjo4D7gaWErhPMZwB0VXAl0RsTu7rchGYHVxh4g4EBHbgL4xxG5mZlVUzgOUVkTEYQBJNwJfjYirRhgzD9hbNN8NnFtGbAHcKymA2yJiQ6lOktYCawHa2tro7Ows40fkr7e3t+5iqgTn1XiaNTfnVT2pBWMBcLxo/jiwcJQxKtFWzqfDL4iIfZLOAO6T9EREPHDKCguFZANAe3t7dHR0lPEj8tfZ2Um9xVQJzqvxNGtuzqt6UgvGXcDDkjZR+KP/duDOUcZ0A2cWzc8H9qUGFhH7su8Hsp+7EjilYJiZWXWk3kvqb4D3AYeA54H3RcTfjjJsG7BE0iJJU4ArgM0pP0/STEmzB6eBi4FHU8aamVk+UvcwAGYAL0TElyS9XNKiiNgzXOeI6Je0DthK4bLaOyLiMUlXZ8vXS3oFsB14GTAg6RoKJ9PnApskDcZ4d0R8ZywJmplZZaReVvvfKVwp9RoKn5uYTOFeUheMNC4itlC4FXpx2/qi6WcpHKoa6gXg7JTYzMysOlIvq307cBnwIvz2/ILvMWBmNoGkFozj2fMvAn57XsHMzCaQ1ILxFUm3AadJ+gvgfuCL+YVlZmb1JvVeUp/OnuX9AoXzGDdExH25RmZmZnVl1IKR3RNqa0S8FXCRMDOboEY9JBURJ4CXJP1eFeIxM7M6lfo5jKPALkn3kV0pBRARH84lKjMzqzupBeNb2ZeZmU1QIxYMSQsi4umI+HK1AjIzs/o02jmMbwxOSLon51jMzKyOjVYwim9R/qo8AzEzs/o2WsGIYabNzGyCGe2k99mSXqCwpzE9myabj4h4Wa7RmZlZ3RixYEREa7UCMTOz+pZ6LykzM5vgXDDMzCyJC4aZmSVxwTAzsyQuGGZmliTXgiHpEklPSuqSdH2J5a+V9GNJxyT953LGmplZdeVWMLLnaNwKXAosBd4paemQbr8GPgx8egxjzcysivLcw1gJdEXE7og4DmwEVhd3iIgDEbEN6Ct3rJmZVVfq7c3HYh6wt2i+Gzi30mMlrQXWArS1tdHZ2Vl2oHnq7e2tu5gqwXk1nmbNzXlVT54FQyXaUu9HlTw2IjYAGwDa29ujo6Mj8UdUR2dnJ/UWUyU4r8bTrLk5r+rJ85BUN3Bm0fx8YF8VxpqZWQ7yLBjbgCWSFkmaAlwBbK7CWDMzy0Fuh6Qiol/SOmAr0ArcERGPSbo6W75e0iuA7cDLgAFJ1wBLI+KFUmPzitXMzEaX5zkMImILsGVI2/qi6WcpHG5KGmtmZrXjT3qbmVkSFwwzM0vigmFmZklcMMzMLIkLhpmZJXHBMDOzJC4YZmaWxAXDzMySuGCYmVkSFwwzM0vigmFmZklcMMzMLIkLhpmZJXHBMDOzJC4YZmaWxAXDzMySuGCYmVkSFwwzM0vigmFmZklyLRiSLpH0pKQuSdeXWC5Jt2TLd0paUbTsKUm7JO2QtD3POM3MbHST8lqxpFbgVuAioBvYJmlzRDxe1O1SYEn2dS7whez7oFURcTCvGM3MLF2eexgrga6I2B0Rx4GNwOohfVYDd0bBQ8Bpkl6ZY0xmZjZGue1hAPOAvUXz3Zy89zBcn3nAr4AA7pUUwG0RsaHUD5G0FlgL0NbWRmdnZ0WCr5Te3t66i6kSnFfjadbcnFf15FkwVKItyuhzQUTsk3QGcJ+kJyLigVM6FwrJBoD29vbo6OgYR8iV19nZSb3FVAnOq/E0a27Oq3ryPCTVDZxZND8f2JfaJyIGvx8ANlE4xGVmZjWSZ8HYBiyRtEjSFOAKYPOQPpuBK7Orpc4DfhMRv5I0U9JsAEkzgYuBR3OM1czMRpHbIamI6Je0DtgKtAJ3RMRjkq7Olq8HtgBvA7qAl4D3ZcPbgE2SBmO8OyK+k1esZmY2ujzPYRARWygUheK29UXTAXyoxLjdwNl5xmZmZuXxJ73NzCyJC4aZmSVxwTAzsyQuGGZmlsQFw8zMkrhgmJlZEhcMMzNL4oJhZmZJXDDMzCyJC4aZmSVxwTAzsyQuGGZmlsQFw8zMkrhgmJlZEhcMMzNL4oJhZmZJXDDMzCyJC4aZmSVxwTAzsyS5FgxJl0h6UlKXpOtLLJekW7LlOyWtSB1rVks9vcf42d7n6ek9Nmrf7Xt6+Oy9T7J9T09F11tO3679hzn0Uh9d+w+P2jdVXrGWG8ORvhPJr8HXtu9tytcgj/WWMimvFUtqBW4FLgK6gW2SNkfE40XdLgWWZF/nAl8Azk0ca1YT39zxDNfds5PJLS30DQxw85plXLZ8Xsm+7779IR7sKhSKW77XxZsWz+Guq84b93rL6XvDN3Zx50NP85Gz+rn2cw9w5fkL+MTqs8aQef6xjiWGD7+uj2tv+l7SazCo2V6DSq93OHnuYawEuiJid0QcBzYCq4f0WQ3cGQUPAadJemXiWLOq6+k9xnX37ORo3wCHj/VztG+Aj96zs+R/d9v39Py2WAz6p66eknsa5ay3nL5d+w+f9IcS4M4fPz2u/7LzinWsMZyImPCvQSXXOxJFRD4rlv4UuCQirsrm3wOcGxHrivr8X+BTEfFgNv9d4Dpg4Whji9axFlgL0NbWds7GjRtzyWesent7mTVrVq3DqLiJmteRvhPsee5FThS9b1olFr18JtMnt57Ud/8Lxzhw+Ogp6zhj9jTaXjZ1zOstp++hl/roPvQSAG3TYf+RQvv802dw+ozJw+Y5krxiHWsMg3mlvAbF6v01qOTv4khWrVr1SES0p/TN7ZAUoBJtQ6vTcH1SxhYaIzYAGwDa29ujo6OjjBDz19nZSb3FVAkTNa+e3mNce9P3ONo38Nu2aZNb+OFlFzJn1slFYPueHq677aFT1vG1D7TTvmjOmNdbTt+u/Ye59nMPAPCRs/r5zK7CW/7+a89jcdvsYfMcSV6xjjWGwbxSXoNi9f4aVPJ3sVLyPCTVDZxZND8f2JfYJ2WsWdXNmTWVm9csY9rkFmZPLfyRunnNspJv0PZFc3jT4pMLw5sWzzmlWJS73nL6Lm6bzZXnLzip7crzF4z5D2WesY41hlZpwr8GlVzviCIily8Key+7gUXAFOBnwL8e0uffAd+msEdxHvBw6thSX+ecc07Um+9///u1DiEXEz2vg4ePxo6nD8XBw0dH7btt98H4zNYnYtvugxVdbzl9f/HsC/H1b90bv3j2hVH7psor1nJj2HLv/cmvwVe3Pd0wr0Eev4ulANsj8e96boekIqJf0jpgK9AK3BERj0m6Olu+HtgCvA3oAl4C3jfS2LxiNSvXnFlTk/+Ta19Ueq9ivOstp+/ittl0z5g8rv+qx/Pzy+lbbgzTJ7cmrXtx2+yK5j/48+vhNch1r6JInucwiIgtFIpCcdv6oukAPpQ61szMasef9DYzsyQuGGZmlsQFw8zMkrhgmJlZEhcMMzNL4oJhZmZJXDDMzCyJC4aZmSVxwTAzsyQuGGZmliS352HUgqTngF/WOo4h5gIHax1EDpxX42nW3JzX+PxBRLw8pWNTFYx6JGl7JD6cpJE4r8bTrLk5r+rxISkzM0vigmFmZklcMPK3odYB5MR5NZ5mzc15VYnPYZiZWRLvYZiZWRIXjAqS9JSkXZJ2SNqetd0o6ZmsbYekt9U6znJJOk3S1yQ9Iennks6X9K8k3SfpF9n302sd51gMk1tDbzNJrymKfYekFyRd0+jbbIS8Gnp7AUi6VtJjkh6V9A+SptXj9vIhqQqS9BTQHhEHi9puBHoj4tO1imu8JH0Z+KeIuF3SFGAG8DHg1xHxKUnXA6dHxHU1DXQMhsntGhp8mw2S1Ao8A5xL4XHIDb/N4JS83kcDby9J84AHgaURcUTSVyg8nnopdba9vIdhI5L0MuDNwP8CiIjjEfE8sBr4ctbty8Cf1CbCsRsht2byFuD/RcQvaYJtVqQ4r2YwCZguaRKFf1r2UYfbywWjsgK4V9IjktYWta+TtFPSHfWwW1mmVwHPAV+S9M+Sbpc0E2iLiF8BZN/PqGWQYzRcbtDY26zYFcA/ZNPNsM0GFecFDby9IuIZ4NPA08CvgN9ExL3U4fZywaisCyJiBXAp8CFJbwa+ALwaWE7hl+EzNYxvLCYBK4AvRMQbgBeB62sbUsUMl1ujbzMAskNslwFfrXUslVQir4beXlmBWw0sAn4fmCnp3bWNqjQXjAqKiH3Z9wPAJmBlROyPiBMRMQB8EVhZyxjHoBvojoifZPNfo/BHdr+kVwJk3w/UKL7xKJlbE2yzQZcCP42I/dl8M2wzGJJXE2yvtwJ7IuK5iOgDvg78IXW4vVwwKkTSTEmzB6eBi4FHBzd45u3Ao7WIb6wi4llgr6TXZE1vAR4HNgPvzdreC3yzBuGNy3C5Nfo2K/JOTj5s0/DbLHNSXk2wvZ4GzpM0Q5Io/B7+nDrcXr5KqkIkvYrCXgUUDnXcHRF/I+kuCrvKATwFfGDwuGSjkLQcuB2YAuymcFVKC/AVYAGFX/h3RMSvaxbkGA2T2y00/jabAewFXhURv8na5tDg22yYvJrhPfZXwOVAP/DPwFXALOpse7lgmJlZEh+SMjOzJC4YZmaWxAXDzMySuGCYmVkSFwwzM0syqdYBmOUtu5z0u9nsK4ATFG4JAoUPVx6vSWAjkPR+YEv2WRGzuuDLam1Cqae7B0tqjYgTwyx7EFgXETvKWN+kiOivWIBmQ/iQlE1okt4r6eHsOQqfl9QiaZKk5yX9D0k/lbRV0rmSfiBp9+DzFiRdJWlTtvxJSR9PXO8nJT0MrJT0V5K2Zc9BWK+Cyyl8EO0fs/FTJHVLOi1b93mS7s+mPynpNkn3UbiJ4iRJn81+9k5JV1X/VbVm5YJhE5ak11O4lcQfRsRyCodor8gW/x5wb3YzyePAjRRu2fAO4BNFq1mZjVkB/AdJyxPW+9OIWBkRPwb+LiLeCJyVLbskIv4R2AFcHhHLEw6ZvQH444h4D7AWOBARK4E3UrgJ5oKxvD5mQ/kchk1kb6XwR3V74RY+TKdw2wmAIxFxXza9i8Itp/sl7QIWFq1ja0QcApD0DeBCCu+r4dZ7nN/dQgbgLZL+CzANmAs8Any7zDy+GRFHs+mLgddJKi5QSyjcWsJsXFwwbCITcEdE/LeTGgsPsSn+r34AOFY0Xfy+GXoSMEZZ75HIThxm90X6ewp3yH1G0icpFI5S+vndEYGhfV4cktMHI+K7mFWYD0nZRHY/8GeS5kLhaqoxHL65WIXngs+g8EyDH5ax3ukUCtDB7E7Ha4qWHQZmF80/BZyTTRf3G2or8MGsOA0+B3t6mTmZleQ9DJuwImJXdpfQ+yW1AH3A1RQej5nqQeBuCg/wuWvwqqaU9UZEjwrPFH8U+CXwk6LFXwJul3SEwnmSG4EvSnoWeHiEeG6jcHfTHdnhsAMUCpnZuPmyWrMxyq5Aen1EXFPrWMyqwYekzMwsifcwzMwsifcwzMwsiQuGmZklccEwM7MkLhhmZpbEBcPMzJK4YJiZWZL/D+TWzOgkFsXDAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYwAAAEKCAYAAAAB0GKPAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMi4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvIxREBQAAHOxJREFUeJzt3X+8VXWd7/HXm18KHEyTPJpgUJBexh8EJ5SrU+fkj7CZpK422e1qv4hLE3Z1xjvi2INpmmsPrWbsNmMRld1sps7DfmiYTIg8PDGZGOAQvxQ7gcoJwR9RcAAF5HP/WIvcHM+P7+GctfdC38/HYz/OWt/1/e713pyz94f1Y6+liMDMzKwnA2odwMzMjgwuGGZmlsQFw8zMkrhgmJlZEhcMMzNL4oJhZmZJXDDMzCyJC4aZmSVxwTAzsySDah2gP40cOTLGjBlT6xgA7Nq1i+HDh9c6RrfKnrHs+aD8GcueD5yxP/Ql38qVK5+NiNcldY6IV8xj8uTJURb3339/rSP0qOwZy54vovwZy54vwhn7Q1/yASsi8TPWu6TMzCxJoQVD0jRJGyS1SprTyfLpklZLWiVphaTzKpY9LmnNwWVF5jQzs54VdgxD0kDgVuBCoA1YLmlBRKyv6LYEWBARIelM4A7gtIrlTRHxbFEZzcwsXZFbGFOA1ojYGBF7gWZgemWHiGjP96EBDAd8rXUzs5IqsmCcDGyumG/L2w4h6b2SHgXuAT5asSiAeyWtlDSzwJxmZpZAL/0Hv5+fWHof8M6ImJHPXwFMiYiruuj/NmBuRFyQz78+IrZIOgFYDFwVEUs7GTcTmAlQX18/ubm5uZDX01vt7e3U1dXVOka3yp6x7Pmg/BnLng+csT/0JV9TU9PKiGhI6px6OlVvH8BUYFHF/PXA9T2M2QSM7KT9M8C1Pa3Tp9X2Ttkzlj1fRPkzlj1fhDP2h1fCabXLgfGSxkoaAlwOLKjsIGmcJOXTk4AhwHOShksakbcPBy4C1haY1czMelDYWVIRsV/SbGARMBC4LSLWSZqVL58HXApcKWkfsAd4f0SEpHrgzryWDAK+GxE/LSqrmZn1rNBLg0TEQmBhh7Z5FdM3Azd3Mm4jcFaR2czMrHf8TW8zM0vigmFmdgRr3baT7bv30bptZ+HrcsEwMztCzb1rDRfcspS27bu54JalzP3xmkLX54JhZnYEat22k9uXPXlI2+0PPlnoloYLhpnZEWjV5t/3qr0/uGCYmR2BJo4+tlft/cEFw8zsCDSufgRXTj3lkLYrp57CuPoRha3zFXWLVjOzV5PPTj+DK88Zw5qVy7jvmnMKLRbgLQwzsyPauPoRHDdscOHFAlwwzMwskQuGmZklccEwM7MkLhhmZpbEBcPMzJK4YJiZWRIXDDMzS+KCYWZmSVwwzMwsiQuGmZklKbRgSJomaYOkVklzOlk+XdJqSaskrZB0XupYMzOrrsIKhqSBwK3AxcAE4AOSJnTotgQ4KyImAh8FvtGLsWZmVkVFbmFMAVojYmNE7AWagemVHSKiPSIinx0OROpYMzOrriILxsnA5or5trztEJLeK+lR4B6yrYzksWZmVj1F3g9DnbTFyxoi7gTulPQ24B+AC1LHAkiaCcwEqK+vp6Wl5XDz9qv29vbSZOlK2TOWPR+UP2PZ84Ez9odq5SuyYLQBoyvmRwFbuuocEUslvUnSyN6MjYj5wHyAhoaGaGxs7GPs/tHS0kJZsnSl7BnLng/Kn7Hs+cAZ+0O18hW5S2o5MF7SWElDgMuBBZUdJI2TpHx6EjAEeC5lrJmZVVdhWxgRsV/SbGARMBC4LSLWSZqVL58HXApcKWkfsAd4f34QvNOxRWU1M7OeFXpP74hYCCzs0DavYvpm4ObUsWZmVjv+preZmSVxwTAzsyQuGGZmlsQFw8zMkrhgmJlZEhcMMzNL4oJhZmZJXDDMzCyJC4aZmSVxwTAzsyQuGGZmlsQFw8zMkrhgmJlZEhcMMzNL4oJhZmZJXDDMzCyJC4aZmSVxwTAzsyQuGGZmlqTQgiFpmqQNklolzelk+Qclrc4fv5B0VsWyxyWtkbRK0ooic5qZWc8GFfXEkgYCtwIXAm3AckkLImJ9RbdNwNsjYruki4H5wNkVy5si4tmiMpqZWboitzCmAK0RsTEi9gLNwPTKDhHxi4jYns8uA0YVmMfMzPqgyIJxMrC5Yr4tb+vKx4B/r5gP4F5JKyXNLCCfmZn1giKimCeW3ge8MyJm5PNXAFMi4qpO+jYBXwHOi4jn8rbXR8QWSScAi4GrImJpJ2NnAjMB6uvrJzc3Nxfyenqrvb2durq6WsfoVtkzlj0flD9j2fOBM/aHvuRrampaGRENSZ0jopAHMBVYVDF/PXB9J/3OBH4DvLmb5/oMcG1P65w8eXKUxf3331/rCD0qe8ay54sof8ay54twxv7Ql3zAikj8XC9yl9RyYLyksZKGAJcDCyo7SDoF+BFwRUQ8VtE+XNKIg9PARcDaArOamVkPCjtLKiL2S5oNLAIGArdFxDpJs/Ll84C5wPHAVyQB7I9s06geuDNvGwR8NyJ+WlRWMzPrWWEFAyAiFgILO7TNq5ieAczoZNxG4KyO7WZmVjv+preZmSVxwTAzsyQuGGZmlsQFw8zMkrhgmJlZEhcMMzNL4oJhZmZJXDDMzCyJC4aZmSVxwTAzsyQuGGZmlsQFw8zMkrhgmJlZEhcMMzNL4oJhZmZJXDDMzCyJC4aZmSVJKhiSTi86iJmZlVvqFsY8Sb+U9JeSji00kZmZlVJSwYiI84APAqOBFZK+K+nCnsZJmiZpg6RWSXM6Wf5BSavzxy8knZU61szMqiv5GEZE/Br4NHAd8Hbgy5IelfTfOusvaSBwK3AxMAH4gKQJHbptAt4eEWcC/wDM78VYMzOrotRjGGdKugV4BHgH8O6I+C/59C1dDJsCtEbExojYCzQD0ys7RMQvImJ7PrsMGJU61szMqksR0XMnaSnwdeAHEbGnw7IrIuI7nYy5DJgWETMO9gPOjojZXazjWuC0iJjRm7GSZgIzAerr6yc3Nzf3+Hqqob29nbq6ulrH6FbZM5Y9H5Q/Y9nzgTP2h77ka2pqWhkRDSl9ByU+57uAPRHxIoCkAcDREbG7s2KRUydtnVYnSU3Ax4Dzejs2IuaT78pqaGiIxsbGrl5DVbW0tFCWLF0pe8ay54PyZyx7PnDG/lCtfKnHMO4DhlbMD8vbutNGdpD8oFHAlo6dJJ0JfAOYHhHP9WasmZlVT2rBODoi2g/O5NPDehizHBgvaaykIcDlwILKDpJOAX4EXBERj/VmrJmZVVfqLqldkiZFxMMAkiYDe7obEBH7Jc0GFgEDgdsiYp2kWfnyecBc4HjgK5IA9kdEQ1djD+P1mZlZP0ktGFcD35d0cLfQScD7exoUEQuBhR3a5lVMzwBmpI41M7PaSSoYEbFc0mnAqWQHpB+NiH2FJjMzs1JJ3cIAeCswJh/zFklExO2FpDIzs9JJKhiSvgO8CVgFvJg3B+CCYWb2KpG6hdEATIiUb/mZmdkrUupptWuBE4sMYmZm5Za6hTESWC/pl8ALBxsj4pJCUpmZWemkFozPFBnCzMzKL/W02p9JegMwPiLukzSM7At1Zmb2KpF6efOPAz8AvpY3nQzcVVQoMzMrn9SD3p8EzgV2wB9vpnRCUaHMzKx8UgvGC/mNjACQNIguLjd+JHqu/QV+tfn3PNf+Qs+dzcxKZMn6rfx2+x6WrN9a+LpSD3r/TNLfAkPze3n/JXB3cbGq58erfst1P1zN4AED2HfgAJ+/9EwumXhyrWOZmfXooltaeGzbLv76jP3ccPtKTq0fzqJrGgtbX+oWxhzgGWAN8D/JLgr46aJCVctz7S9w3Q9X8/y+A+x8YT/P7zvA3/xwtbc0zKz0lqzfymPbdh3StmHbrkK3NFLPkjpAdovWrxeWpAbatu9h8IABPM+BP7YNHjCAtu17OL7uqBomMzPr3r3rt3XZfv6EYr5nnXqW1CZJGzs+CklURaOOG8q+AwcOadt34ACjjhvaxQgzs3K4aEJ9r9r7Q+ouqQayq9W+FfhT4MvAvxYVqlqOrzuKz196JkcPHsCIowZx9OABfP7SM711YWald/6EEzm1fvghbafWDy9s6wLSd0k916HpS5J+TnbHvCPaJRNP5txxI2nbvodRxw11sTCzI8aiaxpZsn4rWzf8J9+88i2FFgtIv7z5pIrZAWRbHCMKSVQDx9cd5UJhZkek8yecSMvTQ2ksuFhA+mm1/1gxvR94HPiLfk9jZmallbpLqulwnlzSNOD/kl136hsRcVOH5acB3wImATdExBcrlj0O7CS7YdP+iGg4nAxmZtY/UndJ/VV3yyPinzoZMxC4FbgQaAOWS1oQEesruv0O+BTwni6euikink3JaGZmxerNWVKfILvo4MnALGAC2XGMro5lTAFaI2JjflmRZmB6ZYeIeDoilgP7DiO7mZlVkVLuuirpXuDSiNiZz48Avh8R07oZcxkwLSJm5PNXAGdHxOxO+n4GaO+wS2oTsJ3smlVfi4j5XaxnJjAToL6+fnJzc3OPr6ca2tvbqaurq3WMbpU9Y9nzQfkzlj0fOGN/6Eu+pqamlam7/FMPep8C7K2Y3wuM6WGMOmnrzQULz42ILZJOABZLejQilr7sCbNCMh+goaEhGhsbe7GK4rS0tFCWLF0pe8ay54PyZyx7PnDG/lCtfKkF4zvALyXdSfah/17g9h7GtAGjK+ZHAVtSg0XElvzn0/l6pwAvKxhmZlYdSccwIuJG4CNku4h+D3wkIj7Xw7DlwHhJYyUNAS4HFqSsT9LwfLcXkoYDFwFrU8aamVkxUrcwAIYBOyLiW5JeJ2lsRGzqqnNE7Jc0G1hEdlrtbRGxTtKsfPk8SScCK4BjgAOSriY7mD4SuFPSwYzfjYifHs4LNDOz/pF6Wu3fkZ0pdSrZ9yYGk11L6tzuxkXEQrJLoVe2zauY3kq2q6qjHcBZKdnMzKw6Uk+rfS9wCbAL/nh84RVzaRAzM+tZasHYG9n5twF/PK5gZmavIqkF4w5JXwOOlfRx4D5eYTdTMjOz7qVeS+qL+b28d5Adx5gbEYsLTWZmZqXSY8HIrwm1KCIuAFwkzMxepXrcJRURLwK7Jb2mCnnMzKykUr+H8TywRtJi8jOlACLiU4WkMjOz0kktGPfkDzMze5XqtmBIOiUinoyIb1crkJmZlVNPxzDuOjgh6YcFZzEzsxLrqWBUXqL8jUUGMTOzcuupYEQX02Zm9irT00HvsyTtINvSGJpPk89HRBxTaDozMyuNbgtGRAysVhAzMyu31GtJmZnZq5wLhpmZJXHBMDOzJC4YZmaWxAXDzMySFFowJE2TtEFSq6Q5nSw/TdKDkl6QdG1vxpqZWXUVVjDy+2jcClwMTAA+IGlCh26/Az4FfPEwxpqZWRUVuYUxBWiNiI0RsRdoBqZXdoiIpyNiObCvt2PNzKy6Ui9vfjhOBjZXzLcBZ/f3WEkzgZkA9fX1tLS09DpoEdrb20uTpStlz1j2fFD+jGXPB87YH6qVr8iCoU7aUq9HlTw2IuYD8wEaGhqisbExcRXFamlpoSxZulL2jGXPB+XPWPZ84Iz9oVr5itwl1QaMrpgfBWypwlgzMytAkQVjOTBe0lhJQ4DLgQVVGGtmZgUobJdUROyXNBtYBAwEbouIdZJm5cvnSToRWAEcAxyQdDUwISJ2dDa2qKxmZtazIo9hEBELgYUd2uZVTG8l292UNNbMzGrH3/Q2M7MkLhhmZpbEBcPMzJK4YJiZWRIXDDMzS+KCYWZmSVwwzMwsiQuGmZklccEwM7MkLhhmZpbEBcPMzJK4YJiZWRIXDDMzS+KCYWZmSVwwzMwsiQuGmZklccEwM7MkLhhmZpbEBcPMzJIUWjAkTZO0QVKrpDmdLJekL+fLV0uaVLHscUlrJK2StKLInGZm1rNBRT2xpIHArcCFQBuwXNKCiFhf0e1iYHz+OBv4av7zoKaIeLaojGZmlq7ILYwpQGtEbIyIvUAzML1Dn+nA7ZFZBhwr6aQCM5mZ2WFSRBTzxNJlwLSImJHPXwGcHRGzK/r8BLgpIn6ezy8BrouIFZI2AduBAL4WEfO7WM9MYCZAfX395Obm5kJeT2+1t7dTV1dX6xjdKnvGsueD8mcsez5wxv7Ql3xNTU0rI6IhpW9hu6QAddLWsTp11+fciNgi6QRgsaRHI2LpyzpnhWQ+QENDQzQ2NvYhcv9paWmhLFm6UvaMZc8H5c9Y9nzgjP2hWvmK3CXVBoyumB8FbEntExEHfz4N3Em2i8vMzGqkyIKxHBgvaaykIcDlwIIOfRYAV+ZnS50D/CEinpI0XNIIAEnDgYuAtQVmNTOzHhS2Syoi9kuaDSwCBgK3RcQ6SbPy5fOAhcC7gFZgN/CRfHg9cKekgxm/GxE/LSqrmZn1rMhjGETEQrKiUNk2r2I6gE92Mm4jcFaR2czMrHf8TW8zM0vigmFmZklcMMzMLIkLhpmZJXHBMDOzJC4YZmaWxAXDzMySuGCYmVkSFwwzM0vigmFmZklcMMzMLIkLhpmZJXHBMDOzJC4YZmaWxAXDzMySuGCYmVkSFwwzM0vigmFmZklcMMzMLEmhBUPSNEkbJLVKmtPJckn6cr58taRJqWP7010Pb2bGt5dz18Obi1yN9ULrtp1s372P1m07ax3FrNRuvHstj27dyY13ry18XYUVDEkDgVuBi4EJwAckTejQ7WJgfP6YCXy1F2P7xTmfW8zVd6zmvkee5uo7VjP1c4uLWI31wty71nDBLUtp276bC25Zytwfr6l1JLNSeuOce/j6A0+w78UDfP2BJ3jjnHsKXV+RWxhTgNaI2BgRe4FmYHqHPtOB2yOzDDhW0kmJY/vsroc3s3XH3kPantqx11saNdS6bSe3L3vykLbbH3zSWxpmHdx491oOdGg7kLcXRRFRzBNLlwHTImJGPn8FcHZEzK7o8xPgpoj4eT6/BLgOGNPT2IrnmEm2dUJ9ff3k5ubm5IxPPLebHc/ve1n7MUcP5g3HD0t+ns60t7dTV1fXp+coWhkzbt+9j7btuwGoHwrb9mTto44bxnHDBtcwWefK+G9Yqez5wBkP16Nbd7LvxaxkVL5XBg8cwGknjkh+nqamppUR0ZDSd1DvYyZTJ20dq1NXfVLGZo0R84H5AA0NDdHY2Jgc8K6HN/N3d6x+WfuX/uIMGieNTn6ezrS0tNCbLLVQxoyt23ZyzS1LAfjrM/bzj2uyP9H7rjmHcfXpb4JqKeO/YaWy5wNnPFwP3L2Wrz/wBHDoe+Xj576BxsbTC1lnkbuk2oDKT91RwJbEPilj++w9k0Zz0jFDDmk76ZghvKePxcIO37j6EVw59ZRD2q6cekopi4VZLd3w7tNf9gE+IG8vSpFbGMuB8ZLGAr8FLgf+e4c+C4DZkpqBs4E/RMRTkp5JGNsvHvzbC7nr4c38ZM1W/vyME10sSuCz08/gynPGsGblstJuWZiVwcab/owb717L4F2b+Pi5byi0WECBBSMi9kuaDSwCBgK3RcQ6SbPy5fOAhcC7gFZgN/CR7sYWlfU9k0a7UJTMuPoRtA0b7GJh1oMb3n06LS3PFrYbqlKRWxhExEKyolDZNq9iOoBPpo41M7Pa8Te9zcwsiQuGmZklccEwM7MkLhhmZpbEBcPMzJK4YJiZWRIXDDMzS+KCYWZmSVwwzMwsiQuGmZklKex+GLWQX7TwiVrnyI0Enq11iB6UPWPZ80H5M5Y9Hzhjf+hLvjdExOtSOr6iCkaZSFqRelOSWil7xrLng/JnLHs+cMb+UK183iVlZmZJXDDMzCyJC0Zx5tc6QIKyZyx7Pih/xrLnA2fsD1XJ52MYZmaWxFsYZmaWxAWjH0g6VtIPJD0q6RFJUyW9VtJiSb/Ofx5X44zXSFonaa2k70k6utYZJd0m6WlJayvauswk6XpJrZI2SHpnjfJ9If89r5Z0p6Rja5Wvq4wVy66VFJJGljGjpKvyHOskfb5WGbv4PU+UtEzSKkkrJE2pYb7Rku7PP1vWSfpfeXv13ysR4UcfH8C3gRn59BDgWODzwJy8bQ5wcw3znQxsAobm83cAH651RuBtwCRgbUVbp5mACcCvgKOAscBvgIE1yHcRMCifvrmW+brKmLePBhaRfS9pZNkyAk3AfcBR+fwJJfs93wtcnE+/C2ipYb6TgEn59AjgsTxH1d8r3sLoI0nHkP3BfRMgIvZGxO+B6WSFhPzne2qT8I8GAUMlDQKGAVuoccaIWAr8rkNzV5mmA80R8UJEbAJagSkUqLN8EXFvROzPZ5cBo2qVr6uMuVuAvwEqD1KWKeMngJsi4oW8z9O1ythFvgCOyadfQ/Z+qVW+pyLi4Xx6J/AI2X8Cq/5eccHouzcCzwDfkvSfkr4haThQHxFPQfYLB06oVcCI+C3wReBJ4CngDxFxb5kyVugq08nA5op+bXlbLX0U+Pd8ujT5JF0C/DYiftVhUWkyAm8G/lTSQ5J+JumteXtZMl4NfEHSZrL3zvV5e03zSRoDvAV4iBq8V1ww+m4Q2ebsVyPiLcAuss3D0sj3bU4n2zx9PTBc0v+obapeUydtNTvFT9INwH7g3w42ddKt6vkkDQNuAOZ2triTtlr9Gw4CjgPOAf43cIckUZ6MnwCuiYjRwDXkexCoYT5JdcAPgasjYkd3XTtp65eMLhh91wa0RcRD+fwPyArINkknAeQ/n+5ifDVcAGyKiGciYh/wI+C/lizjQV1laiPbL3/QKF7aTVBVkj4E/Dnwwch3GlOefG8i+4/BryQ9nud4WNKJlCcjeZYfReaXwAGy6yGVJeOHyN4nAN/npV06NcknaTBZsfi3iDiYq+rvFReMPoqIrcBmSafmTecD64EFZH905D9/XIN4Bz0JnCNpWP6/uPPJ9oOWKeNBXWVaAFwu6ShJY4HxwC+rHU7SNOA64JKI2F2xqBT5ImJNRJwQEWMiYgzZh8ek/O+0FBlzdwHvAJD0ZrKTRZ4tUcYtwNvz6XcAv86nq54vf89+E3gkIv6pYlH13ytFHt1/tTyAicAKYDXZG+E44HhgCdkf2hLgtTXO+PfAo8Ba4DtkZ1DUNCPwPbJjKvvIPtg+1l0msl0tvwE2kJ/BUoN8rWT7h1flj3m1ytdVxg7LHyc/S6pMGckKxL/mf48PA+8o2e/5PGAl2dlGDwGTa5jvPLJdSqsr/u7eVYv3ir/pbWZmSbxLyszMkrhgmJlZEhcMMzNL4oJhZmZJXDDMzCzJoFoHMCsbSS8Ca8jeH48AH4pDv3Nh9qrkLQyzl9sTERMj4nRgLzCrcqEyVXvvSBpYrXWZdccFw6x7/wGMkzQmvx/BV8i+aDZa0kWSHpT0sKTv59f6QdJNktbn98z4Yt72PmX3IvmVpKV524cl/cvBFUn6iaTGfLpd0mclPQRMlTQ5v0jfSkmLDl4SwqyaXDDMupBfCv5ist1TAKcCt8dLF5n8NHBBREwi+6b/X0l6LfBe4E8i4kzg/+Rj5wLvjIizgEsSVj+c7P4MZ5N90/ifgcsiYjJwG3Bjf7xGs97wMQyzlxsqaVU+/R9k1/F5PfBERCzL288hu1HNA9mlfhgCPAjsAJ4HviHpHuAnef8HgP8n6Q5euqhdd14ku9gcZIXqdGBxvq6BZJeyMKsqFwyzl9sTERMrG/IP6l2VTcDiiPhAx8H57TzPBy4HZpNdJ2mWpLOBPwNWSZpIdnn0yq38oyumn4+IFyvWtS4ipvbtZZn1jXdJmR2eZcC5ksZBdh8KSW/Oj2O8JiIWkt2EZ2K+/E0R8VBEzCW7KutosgsDTpQ0QNJour4r2gbgdZKm5s81WNKfFPnizDrjLQyzwxARz0j6MPA9SUflzZ8GdgI/lnQ02ZbBNfmyL0gan7ctIbsKKmT3Wl/DS1dt7WxdeyVdBnxZ0mvI3rdfAtb1+wsz64avVmtmZkm8S8rMzJK4YJiZWRIXDDMzS+KCYWZmSVwwzMwsiQuGmZklccEwM7MkLhhmZpbk/wNMXCyfmBbh1QAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXQAAADxCAYAAADbaUyMAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMi4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvIxREBQAAIABJREFUeJzt3XucVWXZ//HPd2YYYAAFREYEFCzymJIimloNeULTh6w0qNQ8UT5iaT0VWWlqByufTB8PiEVqJXQS42coojWeKhUMD6AoIMTI+ajDwJz29ftjrcHNsGfvNTP7sPZwvV+v9ZrZ63RfMy9e19xc6173LTPDOedc8SspdADOOeeywxO6c851EZ7QnXOui/CE7pxzXYQndOec6yI8oTvnXBfhCd0553JE0lBJf5f0mqSFkr6a4hxJuk3SEkkvSzo66dhYSYvDY5MztecJ3TnncqcJ+LqZHQocD1wh6bBW55wBjAi3icBdAJJKgTvC44cBE1JcuwtP6M45lyNmttrMXgy/fxd4DRjc6rRxwP0W+BfQV9IgYDSwxMyWmVkDMCM8t01lWf8JcmDAgAE2bNiwyOdv27aNXr165S6gLIh7jHGPD+IfY9zjg64Z4/z58zeY2b6dafP0Mb1s46bmzG29XL8Q2JG0a6qZTU11rqRhwIeA51odGgysTPpcE+5Ltf+4dPEURUIfNmwY8+bNi3x+dXU1VVVVuQsoC+IeY9zjg/jHGPf4oGvGKGlFZ9vcuKmZ5+cckPG80kFv7jCzURFi6g38GbjKzN5pfTjFJZZmf5uKIqE751w+GZAgkZV7SepGkMx/Z2YPpjilBhia9HkIsAoob2N/m7yG7pxzrRhGozVn3DKRJOBXwGtm9vM2TpsFXBCOdjke2Gpmq4EXgBGShksqB8aH57bJe+jOOZdClnroJwLnA69IWhDuuwY4AMDMpgCzgTOBJUAdcFF4rEnSJGAOUApMM7OF6RrzhO6cc60YRnMWphY3s2dIXQtPPseAK9o4Npsg4UfiCd0551JIpH/+GEue0J1zsfTmi8t488W3oLIBMyMoR+eHAc1FmND9oahzLnZm3TWHqz/6Pe686tesWrKGn37xdvK9uloCy7jFjSd051ys1G+vZ8rV91Jf10B9XT2WMJ558Dlef35J3mIwoNEs4xY3XnJxzsXKu5u3odJd+5olpSVsXrMlbzEY5iUX55zrrP779aVP/94kl8wTTQnef/Tw/AVh0BxhixtP6M65WCkpKeFnj1/LfsMrkURJaQnXPfgNBg4dkLcYgjdFM29x4yUX51zsDD14MPcvuZ3Ghkae/cezjKo6Ks8RiOb0w8djKac9dElXh5O6vyppuqQekvpLmivpzfBrv1zG4JwrXt3KuxWk3eChqDJucZOzhC5pMPAVYJSZHUHw6up4YDLwhJmNAJ4IPzvnXGwE49CVcYubXNfQy4CeksqACoKZwsYB94XH7wM+meMYnHOu3RKmjFvcKJeD9cP1834IbAceM7PPS9piZn2TztlsZruVXSRNJFiOicrKymNmzJgRud3a2lp69+7d6fhzKe4xxj0+iH+McY8PumaMY8aMmR9ljvJ0Djuy3H778H4ZzzvmwJWdbiubcvZQNKyNjwOGA1uAP0r6QtTrw1U/pgKMGjXK2jPBfVectD/f4h4fxD/GuMcHHmNbDNFchIMAcznK5RTgLTNbDyDpQeAEYK2kQWa2Olw3b10OY3DOuQ6JY0klk1wm9P8Ax0uqICi5nAzMA7YBFwI3hV//ksMYnHOu3QzRYKWFDqPdcpbQzew5SX8CXgSagH8TlFB6A3+QdAlB0j83VzE454rXupUbWLFwJQ005r3t4MUiL7nswsyuA65rtbueoLfunHMpPf3gc/zkgv+jrFsJZ107hmnPTOfiH0zIawxxHJaYSfH9CXLOdWmNDY385PzbqK+rZ9vW7VjCePAXD7P0peV5i8FMNFtJxi1u4heRc26PtnX9O7vtKy0rZc1b+R0/kUAZtygkTZO0TtKrbRz/hqQF4faqpGZJ/cNjyyW9Eh6bl6ktT+jOuVjpV9mX8h67vvLf3NjMsCOG5i2G4KFoWcYtonuBsW22ZfYzMxtpZiOBbwNPmtmmpFPGhMczjnf3hO6ci5XSslJ+OPs79O7Xix69uqMS8ZW7LmPw+wflLYaWh6KZtkj3MnsK2JTxxMAEYHoHw/bZFp1z8XPocSP4w+p72PD2Jl5/axFjPl6V9xia8zwOPRziPRaYlLTbgMckGXB3+MJlmzyhO+diafu7O9i0egvNTc15b7sdb4oOaFXbnpop6aZxNvBsq3LLiWa2StJAYK6k18Mef0qe0J1zsfPCnAXc8JmbKSkt5ezrqvj9gof47DfzO49fItoolg1ZnMtlPK3KLWa2Kvy6TtJMYDTQZkL3GrpzLlaaGpu44dz/Zce2eureqcMSxm+u/yPLF67MWwzB9LklGbdskbQ38DGS3pyX1EtSn5bvgdOAlCNlWngP3TkXK1vWbcWad13grbRbKW+/uZphh+dnpIshGrP06r+k6UAVQXmmhuBly24AZjYlPO0cghlptyVdWgnMVLC4ahnwgJk9mq4tT+jOuVjpO3BvSruVBjNAhZobmxl6yOC8xWBG1l4cMrOMr7ia2b0EwxuT9y0D2rX2npdcnHOxUtatjBse+hYVfXpSsVdPVCIu/ekXOCCPCZ0ILxVFfbEon7yH7pyLnaOqDmd6zd2sXraW5euXcvIp+Z3+ycheDz2fii9i59weoaJPT9531DBKywozjW0+H4pmi/fQnXOuFSOea4Zm4gndOedaMaAx+lwtsVF8ETvnXM6pKOdD94TunHOtGJHfFI0VT+jOOZeC99Cdc64LMJP30J1zLlt21O1g6b+Xk2g1DUA+BA9FCzNcsjOK70+Qc67Le+y+v3N2n/O56iPfY+mC5fzi8rvzHIGvKeqcc53WsKOBmy++M+gmh/569+P8+2+v5C2G4KGoMm5x4wndORcry1/9D2a7758/9+W8xuFvijrnXCcN/sD+Kfe/b+SwvMVQrG+Kxu9PjHNuj9Zrrwo+/bWzdtl3yOj3M+azJ+Y1jmwtEp1P3kN3zsXOl2++kI+ccxzz575E//f15rJ/XZzX9s2gMRG/hJ1JzhK6pIOB3yftOgi4FugLXAasD/dfY2azcxWHc644HfrhD3DAoXsx/+U38t52UHIpvoSes4jNbLGZjTSzkcAxQB0wMzx8S8sxT+bOudZqFlWz6dUjKa8dQ6JhIYuevDPvMTSH87mk2+ImX3+CTgaWmtmKPLXnnCtSzU0NVCSuoO+ABrqVGyUlxrADbmPdipfyFoMPW0xvPDA96fMkSS9LmiapX55icM4VgS3r3qJnryZKkrJTc5PYuOIfeYwiKLlk2iLdKchz6yS92sbxKklbJS0It2uTjo2VtFjSEkmTM7ZlqQZ8ZpGkcmAVcLiZrZVUCWwg+CN4IzDIzHZ74iFpIjARoLKy8pgZM2ZEbrO2tpbevXtnI/yciXuMcY8P4h9j3OODeMZoloDG11BJkJtqt1dS0X0dTc1DKe+5V8brx4wZM9/MRnUmhn0PG2Dj7j8r43m/Ova+jG1J+ihQC9xvZkekOF4F/I+ZndVqfynwBnAqUAO8AEwws0VttZWPUS5nAC+a2VqAlq8Aku4BHk51kZlNBaYCjBo1yqqqqiI3WF1dTXvOL4S4xxj3+CD+McY9PohvjC8/9k9GvO83NDfD88uuZJ+Gao466495az8Y5ZKduVzM7ClJwzpw6WhgiZktA5A0AxgHtJnQ81FymUBSuUXSoKRj5wAp/xvinNtzHXnad9hUfw/LV1xKY9PgvCZzeO/Fogg19AGS5iVtEzvY5IclvSTpEUmHh/sGAyuTzqkJ97Uppz10SRUE/134UtLun0oaSVByWd7qmHPOATD4kI8w+JCPsKG6uiDtJ6KNYtnQ2fIO8CJwoJnVSjoTeAgYASkDSFsjz2lCN7M6YJ9W+87PZZvOOddZLaNc8tKW2TtJ38+WdKekAQQ98qFJpw4heB7ZJn9T1DnnUsjXi0WS9gPWmplJGk1QCt8IbAFGSBoOvE0wWvBz6e7lCd0551oxE01ZSuiSpgNVBPX2GuA6oFvQjk0BPgNcLqkJ2A6Mt2D4YZOkScAcoBSYZmYL07XlCd0551LIVsnFzCZkOH47cHsbx2YDkd+m94TunHOt5LOGnk2e0J1zLgVP6M451wX4AheuS3jur/M5b9ClvPniW3zj5OvZuuGdzBc51wUlUMYtbjyhu51WLFrJjef9nM1rt4IZrz7zGted87NCh+Vc3plBU6Ik4xY3XnJxO7385K5TRDQ1NrPon4tJJBKUlMTvH69zuVSMJRdP6G6nvfbpg0p3/Ufco6K7J3O3x/Eauit6J3zyWA48dAg9enUHQfeKcq6849JCh+VcQZgp4xY33kN3O3Ur78YtT9/I36c/y7YeW/nZE9/n0ONGFDos5woijg89M/EeuttFt/JunHZhFf0q9/Zk7vZYZr4EXcHN/c2TjB8ykaUvLee2Sb+kqbGp0CE554qSaE6UZNzipsuUXF58/GVuvXwq9XUNJJoSPHbv3+nWrZTLb7mo0KE554pQHGvkmcTvT0wHPf3gc9TXNez8XF/XwNN/fq6AETnnilXLXC7FVnLpMj30vfbpTWlZCc1NiZ37KvbqWcCInHNFy4I6erHpMj30cVeMpXe/3pSVB3+juleUc/ktXyxsUM65olWMr/53mR56//36cc8rP2fufdWU7p/gF0//gPd/aHihw3LOFSELH4oWm+KLOI1+A/fmvG+MY59B/TyZO+c6xSzzFjeRErqkIZLGhN93l9Qrt2E551xhFeObohkTuqSLgVnAL8NdBwJ/yWVQzjlXSEEPPDsJXdI0SeskvdrG8c9Lejnc/iHpqKRjyyW9ImmBpHmZ2orSQ/8KcDzwTvCD2hvAwEg/iXPOFaksDlu8Fxib5vhbwMfM7EjgRmBqq+NjzGykmY3K1FCUh6I7zKxBCoKXVAoxfLzrnHNZlK0auZk9JWlYmuP/SPr4L2BIR9uK0kN/VtI3gR5hHf33wMMdbdA55+LOEIlEScYNGCBpXtI2sZNNXwI8skso8Jik+VHuHaWH/k1gIvA68FVgDnB3BwJ1zrmiEbGDviFKKSSKsMN8CXBS0u4TzWyVpIHAXEmvm9lTbd0jbUIPyyvTzOxC4K5sBO2cc7Fn+Z3LRdKRBANPzjCzjTvDMFsVfl0naSYwGmgzoactuZhZMzBIUrcOBHhw+GS2ZXtH0lWS+kuaK+nN8Gu/9t67Lc/8ZwXjZvyWNzdt5O75z2NxHCgac2+9soKvVV3H8oUrufOqX9NQ31jokNweaP3GrZx3+jWcMuhClry2kgWvLM1/EBZhywJJBwAPAueHg05a9veS1Kfle+A0IOVImRZRSi7LgKcl/QXY1rLTzG5Ld5GZLQZGhsGUAm8DM4HJwBNmdpOkyeHnb0WII62X1qxm4sMPsaOpiVP3G8Jdz/2TpkSCK449vrO33mOsr9nIVSd9j+212zn07AOZfc8TbF67he9Mv7rQobk9SCKR4PyTvkXj0k2UNBmJuka+ftL3+O2bdzBoYNb6fxllq4cuaTpQRVBvrwGuA7oFbdgU4FpgH+DOcPBJU1jGqQRmhvvKgAfM7NF0bUVJ6OuBuUBFuHXEycBSM1shaRzBDwdwH1BNFhL6Q4tfY0fTe/Ofb29qYsarL3tCb4cXHvk3zc3NO5/u128PZqz0RaJdPq1YuY6mMJlDOKSuMcHDf/0nl110Zl5iMCCRyE5CN7MJGY5fCuy21qOZLQOO2v2KtikfZQlJ04AXzex2SVvMrG/Ssc1mttuf3fCJ7kSAysrKY2bMmJG2jTW1tWyoC/4DUdmtnLWNDZSXlvKBfQZk80fJmtraWnr37l3oMHbxzsZa1v1nPZYw+g3Zm801W5HE+4+O5zQKcfwdJot7fBDPGBubmnnr5RUoTE39huzNpre3ss+wfdmnf5+M148ZM2Z+Zx9Udj9osA350RUZz1s24TudbiubMvbQJc0lRbXIzE6L0oCkcuC/gG+3JzAzm0o4wH7UqFFWVVWV9vwVW7Zw9vTfsK2xgav3G8KdG9Zw08mnUXXwoe1pNm+qq6vJ9DPl27at27j0iK+xZd1WPvWj05j1/b/xqavOil2cLeL4O0wW9/ggvjGe+5PJbHr6LUoaEpx78xn88RfP8NDiu+jVs3veYijGR3BRSi7fTfq+B/BpoL4dbZxB0DtfG35eK2mQma2WNAhY1457tenAvn35y4Qv8MsX57F3bR13nflffGxYPHuWcdVr717c9eJPmf7jmfTu15NJt1/CaRdUFTostwf63ewfcN13p/H6P9+gR/9eTF/wi7wmcyBrDz3zKWNCN7PWy/48KenJdrQxAZie9HkWcCFwU/g1a/PCDO/bjx9+/FSqq6s9mXdQ33335vKffzG2PTe3ZyjvVsaPfxK8R1NdXc2A/nvlOYJ4Tr6VSZTJufZK2vpKOhkYFOXmkiqAUwmG5LS4CThV0pvhsZs6EHdKi/65mG+eegMrX3+bv94z14ctOuc6Lk/DFrMpSsllIUHoApoIJpK5LMrNzayOYDhO8r6NBKNesmrJgrf45qk3Ul9Xz4ixQ5hy3X3U1zXwqa9+IttNOee6OgPL0iiXfIoyFu0gMzvAzIaa2XAz+zjwbK4Da6/H7qumvu690v6Ounpm3jq7gBE554qbImzxEiWht66hAzyf7UA6q6SkZLffr0ri9wt3zhWJIiy5tJnQJQ0MJ1rvKemDko4Mt5Po+AtGOXPGpSfTo+K9p+DdK7ozfvInCxiRc66oFWFCT1dD/wRwMcHcvHcm7X8X+F4ug+qIAw8dwq3P/pAHfvhnKvau4BvT/puPnXdCocNyzhUjA4pwlEubCd3Mfg38WtJ5ZvaHPMbUYQcdeSDf/f3XgmGLVZ7MnXMdV4yD5KKMQ/+DpNOBwwleLGrZ/6NcBtYRS19azu9+8Gfed9r+VK99lqrPnljokJxzxaoIR7lEefX/TqAv8FHg1wRviv4rx3G124rXarjqpO+yY1s9g44/g5uvvou6d7dz5qWnFDo051wRUhH20KOMcjnJzD4HbDSz7wHH0Yk173LlkV8+wY5t7w1brK+rZ8ZNMwsYkXOuaEV5IBrDhB9pkeiWr5L2AzYCw3IWUQclEond9hVjDcw5FwcqyoeiUXrosyX1BW4GFgDLgT/lMqiOOP2LY+jeatjiOV/Jz9zJzrkuqAh76GkTuqQS4BEz22JmfwSGAx80s2vyEl07vO+oYfz08Ws5+tQj6dG7B//9iy96QnfOdVwiwhYzmdYUTQC3Jn3ebmabch5VBx12/Af4yZzvMfTg/Tnz0lMIl25yzrn2aRmHnmmLmSgll7nhsnGxV/PGKn4+cQpr3lrH84/8u9DhOOc6qKGhkes//TPGD5nIytffZsOq/PcjZZm3SPeRpklaJynlAs8K3CZpiaSXJR2ddGyspMXhscmZ2oqS0CcRLFS6XdImSZslxa6XvmrpGq44djKP/upvvLuplhvOvZm/PfB0ocNyznXAJYdexTMzn2fjqs3s2FbPBe+7gm3v1OU3iOzV0O8FxqY5fgYwItwmAncBSCoF7giPHwZMkHRYuoaiJPQBBCtU9wb2DT/vG+G6vHr47sfYXrtj5xzo9XUN3H99Ubzg6pxLsnbFOta8tetCZo31TTw85bECRdQ5ZvYUkK4TPA643wL/AvqGq7mNBpaY2TIzawBmhOe2KWNCN7Nm4FzgW+H3g4CR0X6U/GnY0bTbghaNDU0FisY511E76hpS7q/fnnp/rkQsuQyQNC9pm9iBpgYDK5M+14T72trfpigrFt0OjAHOD3fVAVPaEWxenPKFj9C9onzn5x4V3fnEZacWMCLnXEcMPXh/Kvbqucs+SYy9aEz+gjCCV/8zbbDBzEYlbVM70Fqqp6uWZn+bopRcTjCzLxG+YBSOcilPf0n+HTJ6BDf8ZTIHj34/5T3LOf/75zHh2+cUOiznXDuVlJTwq0W/YPCI/SjrVkppWSm3PHUDAw/Ic6U3f+PQa4ChSZ+HAKvS7G9TlDdFG8Px6AYgaR9iOQITjj75gxx98o99gWPnityA/ftz7+L/A4JFog8/8ZC8x5DHuVxmAZMkzSCYWmWrma2WtB4YIWk48DYwHvhcuhtFSeh3AH8G9pV0PXAecH1nos+Vdf9Zz6w757DX4T145enX+OBHDi10SM65YpWlhC5pOlBFUG+vAa4jGGiCmU0BZgNnAksIStoXhceaJE0C5gClwDQzW5iurSjT594vaT7QMm3huWaWcjxlIa1buYEvjfwGde9u5zM3nc63L/8B1/zuKk4Yd2yhQ3POFaMsJXQzm5DhuAFXtHFsNkHCjyRKDR2Cvw6NQEM7rsmrWXfOoe7d7SSag2pQfV0Dv/z2bwsclXOuGEUZ4RLH6XWjjHL5DjAd2J+gKP+ApG/nOrD22l77XjJvUd/G8CfnnMso2iiXWInS2/4CcKyZfdfMvkMw2P2C3IbVfmM+eyLde743+KZ7RXdOOf+jBYzIOVfMumQPHVjBrrX2MmBZbsLpuCNOOpRrHriKoYfsT1l5GZ+++hNc8P3zCh2Wc65YFeH0uVFGudQBCyXNIfgRTgOekfRzADP7Wg7ja5cTxh3LCeOO9WGLzrnOiWkPPJMoCf2v4dYi8nqi4cIYvwSOIPhjcDFwOnAZsD487ZrwSa5zzsVHV0zoZvarTtz/VuBRM/uMpHKggiCh32JmN3fivs45l1OK5euT6UUZ5TJW0gvhfL6Rp8+VtBfwUeBXAGbWYGZbOh+yc865VKKUXG4neDv0Fdr3yv9BBGWVX0s6CpgPfDU8NknSBcA84Otmtrn1xeGsZRMBKisrqa6ujtxwbW1tu84vhLjHGPf4IP4xxj0+8BjT6oolF4IJYhaEy9G1995HA1ea2XOSbgUmE/yBuJHg13Uj8L8EtfVdhLOWTQUYNWqUtechZzE8FI17jHGPD+IfY9zjA4+xTV34oeg3gf8nqRqob9lpZrdluK4GqDGz58LPfwImm9nalhMk3QM83K6InXMuH4owoUcZh3490Az0JVipqGVLy8zWACslHRzuOhlYFK7E0eIcIHbzwjjnXFcdhz7QzI7p4P2vBH4XjnBZRjCL2G2SRhL8OpYDX+rgvZ1zLidEcY5yiZLQn5D0cTP7W3tvbmYLgFGtdp+f6lznnIuNIq2hRym5XAY8Lqm2PcMWnXOuqHXRksuAnEfhnHNxE8OEnUnGHrqZNQPnAt8Kvx8EjMx1YM45V0hdcrZFSbcDY3iv9l0HTMllUM45V3BFWHKJUkM/wcy+BOwAMLNNQHn6S5xzrohZMMol0xZFOH3KYklLJE1OcfwbkhaE26uSmiX1D48tl/RKeGxeprai1NAbJZUQ/j2StA/tmwLAOefa7aU1q3ll3Vr2bmjAzJDyvEJQFnrgkkqBO4BTCV62fEHSLDNbtLMZs58BPwvPPxu4Ouw4txhjZhuitNdmD11SS7K/A/gzsK+k64FngJ9E/5Gcc6597l3wIp978A/86OknWbl1C1fNmU2wlnL+ZKmGPhpYYmbLzKwBmAGMS3P+BIIlPzskXcnleQAzux/4LnAzsBk418xmdLRB55xLZ0dTIz9+5km2NzWxo7mJhBlPLFvKgjWr8xtItBr6AEnzkraJre4yGFiZ9Lkm3LcbSRXAWIIOdHIUj0man+Leu0lXctn5/xszWwgszHQz55zrrK076ilpVV4pKREb6uryF0T0h54bzKz1y5PJUtWJ2rrz2cCzrcotJ5rZKkkDgbmSXjezp9pqLF1C31dSm8vLmdnP01zrnHMdsm+vXvTvWcGa2nd3Zr7mRIIjBlbmLQaRtWGJNcDQpM9DgFVtnDueVuUWM1sVfl0naSZBCafNhJ6u5FIK9Ab6tLE551zWlUj89lPnMrxfPwSUlpQw9exPMqhPftNOlmroLwAjJA0P57QaD8zarS1pb+BjwF+S9vWS1Kfle4L1nNNOZpiuh77azG6IFLJzzmXR8L79ePz8i2lOJHj6qac4ceiB+Q8iCz10M2uSNAmYQ9BJnmZmCyV9OTze8k7POcBjZrYt6fJKYGY4uqcMeMDMHk3XXqQaunPOFUJpSZRXZXIkS4NqzGw2MLvVvimtPt8L3Ntq3zLgqPa0lS6hn9yeGznnXJcR01f7M2kzobd60uqcc3uWrpTQnXNuT9ZVF7hwzrk9TpcquTjn3B4rprMpZuIJ3TnnUvGE7pxzxS+Lb4rmlSd055xLQYniy+ie0J1zrjWvoTvnXNfhJRfnnOsqPKE751zX4D1055zrKoowoed0KjNJfSX9SdLrkl6T9GFJ/SXNlfRm+LVfLmNwzrl2s+DV/0xb3OR6bspbgUfN7BCCaSBfAyYDT5jZCOCJ8LNzzsVGyzj0LCxwkVc5S+iS9gI+CvwKwMwazGwLwYrX94Wn3Qd8MlcxOOdch5ll3mImlzX0g4D1wK8lHQXMB74KVJrZagAzWx0ufrqbcIXriQCVlZVUV1dHbri2trZd5xdC3GOMe3wQ/xjjHh94jOnEsQeeSS4TehlwNHClmT0n6VbaUV4xs6nAVIBRo0ZZVVVV5Iarq6tpz/mFEPcY4x4fxD/GuMcHHmObivTFolzW0GuAGjN7Lvz8J4IEv1bSIIDw67ocxuCccx2SrYeiksZKWixpiaTdOrWSqiRtlbQg3K6Nem1rOUvoZrYGWCnp4HDXycAighWvLwz3XUjSKtfOORcX2UjokkqBO4AzgMOACZIOS3Hq02Y2MtxuaOe1O+V6HPqVwO8klQPLgIsI/oj8QdIlwH+Ac3Mcg3POtY+RrYeeo4El4YLPSJpBMDBkUS6uzWlCN7MFwKgUh3wBaudcrEV8KDpA0rykz1PD538tBgMrkz7XAMeluM+HJb0ErAL+x8wWtuPanfxNUeecSyVaQt9gZqk6rS0U4c4vAgeaWa2kM4GHgBERr91Frl8scs65opPFF4tqgKFJn4cQ9MJ3MrN3zKw2/H420E3SgCjXtuYJ3TnnWjNDicxbBC8AIyQND58ljicYGLKTpP0kKfxCpWOpAAAHcUlEQVR+NEFe3hjl2ta85OKcc6lk4ZmomTVJmgTMAUqBaWa2UNKXw+NTgM8Al0tqArYD483MgJTXpmvPE7pzzqWQrTdFwzLK7Fb7piR9fztwe9Rr0/GE7pyLpfqmJtbU1pIoxJwpBviaos4513n/Xr2Ki2Y9SFNzgsv33Y9Ni17lM4cdkd8gii+f+0NR51y8NCUSXDzrQd6pr6euqREz49rqJ1i+ZXNe4/Dpc51zrpM21G2jvrl5l31lJSW8sXFDXuPI0iiXvPKSi3MuVvr3rNhtX1MiwdC9++YvCJ9t0TnnOq+8tJRbTjuTnmVl9CkvRxKXfOgYDh2wb95iCF4ssoxb3HgP3TkXO6e/fwSPV17Mkk0bqVuylNM/fFL+g4jhmqGZeA/dORdLlb17c9i+A+leWlqQ9r2H7pxzWfDa+nVc8NCfebehnisH7s+ON99g7IgP5C8Ar6E751znNScSXPDQn9i4vY6G5mbMjK/NfYSVW7fmMYqszeWSV57QnXOxsr5uG7UNjbvs61ZSwusb1uc3ELPMW8x4ycU5Fyv9evSkdb2jKZFgUJ8++QvCoq8ZGifeQ3fOxUr3sjJ++PFT6VFWRq9uwbDF8w77IEcMrMxvIN5Dd865zvvUoYczcr9BvL5hA1qxgiuqPp7/IOKXrzPyhO6ci6WD+vXnoH79qX477SI9OaNE8dVcPKE751xrRlG+WOQJ3TnnWhHxfHEoE0/ozjmXShEmdB/l4pxzqWRplIuksZIWS1oiaXKK45+X9HK4/UPSUUnHlkt6RdICSfMyteU9dOecay1LNXRJpcAdwKlADfCCpFlmtijptLeAj5nZZklnAFOB45KOjzGzSJPBe0J3zrkUsjTKZTSwxMyWAUiaAYwDdiZ0M/tH0vn/AoZ0tDEvuTjn3G4ilFuilVwGAyuTPteE+9pyCfDIroHwmKT5kiZmaiynPXRJy4F3gWagycxGSfo+cBnQMjHDNWY2O5dxOOdcuxhRE/aAVrXtqWY2Nemz2rj7biSNIUjoyZO/n2hmqyQNBOZKet3MnmormHyUXFLVf24xs5vz0LZzznVMtIrLBjMbleZ4DTA06fMQYLc3pSQdCfwSOMPMNrbsN7NV4dd1kmYSlHDaTOhecnHOuRSytMDFC8AIScMllQPjgVm7tCMdADwInG9mbyTt7yWpT8v3wGnAq+kay3UPvaX+Y8DdSf8VmSTpAmAe8HUz29z6wrBeNBGgsrKS6urqyI3W1ta26/xCiHuMcY8P4h9j3OMDjzGtLIxDN7MmSZOAOUApMM3MFkr6cnh8CnAtsA9wpyQIy9NAJTAz3FcGPGBmj6ZrT5bDwfOS9k+u/wBXAouBDQTJ/kZgkJldnO4+o0aNsnnzMg7B3Km6upqqqqoOx50PcY8x7vFB/GOMe3zQNWOUND9DGSSjvXvsZycccGHG8x5986edbiubclpySa7/ADOB0Wa21syazSwB3ENQE3LOuXgpwulzc5bQ26r/SBqUdNo5ZKgJOedcQRRhQs9lDT1l/UfSbySNJCi5LAe+lMMYnHOu/QyI4ZqhmeQsoYdvRh2VYv/5uWrTOeeyw8CKb/5cf/XfOedaM6DZE7pzznUNMayRZ+IJ3TnnUvGE7pxzXUE8R7Fk4gndOedaM6AIF4n2uVycc7FjiU0kNk0kse7D0LQUa3yzAEEU3zh0T+jOuVgxM2zTRdDwDCQ2AjuwTROwxG5TPuUyimCUS6YtZjyhO+fiJbERmpYCTeGOcD24hgX5i8HALJFxixuvoTvn4kU92H0y8gSoIr9xFOGbot5Dd87Fikp6Q8XngZ7hnhIo+wCUH5PfQIqwhu49dOdc7KjPNdDtKKxxAZRUov6/RcpjujIrylEuntCdc7EjCXqehXqeBSXVBIv95FkMe+CZeEJ3zrndGNbcXOgg2s0TunPOtebT5zrnXBcSw2GJmfgoF+eca8UAS1jGLQpJYyUtlrRE0uQUxyXptvD4y5KOjnpta57QnXOuNQsXuMi0ZSCpFLgDOAM4DJgg6bBWp50BjAi3icBd7bh2F15ycc65FLL0UHQ0sCRcwQ1JM4BxwKKkc8YB95uZAf+S1Ddce3lYhGt3URQJff78+RskrWjHJQcA/8lVPFkS9xjjHh/EP8a4xwddM8YDO9vgu2ye87j9aUCEU3tImpf0eaqZTU36PBhYmfS5Bjiu1T1SnTM44rW7KIqEbmb7tud8SevNbFSu4smGuMcY9/gg/jHGPT7wGNtiZmOzdCulun3Ec6Jcu4uiSOgdsKXQAUQQ9xjjHh/EP8a4xwceY67VAEOTPg8BVkU8pzzCtbvoqg9FtxY6gAjiHmPc44P4xxj3+MBjzLUXgBGShit43XU8MKvVObOAC8LRLscDW81sdcRrd9FVe+hTM59ScHGPMe7xQfxjjHt84DHmlJk1SZoEzAFKgWlmtlDSl8PjU4DZwJnAEqAOuCjdtenakxXhfAXOOed211VLLs45t8fxhO6cc12EJ3TnnOsiPKE751wX4QndOee6CE/ozjnXRXhCd865LuL/Axwiu7MmG3OGAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 432x288 with 2 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "%matplotlib inline\n",
    "pd.set_option('mode.chained_assignment',None) # this removes a useless warning from pandas\n",
    "import matplotlib.pyplot as plt\n",
    "\n",
    "data[\"Frequency\"]=data.Malfunction/data.Count\n",
    "data.plot(x=\"Temperature\",y=\"Frequency\",kind=\"scatter\")\n",
    "plt.grid(True)\n",
    "\n",
    "data.plot(x=\"Pressure\", y=\"Frequency\", kind=\"scatter\")\n",
    "plt.grid(True)\n",
    "\n",
    "data.plot(x=\"Pressure\", y=\"Temperature\", kind=\"scatter\", color=data[\"Malfunction\"], colormap='viridis')\n",
    "plt.grid(True)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "À première vue, ce n'est pas flagrant mais bon, essayons quand même\n",
    "d'estimer l'impact de la température $t$ sur la probabilité de\n",
    "dysfonctionnements d'un joint. \n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Estimation de l'influence de la température\n",
    "\n",
    "Supposons que chacun des 6 joints toriques est endommagé avec la même\n",
    "probabilité et indépendamment des autres et que cette probabilité ne\n",
    "dépend que de la température. Si on note $p(t)$ cette probabilité, le\n",
    "nombre de joints $D$ dysfonctionnant lorsque l'on effectue le vol à\n",
    "température $t$ suit une loi binomiale de paramètre $n=6$ et\n",
    "$p=p(t)$. Pour relier $p(t)$ à $t$, on va donc effectuer une\n",
    "régression logistique."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<table class=\"simpletable\">\n",
       "<caption>Generalized Linear Model Regression Results</caption>\n",
       "<tr>\n",
       "  <th>Dep. Variable:</th>      <td>Frequency</td>    <th>  No. Observations:  </th>  <td>    23</td>  \n",
       "</tr>\n",
       "<tr>\n",
       "  <th>Model:</th>                 <td>GLM</td>       <th>  Df Residuals:      </th>  <td>    21</td>  \n",
       "</tr>\n",
       "<tr>\n",
       "  <th>Model Family:</th>       <td>Binomial</td>     <th>  Df Model:          </th>  <td>     1</td>  \n",
       "</tr>\n",
       "<tr>\n",
       "  <th>Link Function:</th>        <td>logit</td>      <th>  Scale:             </th> <td>  1.0000</td> \n",
       "</tr>\n",
       "<tr>\n",
       "  <th>Method:</th>               <td>IRLS</td>       <th>  Log-Likelihood:    </th> <td> -3.9210</td> \n",
       "</tr>\n",
       "<tr>\n",
       "  <th>Date:</th>           <td>Fri, 26 Jun 2020</td> <th>  Deviance:          </th> <td>  3.0144</td> \n",
       "</tr>\n",
       "<tr>\n",
       "  <th>Time:</th>               <td>12:44:47</td>     <th>  Pearson chi2:      </th>  <td>  5.00</td>  \n",
       "</tr>\n",
       "<tr>\n",
       "  <th>No. Iterations:</th>         <td>6</td>        <th>  Covariance Type:   </th> <td>nonrobust</td>\n",
       "</tr>\n",
       "</table>\n",
       "<table class=\"simpletable\">\n",
       "<tr>\n",
       "       <td></td>          <th>coef</th>     <th>std err</th>      <th>z</th>      <th>P>|z|</th>  <th>[0.025</th>    <th>0.975]</th>  \n",
       "</tr>\n",
       "<tr>\n",
       "  <th>Intercept</th>   <td>    5.0850</td> <td>    7.477</td> <td>    0.680</td> <td> 0.496</td> <td>   -9.570</td> <td>   19.740</td>\n",
       "</tr>\n",
       "<tr>\n",
       "  <th>Temperature</th> <td>   -0.1156</td> <td>    0.115</td> <td>   -1.004</td> <td> 0.316</td> <td>   -0.341</td> <td>    0.110</td>\n",
       "</tr>\n",
       "</table>"
      ],
      "text/plain": [
       "<class 'statsmodels.iolib.summary.Summary'>\n",
       "\"\"\"\n",
       "                 Generalized Linear Model Regression Results                  \n",
       "==============================================================================\n",
       "Dep. Variable:              Frequency   No. Observations:                   23\n",
       "Model:                            GLM   Df Residuals:                       21\n",
       "Model Family:                Binomial   Df Model:                            1\n",
       "Link Function:                  logit   Scale:                          1.0000\n",
       "Method:                          IRLS   Log-Likelihood:                -3.9210\n",
       "Date:                Fri, 26 Jun 2020   Deviance:                       3.0144\n",
       "Time:                        12:44:47   Pearson chi2:                     5.00\n",
       "No. Iterations:                     6   Covariance Type:             nonrobust\n",
       "===============================================================================\n",
       "                  coef    std err          z      P>|z|      [0.025      0.975]\n",
       "-------------------------------------------------------------------------------\n",
       "Intercept       5.0850      7.477      0.680      0.496      -9.570      19.740\n",
       "Temperature    -0.1156      0.115     -1.004      0.316      -0.341       0.110\n",
       "===============================================================================\n",
       "\"\"\""
      ]
     },
     "execution_count": 4,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "import statsmodels.api as sm\n",
    "\n",
    "data[\"Success\"]=data.Count-data.Malfunction\n",
    "data[\"Intercept\"]=1\n",
    "\n",
    "logmodel=sm.GLM(data['Frequency'], data[['Intercept','Temperature']], family=sm.families.Binomial(sm.families.links.logit)).fit()\n",
    "\n",
    "logmodel.summary()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "L'estimateur le plus probable du paramètre de température est 0.0014\n",
    "et l'erreur standard de cet estimateur est de 0.122, autrement dit on\n",
    "ne peut pas distinguer d'impact particulier et il faut prendre nos\n",
    "estimations avec des pincettes.\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Estimation de la probabilité de dysfonctionnant des joints toriques\n",
    "La température prévue le jour du décollage est de 31°F. Essayons\n",
    "d'estimer la probabilité de dysfonctionnement des joints toriques à\n",
    "cette température à partir du modèle que nous venons de construire:\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXcAAAEKCAYAAADpfBXhAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMi4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvIxREBQAAIABJREFUeJzt3Xl4VOXd//H3dyb7QmLYISA7yA5hEXEBrYK2KiriinVBpHWp7SNVn199tE+16oNt1VZxQ3GpgisupYJa44JbQBBkX8UEkJ0kkD33748ZMGAgQzLJLPm8rivXzDlzn3O+dwY+c3LmnPuYcw4REYkunlAXICIiwadwFxGJQgp3EZEopHAXEYlCCncRkSikcBcRiUI1hruZPW1mW83s28O8bmb2sJmtMbPFZjYw+GWKiMjRCGTPfTow+givnwl09f9MBKbWvSwREamLGsPdOfcxsPMITc4FnnM+XwDpZtY6WAWKiMjRiwnCOtoC31eZzvXP23xoQzObiG/vnsTExKx27drVaoOVlZV4PNHxdYH6Ep6ipS/R0g9QX/ZbtWrVdudc85raBSPcrZp51Y5p4Jx7AngCYNCgQW7+/Pm12mB2djYjRoyo1bLhRn0JT9HSl2jpB6gv+5nZd4G0C8bHYC5QdRc8E9gUhPWKiEgtBSPc3wKu8J81czywxzn3k0MyIiLScGo8LGNmLwEjgGZmlgvcCcQCOOceA2YDZwFrgH3AVfVVrIiIBKbGcHfOXVLD6w64PmgViUhEKCsrIzc3l+Li4gbZXlpaGsuXL2+QbdW3QPqSkJBAZmYmsbGxtdpGML5QFZFGKDc3l9TUVDp06IBZdedVBFdBQQGpqan1vp2GUFNfnHPs2LGD3NxcOnbsWKttRMd5RSLS4IqLi2natGmDBHtjY2Y0bdq0Tn8VKdxFpNYU7PWnrr9bhbuISBTSMXcRiVher5c+ffocmJ41axYdOnQIXUFhROEuIhErMTGRRYsWHfb18vJyYmIaZ8zpsIyIRJXp06dz4YUXcvbZZ3PGGWcAMGXKFAYPHkzfvn258847D7S955576N69Oz/72c+45JJLeOCBBwAYMWIE+4dH2b59+4G/BioqKpg8efKBdT3++OPAj8MJjB07lh49enDZZZfhO0sccnJyOOGEE+jXrx9DhgyhoKCAUaNGHfShNHz4cBYvXhzU30Pj/EgTkaD649tLWbYpP6jr7NmmCXee3euIbYqKiujfvz8AHTt25I033gDg888/Z/HixWRkZDB37lxWr17NV199hXOOc845h48//pjk5GRmzJjBwoULKS8vZ+DAgWRlZR1xe9OmTSMtLY2cnBxKSkoYPnz4gQ+QhQsXsnTpUtq0acPw4cOZN28eQ4YM4aKLLmLmzJkMHjyY/Px8EhMTueKKK5g+fToPPvggq1atoqSkhL59+wbht/YjhbuIRKzDHZY5/fTTycjIAGDu3LnMnTuXAQMGAFBYWMjq1aspKCjgvPPOIykpCYBzzjmnxu3NnTuXxYsX8+qrrwKwZ88eVq9eTVxcHEOGDCEzMxOA/v37s2HDBtLS0mjdujWDBw8GoEmTJgCcd955DB8+nClTpvD0009z5ZVX1u0XUQ2Fu4jUWU172A0tOTn5wHPnHLfffjvXXXfdQW0efPDBw55uGBMTQ2VlJcBB55o75/j73//OqFGjDmqfnZ1NfHz8gWmv10t5eTnOuWq3kZSUxOmnn86bb77Jyy+/TG1HyD0SHXMXkag2atQonn76aQoLCwHIy8tj69atnHzyybzxxhsUFRVRUFDA22+/fWCZDh06sGDBAoADe+n71zV16lTKysoAWLVqFXv37j3stnv06MGmTZvIyckBfFemlpeXAzBhwgRuuukmBg8efOCvjGDSnruIRLUzzjiD5cuXM2zYMABSUlJ44YUXGDhwIBdddBH9+/fn2GOP5aSTTjqwzC233MK4ceN4/vnnOfXUUw/MnzBhAhs2bGDgwIE452jevDmzZs067Lbj4uKYOXMmN954I0VFRSQmJvL+++8DkJWVRZMmTbjqqnoaa9E5F5KfrKwsV1sffvhhrZcNN+pLeIqWvtRnP5YtW1Zv665Ofn5+va7/zjvvdFOmTKnXbeyXn5/v8vLyXNeuXV1FRcVh21X3OwbmuwAyVodlREQa2IsvvsjQoUO555576u3WgTosIyIC3HXXXQ22rUsvvfQnX/AGm/bcRaTWnKv2dskSBHX93SrcRaRWEhIS2LFjhwK+Hjj/eO4JCQm1XocOy4hIrWRmZpKbm8u2bdsaZHvFxcV1CrtwEkhf9t+JqbYU7iJSK7GxsbW+S1BtZGdnH7jKNNI1RF90WEZEJAop3EVEopDCXUQkCincRUSikMJdRCQKKdxFRKKQwl1EJAop3EVEopDCXUQkCincRUSiUMSF+77Sct7bUEZ5RWWoSxERCVsRF+7vLN7MP1eUMu7xz/lux+HvXSgi0phFXLiPG9SOSX3jWbO1kDMf+oSZORs15KiIyCEiLtwBjm8Tw7s3n0z/dunc+toSbnhxIXuKykJdlohI2IjIcAdok57IC9cM5dbRPZizdAtnPfQJX2/cFeqyRETCQsSGO4DHY/xqRGde/dUJeDww7rHPefLjdTpMIyKNXkDhbmajzWylma0xs9uqeT3NzN42s2/MbKmZXRX8Ug+vf7t03rnxJE47rgX3zF7OxOcX6DCNiDRqNYa7mXmBR4AzgZ7AJWbW85Bm1wPLnHP9gBHAX8wsLsi1HlFaYiyPXZ7F//yiJx+u2Mq5//iUFVvyG7IEEZGwEcie+xBgjXNunXOuFJgBnHtIGwekmpkBKcBOoDyolQbAzLj6xI7MmHg8+0orGPPIPN76ZlNDlyEiEnJW0/FpMxsLjHbOTfBPjweGOuduqNImFXgL6AGkAhc55/5VzbomAhMBWrZsmTVjxoxaFV1YWEhKSsoR2+wuqeTRRSWs2lXJWR1jGdstFo9ZrbZXnwLpS6RQX8JPtPQD1Jf9Ro4cucA5N6jGhs65I/4AFwJPVZkeD/z9kDZjgb8BBnQB1gNNjrTerKwsV1sffvhhQO1Kyirc/3tjsTv21nfcFdO+dHuKSmu9zfoSaF8igfoSfqKlH86pL/sB810Nue2cC+iwTC7Qrsp0JnDosY6rgNf9217jD/ceAay7XsXFeLh7TB/uPb8P89Zs5/xHP9NVrSLSKAQS7jlAVzPr6P+S9GJ8h2Cq2gicBmBmLYHuwLpgFloXlwxpz/PXDGV7YQnnPjKPL9btCHVJIiL1qsZwd86VAzcAc4DlwMvOuaVmNsnMJvmb/Qk4wcyWAB8AtzrnttdX0bUxrHNT3rx+OE2T4xg/7UveWJgb6pJEROpNTCCNnHOzgdmHzHusyvNNwBnBLS34jm2azOu/Gs6kFxbw25nfsHFHETed1gULwy9aRUTqIqKvUK2NtKRYnr16COcPbMvf3l/Fba8t0fDBIhJ1AtpzjzZxMR7+cmE/MtMTefg/a9hWWMI/Lh1AUlyj/HWISBRqdHvu+5kZvzujO/ec15vslVu59Mkv2bW3NNRliYgERaMN9/0uG3osUy/PYtnmfC58/HM27S4KdUkiInXW6MMdYFSvVjx71RC27Clm7NTPWLetMNQliYjUicLdb1jnpsyYeDzF5ZWMe/xzlm/WoGMiErkU7lX0bpvGy9cNI8bj4aLHP2ehbv4hIhFK4X6ILi1SeGXSMNKT4hg/7StyNuwMdUkiIkdN4V6NdhlJvHzdMFo0ieeKaV/x2ZqwuthWRKRGCvfDaJWWwMyJw2ifkcRV03P4eNW2UJckIhIwhfsRNE+N56WJx9OxWTITnpuvgBeRiKFwr0FGchwvXns8nZunMOG5+XykgBeRCKBwD0BGchwvThhK5+YpTHxuPvN0DF5EwpzCPUDHJMfxzwlD6dA0mWuezdGY8CIS1hTuRyEjOY5/XjuUzGOSuHp6Dgu+02mSIhKeFO5HqVlKPC9eO5SWTRK48ukcvs3bE+qSRER+QuFeCy1SE/jnhKE0SYxl/LQvWbmlINQliYgcROFeS23SE3nx2qHExXi4fNqXuvG2iIQVhXsdHNs0mReuGUp5RSWXPfUlW/YUh7okERFA4V5nXVum8uzVQ9i9r4zLp+mGHyISHhTuQdA3M50nrxjExp37uHJ6DntLykNdkog0cgr3IBnWuSn/uGQAS3J3M+mFBZSW66bbIhI6CvcgOqNXK+47vy+frN7Of73yDZWVLtQliUgjFRPqAqLNuMHt2LG3lPvfXUHzlHju+MVxmFmoyxKRRkbhXg8mndKJrQXFPD1vPS2bxHPdKZ1DXZKINDIK93pgZtzx855sKyjh3n+voEWTeM4bkBnqskSkEVG41xOPx/jLuH7sKCzl968upkVqAsO7NAt1WSLSSOgL1XoUH+PlsfFZdGqWwqTnF7B8c36oSxKRRkLhXs/SEmN55qrBJMfHcNUzOWzeUxTqkkSkEVC4N4A26Yk8c9VgCkvKueqZHAqKy0JdkohEOYV7AzmudRMevWwgq7cWcv2LCymr0EVOIlJ/FO4N6ORuzfnzeb35eNU2/ufNpTini5xEpH7obJkGdtHg9ny3Yx+PZq+lY7MkuoW6IBGJStpzD4FbzujOz/u05t5/r2D+Fg0yJiLBF1C4m9loM1tpZmvM7LbDtBlhZovMbKmZfRTcMqPL/nPg+7dL54nFJSzO3R3qkkQkytQY7mbmBR4BzgR6ApeYWc9D2qQDjwLnOOd6ARfWQ61RJSHWyxPjB5EaZ0x4dr5OkRSRoApkz30IsMY5t845VwrMAM49pM2lwOvOuY0AzrmtwS0zOjVPjee3WQnsK63gmunzNQ68iASN1XTGhpmNBUY75yb4p8cDQ51zN1Rp8yAQC/QCUoGHnHPPVbOuicBEgJYtW2bNmDGjVkUXFhaSkpJSq2XDTWFhIeuKEvjbghL6t/By44B4PBE6imS0vS/R0Jdo6QeoL/uNHDlygXNuUE3tAjlbprqkOfQTIQbIAk4DEoHPzewL59yqgxZy7gngCYBBgwa5ESNGBLD5n8rOzqa2y4ab7OxsbvrFCJq0Wc9dby8jp6Q1t47uEeqyaiXa3pdo6Eu09APUl6MVSLjnAu2qTGcCm6pps905txfYa2YfA/2AVUhAfnlCB1ZvLWRq9lq6NE/hgiyNIikitRfIMfccoKuZdTSzOOBi4K1D2rwJnGRmMWaWBAwFlge31OhmZtx1Ti9O6NyU219fwoLvdoa6JBGJYDWGu3OuHLgBmIMvsF92zi01s0lmNsnfZjnwLrAY+Ap4yjn3bf2VHZ1ivR4evWwgrdMTuO75BeTt1hk0IlI7AZ3n7pyb7Zzr5pzr7Jy7xz/vMefcY1XaTHHO9XTO9XbOPVhfBUe79KQ4pv1yECVllUx4VmfQiEjt6ArVMNSlRSoPXzqAlVvy+a+XdaNtETl6CvcwNbJ7C/77rON4d+kWHv7P6lCXIyIRRgOHhbFrTuzIii0FPPj+arq3TOXMPq1DXZKIRAjtuYcxM+Oe83ozsH06v3v5G5Zt0m36RCQwCvcwt/8+rGmJsVz73Hx2FJaEuiQRiQAK9wjQIjWBJ67IYnthCb/+59e6i5OI1EjhHiH6ZqZz/wV9+XL9Tv749tJQlyMiYU5fqEaQMQPasnxzPo9/vI6erdO4dGj7UJckImFKe+4R5veje3BKt+bc+da3zN+gIQpEpHoK9wjj9RgPXzyAzGOSmPTCAjZpiAIRqYbCPQKlJcXy5BVZFJdVct3zCyguqwh1SSISZhTuEapLi1T+dlF/luTt4fbXl1DTTVdEpHFRuEew03u25Hend+ONhXlM+3R9qMsRkTCicI9wN4zswuherfjz7OV8snpbqMsRkTChcI9wHo/xwLh+dGmRwg0vLmTjjn2hLklEwoDCPQqkxMfw5BW+++Ve+5zGgBcRhXvUOLZpMv+4dACrtxZwyyvf6AtWkUZO4R5FTuranNvO7MG/v93CIx+uCXU5IhJCCvcoc+1JnRjTvw1/eW8V7y/7IdTliEiIKNyjjJlx3wV96dWmCTfPXMSarYWhLklEQkDhHoUSYr08Pn4Q8TEeJj43nz1FZaEuSUQamMI9SrVNT2Tq5Vls3LmP38xYSIVusi3SqCjco9iQjhnceU4vsldu44G5K0Ndjog0II3nHuUuH9qeZZvymZq9lp6tm3B2vzahLklEGoD23KOcmfHHc3ox6NhjmPzqN3ybtyfUJYlIA1C4NwJxMR6mXp7FMUlxTHxuPtt1k22RqKdwbySap8bz5BWD2LG3lF+9sIDSct1kWySaKdwbkd5t0/i/sX3J2bCLO99aqiEKRKKYvlBtZM7t35blmwt47KO19GydyvhhHUJdkojUA+25N0KTR3Xn1B4t+OPby/hs7fZQlyMi9UDh3gh5PcZDF/enQ7Nkrv/n1xoDXiQKKdwbqdSEWJ66YhCVDiY8l0NBsYYoEIkmCvdGrEOzZB69bCBrt+3l5hmLNESBSBRRuDdyw7s0486ze/LBiq1MmaMhCkSihc6WEcYffywrt/jOoOnWMoXzB2aGuiQRqaOA9tzNbLSZrTSzNWZ22xHaDTazCjMbG7wSpb6ZGXed04thnZpy22tLWPDdzlCXJCJ1VGO4m5kXeAQ4E+gJXGJmPQ/T7n5gTrCLlPoX6/Uw9fKBtElPYOJzC8jdpTNoRCJZIHvuQ4A1zrl1zrlSYAZwbjXtbgReA7YGsT5pQOlJcUy7cjClFZVMeHa+zqARiWBW0yXo/kMso51zE/zT44GhzrkbqrRpC7wInApMA95xzr1azbomAhMBWrZsmTVjxoxaFV1YWEhKSkqtlg034diXpdsr+MuCYno383LzwHg8ZgEtF459qa1o6Uu09APUl/1Gjhy5wDk3qKZ2gXyhWt3/7EM/ER4EbnXOVdgRgsA59wTwBMCgQYPciBEjAtj8T2VnZ1PbZcNNOPZlBJCW+R1/mPUtnxS24M6zewW0XDj2pbaipS/R0g9QX45WIOGeC7SrMp0JbDqkzSBghj/YmwFnmVm5c25WUKqUBnf58ceybttenp63nk7NkjUGjUiECSTcc4CuZtYRyAMuBi6t2sA513H/czObju+wjII9wv2/nx/Hdzv2cudbS8nMSGJk9xahLklEAlTjF6rOuXLgBnxnwSwHXnbOLTWzSWY2qb4LlNDxeoyHLxnAca2bcMM/v2bZpvxQlyQiAQroPHfn3GznXDfnXGfn3D3+eY855x6rpu2V1X2ZKpEpOT6Gab8cTGpCLNc8m8OWPcWhLklEAqDhB6RGrdISmHblIPKLyrh6eg6FJeWhLklEaqBwl4D0apPGI5cNZOUPBdzw4teUV+g2fSLhTOEuARvRvQV/Orc32Su3cceb3+o2fSJhTAOHyVG5dGh7cnft49HstbRNT+SGU7uGuiQRqYbCXY7a5FHd2bKnmAfmrqJVWiJjszSKpEi4UbjLUTMz7rugL1sLSrjttcU0T43nlG7Na7WuWQvzmDJnJZt2F9EmPZHJo7ozZkDbIFcs9UXvX/jSMXeplbgY3yiS3Vqm8qsXFvDN97uPeh2zFuZx++tLyNtdhAPydhdx++tLmLUwL/gFS9Dp/QtvCneptdSEWKZfPZimKXFcPT2HLXuP7gyaKXNWUlRWcdC8orIK3REqQuj9C28Kd6mTFqkJPHvVEBzwwPxifsgP/CKnTbuLjmq+hBe9f+FN4S511ql5CtOvGkxhqeOKaV+xZ19g48C3SU88qvkSXvT+hTeFuwRF38x0bhqYwPrte7n62Rz2ldZ8FevkUd1JjPUeNC8x1svkUd3rq0wJIr1/4U3hLkHTs6mXBy/uz8KNu5j0wteUlh/5GPyYAW259/w+tE1PxIC26Ynce34fnW0RIfT+hTedCilBdVaf1tx3fl9+/9pifjtzEQ9fMgCv5/A3cBkzoK3CIILp/QtfCncJunGD25FfXMbd/1pOcryX+87vi+cIAS8iwadwl3ox4aROFBSX89AHq0mM9XLXOb040i0YRSS4FO5Sb27+WVf2lZbz5CfrSYjzctvoHgp4kQaicJd6Y2b891nHsa+0gsc/Wkd8jJffnd4t1GWJNAoKd6lXZsafzu1NaXklD3+wmliPceNpGklSpL4p3KXeeTy+gcYqKh1/eW8VMV4PvxrROdRliUQ1hbs0CK/HmHJhP8orHfe/uwKH49cjuoS6LJGopXCXBuP1GH8d1w8z+L93V+IcXD9SAS9SHxTu0qBivB7+Oq4/hm9UwfIKx02nddFZNCJBpnCXBuf1GH8Z1x+vx8Pf3l9FSXkFk0d1V8CLBJHCXULC6zGmjO1LXIyHR7PXUlxWyR2/OE4BLxIkCncJGY/H+PN5vYmP8fD0vPXsKy3nnvP6HHEsGhEJjMJdQsrMuPPsnqQmxPD3/6yhsKScv47rT1yMBiwVqQuFu4ScmfFfZ3QnJT6Ge/+9goLicqZePpCkOP3zFKkt7R5J2LjulM7cd34fPlm9jUuf/JJde0tDXZJIxFK4S1i5eEh7Hr0si2Wb87nw8c/J0/04RWpF4S5hZ3TvVjx71RB+yC/mvEfmsWxTfqhLEok4CncJS8M6N+WVScPwmDHu8c/5ZPW2UJckElEU7hK2erRqwhvXn0DmMYlc9UwOM3M2hrokkYihcJew1jotkVcmDWNY56bc+toS7n93BZWVLtRliYQ9hbuEvdSEWJ6+cjCXDGnP1Oy1/OqfC9hXWh7qskTCWkDhbmajzWylma0xs9uqef0yM1vs//nMzPoFv1RpzGK9Hv58Xm/u+EVP3lv2A2Onfs4mnUkjclg1hruZeYFHgDOBnsAlZtbzkGbrgVOcc32BPwFPBLtQETPjmhM7Mu3KwXy/cx9n//1Tvlq/M9RliYSlQPbchwBrnHPrnHOlwAzg3KoNnHOfOed2+Se/ADKDW6bIj0Z2b8Eb1w8nLTGWS5/8guc/34BzOg4vUpXV9J/CzMYCo51zE/zT44GhzrkbDtP+FqDH/vaHvDYRmAjQsmXLrBkzZtSq6MLCQlJSUmq1bLhRX2pvX5nj8cUlfLOtguFtYriiVxzx3uAMOhYt70u09APUl/1Gjhy5wDk3qMaGzrkj/gAXAk9VmR4P/P0wbUcCy4GmNa03KyvL1daHH35Y62XDjfpSNxUVle5v7610HW57x43620du/bbCoKw3Wt6XaOmHc+rLfsB8V0O+OucCOiyTC7SrMp0JbDq0kZn1BZ4CznXO7QhgvSJ15vEYN/+sG89cOZjNe4o5+++f8q/Fm0NdlkjIBRLuOUBXM+toZnHAxcBbVRuYWXvgdWC8c25V8MsUObIR3Vvwr5tOpHOLFK5/8WvumPUtxWUVoS5LJGRqDHfnXDlwAzAH3yGXl51zS81skplN8jf7H6Ap8KiZLTKz+fVWschhZB6TxMvXDePakzry/BffMeaReaz+oSDUZYmEREADZjvnZgOzD5n3WJXnE4CffIEq0tBmL9nM7CVbAFj1QwFnPfwJY/q3Zd6a7WzeU0yb9EQmj+rOmAFtg77tWQvzmDJnJZt2F9XrdgLxh1lLeOnL77m5dxnX3D6bS4a24+4xfUJSi4SG7oYgUWPWwjxuf30JRf7DMZUOXIXjlQW5B9rk7S7i9teXAAQ1eA/ddn1tJxB/mLWEF774cRyeCucOTCvgGw8NPyBRY8qclQfCdb/qTvQtKqtgypyV9b7t+thOIF768vujmi/RSeEuUeNohiMI9k1ADrftUAyRUHGYa1cON1+ik8Jdokab9MSA23oM3lyUF7QrWw+37aOpKVi8Vv2FXIebL9FJ4S5RY/Ko7iTGeg+aF+sxYg+5ajU+xkPmMUn8ZsYiJjw7n9xd++pl24mxXiaP6l7ndR+tS4a2O6r5Ep0U7hI1xgxoy73n96FteiIGtE1PZMqF/Zgytt9B8+6/oC8f3jKCP/z8OD5bu4PT//oxT3y8lrKKyqBu+97z+4TkbJm7x/Th8uPbH9hT95px+fHt9WVqI6OzZSSqjBnQttpArW7ehJM6Mbp3K+56ayl/nr2C1xbk8cdze3F8p6ZB3XYo3D2mD3eP6UN2djZrLxsR6nIkBLTnLo1a5jFJPPXLwTwxPovCknIufuILbnppITuLa78XLxIOtOcuApzRqxUndW3O1I/W8thHa3m3spLvvKu47pROJMXpv4lEHu25i/glxnn53end+OB3p9C/hZeHPljNyAeymZmzkQrdt1UijMJd5BDtMpL4df8EXpk0jNZpidz62hLOfOhj3lv2g24KIhFD4S5yGIM7ZPDGr0/g0csGUlbhuPa5+Yx59DM+Xb1dIS9hT+EucgRmxll9WvPeb0/m/gv6sC2/mMunfcm4xz9n3hqFvIQvhbtIAGK8Hi4a3J7/3DKC/z23F9/vLOKyp77kgqmf8Z8VOlwj4UfhLnIUEmK9XDGsA9mTR/Cnc3vxQ34JV0+fz1kPf8obC3PrdCGUSDAp3EVqISHWy3h/yD9wYT/KKir57cxvOPn/PuTxj9ayZ19ZqEuURk4n8IrUQazXw9isTM4f0JbsVVt54uN13PvvFTz4/mouyGrLFcM60K1laqjLlEZI4S4SBB6PcWqPlpzaoyVLN+1h+rwNvDw/lxe+2MjQjhlcfvyxnNGrJfEx3ppXJhIECneRIOvVJo0pF/bj9rOO4+X53/PCF99x40sLyUiO44KBbRk3qB1dtTcv9UzhLlJPMpLjmHRKZyae1IlP12znpa828sy8DTz5yXr6ZaYxNiuTn/dtQ0ZyXKhLlSikcBepZx6PcXK35pzcrTnbC0uYtTCPVxfkcsebS/nj28sY0b05Z/drw8+Oa0lyvP5LSnDoX5JIA2qWEs+Ekzox4aROLN+cz6yFeby5aBPvL99KQqyHU3u04MzerRnZowUpCnqpA/3rEQmR41o34bjWTbh1dA/mf7eLt7/ZxL+/3cLsJVuIi/FwUpdmnN6zJacd15LmqfGhLlcijMJdJMQ8HmNIxwyGdMzgrnN6seC7Xcxespn3lv3AByu2YraEvpnpnNajBSO6N6d3mzQ8Ht0PVY5M4S4SRrxVgv7Os3uyfHMBHyz3hfzf3l/FX99bRUZyHCd2acaJXZtxYpc9/m12AAANA0lEQVRmIbkJt4Q/hbtImDIzerZpQs82TbjxtK5sLyzh09XbyV65lU/XbOetbzYB0LFZMsd3yuD4Tk0Z0jGD1mkKe1G4i0SMZinxB+7T6pxj5Q8FfLp6O5+v3cE732zmpa++B6BdRiKDj81g4LHH4Aoqqah0eHUYp9FRuItEIDOjR6sm9GjVhAkndaK8opLlmwv4asNOvlq/g49Xb+P1hXkA3Jczh76Z6fRrl06/zDT6tkunTVoCZgr8aKZwF4kCMV4PfTLT6JOZxjUndsQ5x8ad+3jh3c8oSWnNwo27mfbpOsoqfEMTH5MUS++2afRs7Tvsc1zrJnRslkysV2MJRguFu0gUMjOObZrM8LaxjBjRG4CS8gpWbC5gce5ulm7K59tNe3hm3gZK/cMUx3k9dGqeTPdWqXRtkULXlql0aZFC+4wkhX4EUriLNBLxMV7foZl26QfmlVVUsm7bXpZvzmfFlgJWbMln/oZdvLlo04E2sV6jfUYSHZul0Kl5Msc2TaJj02TaN02idVqijueHKYW7SCMW6/XQvVUq3VsdPJBZYUk5a7YWsnZrIWu2FbJ+217Wb9/Lx6u3UVpeWWV5o216Iu0yksg8JonMYxJpm55I22MSaZOeSMvUeGK01x8SCncR+YmU+Bj6t0unf5W9fIDKSseW/GI2bN/Lhh37+H7XPjbu2Efurn3M3bSFHXtLD2rvMWiRmkCrtARaNfE9tmgST4vUBFqkxtOiSTzNUuLJSIrThVlBpnAXkYB5PEabdN9e+Qldfvr6vtJyNu0uIndXEZv3FLN5dxGb9hTzQ34xa7YVMm/tdgqKy3+ynNdjZCTH0TQ5jmYp8TRNiSMjOY6MpDgyUuI4JimO73ZU0HJzPsckxZGWGEtCrEdn/ByBwl1EgiYpLoYuLVLp0uLw49XvKy1na34JWwtK2FpQzPaCErYXlrKtoIQde0vZXljCxp372LW3lIKSgz8I7s/55MDzOK+HtKRYmiTEkJYYS5PEWJokxJKaEEPqgccYUuJ//En2//iee0mKi4na7wwCCnczGw08BHiBp5xz9x3yuvlfPwvYB1zpnPs6yLWKRK1ZC/OYMmclm3YX0SY9kcmjuvPK/I3MW7vzQJvhnTO4cFD7n7QDfjJv/nc7eenL77m5dxnX3D6bS4a24+4xfQLa7pgBbQ87P5Dl92+7wjm8Zj/ZdlJcDB2axbDo+9019uXOs3tycrfm7NxXyofzcujQrSe79pXx2drtZK/cxraCEgqKy/B6jLIKx4bte8kvLqeguOzAaZ81iY/xkBTnC/qkOC+JcV4SY70kxXlJiPU9T4jzkhDjJSHWQ0Lsj4/xMR7iY/yPsT8+j9v/4/3xebzXS2yMEev14FxgtdVFjeFuZl7gEeB0IBfIMbO3nHPLqjQ7E+jq/xkKTPU/ikgNZi3M4/bXl1BUVgFA3u4ibp656Cft5q3deVDY5+0uYvKr34CDskp3YN7vZi6isspyFc7xwhcbAQ4K2eq2e/vrS5j/3U5eW5D3k/nAQQFf3fJ12fbkV74B40Ao5+0u4o43l3Lv+X0YM6AtW5p6GdGnNbMW5vHB8q0Hli0uq+T7nUUH2u1XXFZBYUk5hcXlFBSXU1hSzt6ScvaWlrO3pIK9JeXsK61gb2k5+0p9z4tKKw487thbemC6pLyC4rJKisoqqKisezCf2TGWkSPrvJojCmTPfQiwxjm3DsDMZgDnAlXD/VzgOef7OPrCzNLNrLVzbnPQKxaJMlPmrDwQVEerur3TymraAbz05fcHBWx12y0qqziw133o/ClzVh4UntUtX5dtl1UTmoFut7p2vj1sL81SgjtccllFJcVlFZSUV1JS7nte6n/ue6ygpKyS0grf9IHH8krKKnw/nl0bg1pTdaymPw/MbCww2jk3wT89HhjqnLuhSpt3gPucc5/6pz8AbnXOzT9kXROBif7J7sDKWtbdDNhey2XDjfoSnhqsL3GtumTV17or9u3Bm5R2YLp0y5oFddluXZYPwrLNgO1HWrbqNsJcXf59Heuca15To0D23Kv7tuHQT4RA2uCcewJ4IoBtHrkgs/nOuUF1XU84UF/CU7T0xczml+/ZGvH9gOh5T6Bh+hLI1QW5QLsq05nAplq0ERGRBhJIuOcAXc2so5nFARcDbx3S5i3gCvM5Htij4+0iIqFT42EZ51y5md0AzMF3KuTTzrmlZjbJ//pjwGx8p0GuwXcq5FX1VzIQhEM7YUR9CU/R0pdo6QeoL0elxi9URUQk8mhEHxGRKKRwFxGJQmEf7maWYGZfmdk3ZrbUzP7on59hZu+Z2Wr/4zGhrjUQZuY1s4X+awMiuR8bzGyJmS0ys/n+eZHal3Qze9XMVpjZcjMbFol9MbPu/vdj/0++md0coX35rf//+7dm9pI/ByKuHwBm9ht/P5aa2c3+efXel7APd6AEONU51w/oD4z2n5FzG/CBc64r8IF/OhL8BlheZTpS+wEw0jnXv8r5upHal4eAd51zPYB++N6fiOuLc26l//3oD2ThO7nhDSKsL2bWFrgJGOSc643vRI6LibB+AJhZb+BafFf69wN+YWZdaYi+OOci5gdIAr7GN27NSqC1f35rYGWo6wug/kz/G3kq8I5/XsT1w1/rBqDZIfMiri9AE2A9/pMLIrkvh9R/BjAvEvsCtAW+BzLwndH3jr8/EdUPf50X4htscf/0HcDvG6IvkbDnvv9QxiJgK/Cec+5LoKXzn0vvf2wRyhoD9CC+N7bqEByR2A/wXYE818wW+IeVgMjsSydgG/CM/3DZU2aWTGT2paqLgZf8zyOqL865POABYCOwGd91M3OJsH74fQucbGZNzSwJ3ynj7WiAvkREuDvnKpzvT81MYIj/T52IYma/ALY65yJl7IuaDHfODcQ3Iuj1ZnZyqAuqpRhgIDDVOTcA2EsE/Ll/JP6LDc8BXgl1LbXhP/58LtARaAMkm9nloa2qdpxzy4H7gfeAd4FvgJ/eraQeRES47+ec2w1kA6OBH8ysNYD/cWsISwvEcOAcM9sAzABONbMXiLx+AOCc2+R/3IrvuO4QIrMvuUCu/69BgFfxhX0k9mW/M4GvnXM/+KcjrS8/A9Y757Y558qA14ETiLx+AOCcm+acG+icOxnYCaymAfoS9uFuZs3NLN3/PBHfG78C35AHv/Q3+yXwZmgqDIxz7nbnXKZzrgO+P5n/45y7nAjrB4CZJZtZ6v7n+I6HfksE9sU5twX43sy6+2edhm8464jrSxWX8OMhGYi8vmwEjjezJDMzfO/JciKvHwCYWQv/Y3vgfHzvTb33JeyvUDWzvsCz+L4x9wAvO+f+18yaAi8D7fH9Y7jQObfz8GsKH2Y2ArjFOfeLSOyHmXXCt7cOvsMaLzrn7onEvgCYWX/gKSAOWIdv+AwPkdmXJHxfRnZyzu3xz4u498V/yvNF+A5hLAQmAClEWD8AzOwToClQBvzOOfdBQ7wnYR/uIiJy9ML+sIyIiBw9hbuISBRSuIuIRCGFu4hIFFK4i4hEoUBukC3SoPyniX3gn2wFVOAbIgBgiHOuNCSFHYGZXQ3M9p83LxJyOhVSwpqZ3QUUOuceCINavM65isO89ilwg3Nu0VGsL8Y51yCXokvjo8MyElHM7JfmG99/kZk9amYeM4sxs91mNsXMvjazOWY21Mw+MrN1ZnaWf9kJZvaG//WVZvaHANd7t5l9hW9coz+aWY5/fO7HzOcifMNRz/QvH2dmuVWurD7ezN73P7/bzB43s/fwDVYWY2Z/9W97sZlNaPjfqkQjhbtEDP+AcecBJ/gHkovBN5QDQBow1z+YWSlwF77L1i8E/rfKaob4lxkIXGpm/QNY79fOuSHOuc+Bh5xzg4E+/tdGO+dmAouAi5xvPPWaDhsNAM52zo0HJuIbUG4IMBjfIGzta/P7EalKx9wlkvwMXwDO9w05QiK+S+0Bipxz7/mfL8E3TGy5mS0BOlRZxxzn3C4AM5sFnIjv/8Hh1lvKj0MtAJxmZpOBBKAZsAD491H2403nXLH/+RnAcWZW9cOkK75L0kVqTeEukcSAp51zdxw00ywGXwjvV4nvDl77n1f9d37ol0yuhvUWOf8XU/5xW/4BDHTO5ZnZ3fhCvjrl/PiX8aFt9h7Sp1875z5AJIh0WEYiyfvAODNrBr6zampxCOMM890zNQnfmOHzjmK9ifg+LLb7R8W8oMprBUBqlekN+G51xyHtDjUH+LX/g2T/fVATj7JPIj+hPXeJGM65Jf7RAt83Mw++UfYmAZuOYjWfAi8CnYHn95/dEsh6nXM7zOxZfMMbfwd8WeXlZ4CnzKwI33H9u4AnzWwL8NUR6nkc38iAi/yHhLbi+9ARqROdCimNhv9MlN7OuZtDXYtIfdNhGRGRKKQ9dxGRKKQ9dxGRKKRwFxGJQgp3EZEopHAXEYlCCncRkSj0/wHRUJwHFwSFegAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "%matplotlib inline\n",
    "data_pred = pd.DataFrame({'Temperature': np.linspace(start=30, stop=90, num=121), 'Intercept': 1})\n",
    "data_pred['Frequency'] = logmodel.predict(data_pred[['Intercept','Temperature']])\n",
    "data_pred.plot(x=\"Temperature\",y=\"Frequency\",kind=\"line\",ylim=[0,1])\n",
    "plt.scatter(x=data[\"Temperature\"],y=data[\"Frequency\"])\n",
    "plt.grid(True)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "hideCode": false,
    "hidePrompt": false,
    "scrolled": true
   },
   "source": [
    "Comme on pouvait s'attendre au vu des données initiales, la\n",
    "température n'a pas d'impact notable sur la probabilité d'échec des\n",
    "joints toriques. Elle sera d'environ 0.2, comme dans les essais\n",
    "précédents où nous il y a eu défaillance d'au moins un joint. Revenons\n",
    "à l'ensemble des données initiales pour estimer la probabilité de\n",
    "défaillance d'un joint:\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "0.06521739130434782\n"
     ]
    }
   ],
   "source": [
    "data = pd.read_csv(\"shuttle.csv\")\n",
    "print(np.sum(data.Malfunction)/np.sum(data.Count))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 23,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "0.012621523543890678\n"
     ]
    }
   ],
   "source": [
    "p = 1 - (1-0.065**2)**3\n",
    "print(p)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Cette probabilité est donc d'environ $p=0.065$, sachant qu'il existe\n",
    "un joint primaire un joint secondaire sur chacune des trois parties du\n",
    "lançeur, la probabilité de défaillance des deux joints d'un lançeur\n",
    "est de $p^2 \\approx 0.00425$. La probabilité de défaillance d'un des\n",
    "lançeur est donc de $1-(1-p^2)^3 \\approx 1.2%$.  Ça serait vraiment\n",
    "pas de chance... Tout est sous contrôle, le décollage peut donc avoir\n",
    "lieu demain comme prévu.\n",
    "\n",
    "Seulement, le lendemain, la navette Challenger explosera et emportera\n",
    "avec elle ses sept membres d'équipages. L'opinion publique est\n",
    "fortement touchée et lors de l'enquête qui suivra, la fiabilité des\n",
    "joints toriques sera directement mise en cause. Au delà des problèmes\n",
    "de communication interne à la NASA qui sont pour beaucoup dans ce\n",
    "fiasco, l'analyse précédente comporte (au moins) un petit\n",
    "problème... Saurez-vous le trouver ? Vous êtes libre de modifier cette\n",
    "analyse et de regarder ce jeu de données sous tous les angles afin\n",
    "d'expliquer ce qui ne va pas."
   ]
  }
 ],
 "metadata": {
  "celltoolbar": "Hide code",
  "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"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}