diff --git a/module4/.ipynb_checkpoints/challenger-checkpoint.ipynb b/module4/.ipynb_checkpoints/challenger-checkpoint.ipynb new file mode 100644 index 0000000000000000000000000000000000000000..1200c129f650f7c3fc443ec027f30c17b407c0a9 --- /dev/null +++ b/module4/.ipynb_checkpoints/challenger-checkpoint.ipynb @@ -0,0 +1,1092 @@ +{ + "cells": [ + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "# Risk Analysis of the Space Shuttle: Pre-Challenger Prediction of Failure" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "In this document we reperform some of the analysis provided in \n", + "*Risk Analysis of the Space Shuttle: Pre-Challenger Prediction of Failure* by *Siddhartha R. Dalal, Edward B. Fowlkes, Bruce Hoadley* published in *Journal of the American Statistical Association*, Vol. 84, No. 408 (Dec., 1989), pp. 945-957 and available at http://www.jstor.org/stable/2290069. \n", + "\n", + "On the fourth page of this article, they indicate that the maximum likelihood estimates of the logistic regression using only temperature are: $\\hat{\\alpha}=5.085$ and $\\hat{\\beta}=-0.1156$ and their asymptotic standard errors are $s_{\\hat{\\alpha}}=3.052$ and $s_{\\hat{\\beta}}=0.047$. The Goodness of fit indicated for this model was $G^2=18.086$ with 21 degrees of freedom. Our goal is to reproduce the computation behind these values and the Figure 4 of this article, possibly in a nicer looking way." + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Technical information on the computer on which the analysis is run" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "We will be using the python3 language using the pandas, statsmodels, numpy, matplotlib and seaborn libraries." + ] + }, + { + "cell_type": "code", + "execution_count": 1, + "metadata": {}, + "outputs": [], + "source": [ + "#!pip install statsmodels \n", + "#!pip install seaborn" + ] + }, + { + "cell_type": "code", + "execution_count": 2, + "metadata": {}, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "3.7.6 (default, Jan 8 2020, 19:59:22) \n", + "[GCC 7.3.0]\n", + "uname_result(system='Linux', node='H2-SCHMIDER', release='5.3.0-46-generic', version='#38~18.04.1-Ubuntu SMP Tue Mar 31 04:17:56 UTC 2020', machine='x86_64', processor='x86_64')\n", + "IPython 7.13.0\n", + "IPython.core.release 7.13.0\n", + "html 6.0.3\n", + "_csv 1.0\n", + "_ctypes 1.1.0\n", + "_curses b'2.2'\n", + "decimal 1.70\n", + "argparse 1.1\n", + "backcall 0.1.0\n", + "csv 1.0\n", + "ctypes 1.1.0\n", + "cycler 0.10.0\n", + "dateutil 2.8.1\n", + "decimal 1.70\n", + "decorator 4.4.2\n", + "distutils 3.7.6\n", + "ipykernel 5.1.4\n", + "ipykernel._version 5.1.4\n", + "ipython_genutils 0.2.0\n", + "ipython_genutils._version 0.2.0\n", + "jedi 0.16.0\n", + "json 2.0.9\n", + "jupyter_client 6.1.2\n", + "jupyter_client._version 6.1.2\n", + "jupyter_core 4.6.3\n", + "jupyter_core.version 4.6.3\n", + "kiwisolver 1.2.0\n", + "logging 0.5.1.2\n", + "matplotlib 3.2.1\n", + "matplotlib.backends.backend_agg 3.2.1\n", + "numpy 1.18.2\n", + "numpy.core 1.18.2\n", + "numpy.core._multiarray_umath 3.1\n", + "numpy.lib 1.18.2\n", + "numpy.linalg._umath_linalg b'0.1.5'\n", + "pandas 1.0.3\n", + "parso 0.6.2\n", + "patsy 0.5.1\n", + "patsy.version 0.5.1\n", + "pexpect 4.8.0\n", + "pickleshare 0.7.5\n", + "platform 1.0.8\n", + "prompt_toolkit 3.0.4\n", + "ptyprocess 0.6.0\n", + "pygments 2.6.1\n", + "pyparsing 2.4.6\n", + "pytz 2019.3\n", + "re 2.2.1\n", + "scipy 1.4.1\n", + "scipy._lib._uarray 0.5.1+5.ga864a57.scipy\n", + "scipy._lib.decorator 4.0.5\n", + "scipy._lib.six 1.2.0\n", + "scipy.integrate._dop b'$Revision: $'\n", + "scipy.integrate._ode $Id$\n", + "scipy.integrate._odepack 1.9 \n", + "scipy.integrate._quadpack 1.13 \n", + "scipy.integrate.lsoda b'$Revision: $'\n", + "scipy.integrate.vode b'$Revision: $'\n", + "scipy.interpolate._fitpack 1.7 \n", + "scipy.interpolate.dfitpack b'$Revision: $'\n", + "scipy.linalg 0.4.9\n", + "scipy.linalg._fblas b'$Revision: $'\n", + "scipy.linalg._flapack b'$Revision: $'\n", + "scipy.linalg._flinalg b'$Revision: $'\n", + "scipy.ndimage 2.0\n", + "scipy.optimize._cobyla b'$Revision: $'\n", + "scipy.optimize._lbfgsb b'$Revision: $'\n", + "scipy.optimize._minpack 1.10 \n", + "scipy.optimize._nnls b'$Revision: $'\n", + "scipy.optimize._slsqp b'$Revision: $'\n", + "scipy.optimize.minpack2 b'$Revision: $'\n", + "scipy.signal.spline 0.2\n", + "scipy.sparse.linalg.eigen.arpack._arpack b'$Revision: $'\n", + "scipy.sparse.linalg.isolve._iterative b'$Revision: $'\n", + "scipy.special.specfun b'$Revision: $'\n", + "scipy.stats.mvn b'$Revision: $'\n", + "scipy.stats.statlib b'$Revision: $'\n", + "seaborn 0.10.0\n", + "seaborn.external.husl 2.1.0\n", + "six 1.14.0\n", + "statsmodels 0.11.1\n", + "statsmodels.__init__ 0.11.1\n", + "statsmodels.api 0.11.1\n", + "statsmodels.tools.web 0.11.1\n", + "traitlets 4.3.3\n", + "traitlets._version 4.3.3\n", + "urllib.request 3.7\n", + "zlib 1.0\n", + "zmq 18.1.1\n", + "zmq.sugar 18.1.1\n", + "zmq.sugar.version 18.1.1\n" + ] + } + ], + "source": [ + "def print_imported_modules():\n", + " import sys\n", + " for name, val in sorted(sys.modules.items()):\n", + " if(hasattr(val, '__version__')): \n", + " print(val.__name__, val.__version__)\n", + "# else:\n", + "# print(val.__name__, \"(unknown version)\")\n", + "def print_sys_info():\n", + " import sys\n", + " import platform\n", + " print(sys.version)\n", + " print(platform.uname())\n", + "\n", + "import numpy as np\n", + "import pandas as pd\n", + "import matplotlib.pyplot as plt\n", + "import statsmodels.api as sm\n", + "import seaborn as sns\n", + "\n", + "print_sys_info()\n", + "print_imported_modules()" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Loading and inspecting data\n", + "Let's start by reading data." + ] + }, + { + "cell_type": "code", + "execution_count": 3, + "metadata": {}, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "file downloaded on Tue Apr 14 10:06:35 2020\n" + ] + }, + { + "data": { + "text/html": [ + "
\n", + "\n", + "\n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + "
DateCountTemperaturePressureMalfunction
04/12/81666500
111/12/81670501
23/22/82669500
311/11/82668500
44/04/83667500
56/18/82672500
68/30/836731000
711/28/836701000
82/03/846572001
94/06/846632001
108/30/846702001
1110/05/846782000
1211/08/846672000
131/24/856532002
144/12/856672000
154/29/856752000
166/17/856702000
177/2903/856812000
188/27/856762000
1910/03/856792000
2010/30/856752002
2111/26/856762000
221/12/866582001
\n", + "
" + ], + "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/2903/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": 3, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "#data = pd.read_csv(\"https://app-learninglab.inria.fr/moocrr/gitlab/moocrr-session3/moocrr-reproducibility-study/blob/master/data/shuttle.csv\")\n", + "data_url ='https://app-learninglab.inria.fr/moocrr/gitlab/moocrr-session3/moocrr-reproducibility-study/raw/7e01583b99527ad27cbeae0f9d2085fe8f2f1d15/data/shuttle.csv?inline=false'\n", + "data_file='/home/aschmide/Documents/formation_MOOC_RR/mooc-rr/module4/shuttle.csv'\n", + "import os\n", + "import urllib.request\n", + "\n", + "if not os.path.exists(data_file):\n", + " urllib.request.urlretrieve(data_url, data_file)\n", + " \n", + "data = pd.read_csv(data_file)\n", + "\n", + "import time\n", + "modtime=os.path.getmtime(data_file)\n", + "modtime=time.ctime(modtime)\n", + "print('file downloaded on',modtime)\n", + "data" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "We know from our previous experience on this data set that filtering data is a really bad idea. We will therefore process it as such." + ] + }, + { + "cell_type": "code", + "execution_count": 4, + "metadata": {}, + "outputs": [ + { + "data": { + "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYIAAAEKCAYAAAAfGVI8AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4yLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+j8jraAAAYAElEQVR4nO3dfZRddX3v8fdnkiEPJEJMvCMlQbBErlRyEYaHXGqd1IcGu0xq4Sp0IeoV4r2X3Cq6KtTrolyqaxVapVcvVVNKFbQOD7E4rWlRiqPYK5BEQoDQ4BQwTMAAYwgZCMlM5nv/OHv0zHBmZp/J7HPmnN/ntVZWzt77d/b5/s6eM5/ZD+e3FRGYmVm6WupdgJmZ1ZeDwMwscQ4CM7PEOQjMzBLnIDAzS5yDwMwscYUFgaQbJD0j6aExlkvSFyT1SNoq6ZSiajEzs7EVuUfwVWDlOMvPBpZm/9YAXyqwFjMzG0NhQRARPwR+MU6T1cCNUXIPcKSko4qqx8zMKptZx9c+GniybLo3m/f06IaS1lDaa2DOnDmnLlmypCYF5jU0NERLS/OdbmnWfkHz9s39ajy16tujjz76XES8ptKyegZBbhGxDlgH0N7eHps2bapzRSN1d3fT0dFR7zKmXLP2C5q3b+5X46lV3yT9bKxl9YzYnUD5n/aLs3lmZlZD9QyCLuDC7OqhM4E9EfGKw0JmZlaswg4NSfom0AEsktQL/AnQChARXwY2AO8CeoCXgA8VVYuZmY2tsCCIiPMnWB7AJUW9vpmZ5dOcp+HNzCw3B4GZWeIcBGZmiXMQmJklzkFgZpY4B4GZWeIcBGZmiXMQmJklzkFgZpY4B4GZWeIcBGZmiXMQmJklzkFgZpY4B4GZWeIcBGZmiXMQmJklzkFgZpY4B4GZWeIcBGZmiXMQmJklzkFgZpY4B4GZWeIcBGZmiXMQmJklzkFgZpY4B4GZWeIcBGZmiXMQmJklzkFgZpY4B4GZWeIcBGZmiXMQmJklzkFgZpY4B4GZWeIcBGZmiSs0CCStlLRdUo+kyyssP0bS9yXdL2mrpHcVWY+Zmb1SYUEgaQZwHXA2cCJwvqQTRzX7NHBLRLwZOA/4q6LqMTOzyorcIzgd6ImIxyLiANAJrB7VJoBXZY+PAJ4qsB4zM6tAEVHMiqVzgZURcVE2/X7gjIhYW9bmKOC7wALgcODtEbG5wrrWAGsA2traTu3s7Cyk5snq7+9n3rx59S5jyjVrv6B5++Z+NZ5a9W3FihWbI6K90rKZhb/6+M4HvhoRn5O0HLhJ0psiYqi8UUSsA9YBtLe3R0dHR+0rHUd3dzfTraap0Kz9gubtm/vVeKZD34o8NLQTWFI2vTibV+7DwC0AEfFjYDawqMCazMxslCKDYCOwVNJxkg6jdDK4a1SbHcDbACS9kVIQPFtgTWZmNkphQRARg8Ba4A7gEUpXBz0s6SpJq7JmnwAulvQA8E3gg1HUSQszM6uo0HMEEbEB2DBq3hVlj7cBZxVZg5mZjc/fLDYzS5yDwMwscQ4CM7PEOQjMzBLnIDAzS5yDwMwscQ4CM7PEOQjMzBLnIDAzS5yDwMwscQ4CM7PEOQjMzBLnIDAzS5yDwMwscQ4CM7PEOQjMzBLnIDAzS5yDwMwscQ4CM7PEOQjMzBLnIDAzS5yDwMwscQ4CM7PEOQjMzBLnIDAzS5yDwMwscQ4CM7PEOQjMzBLnIDAzS5yDwMwscQ4CM7PEOQjMzBLnIDAzS5yDwMwscYUGgaSVkrZL6pF0+Rht3itpm6SHJf1dkfWYmdkrzSxqxZJmANcB7wB6gY2SuiJiW1mbpcAfA2dFxG5J/6GoeszMrLJcewSSTprEuk8HeiLisYg4AHQCq0e1uRi4LiJ2A0TEM5N4HTMzOwSKiIkbSXcDs4CvAt+IiD05nnMusDIiLsqm3w+cERFry9rcDjwKnAXMAK6MiH+usK41wBqAtra2Uzs7OyfuWQ319/czb968epcx5Zq1X9C8fXO/Gk+t+rZixYrNEdFeaVmuQ0MR8ZbsMM5/BTZLug/424j43iHWNhNYCnQAi4EfSjopIp4f9frrgHUA7e3t0dHRcYgvO7W6u7uZbjVNhWbtFzRv39yvxjMd+pb7ZHFE/BT4NHAZ8FbgC5L+TdLvj/GUncCSsunF2bxyvUBXRAxExOOU9g6W5q3JzMwOXd5zBMskXQs8Avw28O6IeGP2+NoxnrYRWCrpOEmHAecBXaPa3E5pbwBJi4A3AI9V2wkzM5u8vFcNfRG4HvhUROwbnhkRT0n6dKUnRMSgpLXAHZSO/98QEQ9LugrYFBFd2bJ3StoGHAT+KCL6DqE/ZmZWpbxB8LvAvog4CCCpBZgdES9FxE1jPSkiNgAbRs27ouxxAB/P/pmZWR3kPUdwJzCnbHpuNs/MzBpc3iCYHRH9wxPZ47nFlGRmZrWUNwhelHTK8ISkU4F947Q3M7MGkfccwceAWyU9BQh4LfC+wqoyM7OayfuFso2S/iNwQjZre0QMFFeWmZnVSjWDzp0GHJs95xRJRMSNhVRlZmY1kysIJN0E/DqwhdL1/gABOAjMzBpc3j2CduDEyDNCnZmZNZS8Vw09ROkEsZmZNZm8ewSLgG3ZqKP7h2dGxKpCqjIzs5rJGwRXFlmEmZnVT97LR38g6XXA0oi4U9JcSgPJmZlZg8s7DPXFwG3AV7JZR1MaQtrMzBpc3pPFl1C6neQL8Mub1PhG82ZmTSBvEOzPbkAPgKSZlL5HYGZmDS5vEPxA0qeAOZLeAdwK/ENxZZmZWa3kDYLLgWeBB4GPULrZTMU7k5mZWWPJe9XQEPDX2T8zM2sieccaepwK5wQi4vVTXpGZmdVUNWMNDZsN/Bfg1VNfjpmZ1VqucwQR0Vf2b2dE/CWlG9qbmVmDy3to6JSyyRZKewjV3MvAzMymqby/zD9X9ngQeAJ475RXY2ZmNZf3qqEVRRdiZmb1kffQ0MfHWx4Rn5+acszMrNaquWroNKArm343cB/w0yKKMjOz2skbBIuBUyJiL4CkK4HvRMQFRRVmZma1kXeIiTbgQNn0gWyemZk1uLx7BDcC90n6+2z694CvFVOSmZnVUt6rhj4r6Z+At2SzPhQR9xdXlpmZ1UreQ0MAc4EXIuL/AL2SjiuoJjMzq6G8t6r8E+Ay4I+zWa3A14sqyszMaifvHsF7gFXAiwAR8RQwv6iizMysdvIGwYGICLKhqCUdXlxJZmZWS3mD4BZJXwGOlHQxcCe+SY2ZWVOYMAgkCbgZuA1YD5wAXBERX8zx3JWStkvqkXT5OO3OkRSS2sdqY2ZmxZjw8tGICEkbIuIk4Ht5VyxpBnAd8A6gF9goqSsito1qNx/4KHBvVZWbmdmUyHto6CeSTqty3acDPRHxWEQcADqB1RXa/SlwNfByles3M7MpkPebxWcAF0h6gtKVQ6K0s7BsnOccDTxZNt2breeXshveLImI70j6o7FWJGkNsAagra2N7u7unGXXRn9//7SraSo0a7+gefvmfjWe6dC3cYNA0jERsQP4nal+YUktwOeBD07UNiLWAesA2tvbo6OjY6rLOSTd3d1Mt5qmQrP2C5q3b+5X45kOfZtoj+B2SqOO/kzS+og4p4p17wSWlE0vzuYNmw+8CegunY/mtUCXpFURsamK1zEzs0Mw0TkClT1+fZXr3ggslXScpMOA8/jV/QyIiD0RsSgijo2IY4F7AIeAmVmNTRQEMcbjCUXEILAWuAN4BLglIh6WdJWkVdWVaWZmRZno0NB/kvQCpT2DOdlj+NXJ4leN9+SI2ABsGDXvijHaduSq2MzMptS4QRARM2pViJmZ1Uc1w1CbmVkTchCYmSXOQWBmljgHgZlZ4pIJgr7+/Tzw5PP09e+vdylmVoW+/v3sGzjoz26BkgiCb2/ZyVlX38UF19/LWVffRdeWnRM/yczqbviz+/izL/qzW6CmD4K+/v1ctn4rLw8MsXf/IC8PDPHJ9Vv914XZNFf+2T0Y4c9ugZo+CHp376O1ZWQ3W1ta6N29r04VmVke/uzWTtMHweIFcxgYGhoxb2BoiMUL5tSpIjPLw5/d2mn6IFg4bxbXnLOM2a0tzJ81k9mtLVxzzjIWzptV79LMbBzln90Zkj+7Bcp7Y5qGturkoznr+EX07t7H4gVz/INk1iCGP7v3/fhH/Ouq3/RntyBJBAGU/rrwD5FZ41k4bxZzWmf481ugpj80ZGZm43MQmJklzkFgZpY4B4GZWeIcBGZmiXMQmJklzkFgZpY4B4GZWeIcBGZmiXMQmJklzkFgZpY4B4GZWeIcBGZmiXMQmJklzkFgZpY4B4GZWeIcBGZmiXMQmJklzkFgZpY4B4GZWeIcBGZmiSs0CCStlLRdUo+kyyss/7ikbZK2SvoXSa8rsh4zM3ulwoJA0gzgOuBs4ETgfEknjmp2P9AeEcuA24BriqrHzMwqK3KP4HSgJyIei4gDQCewurxBRHw/Il7KJu8BFhdYj5mZVaCIKGbF0rnAyoi4KJt+P3BGRKwdo/3/BX4eEZ+psGwNsAagra3t1M7OzkJqnqz+/n7mzZtX7zKmXLP2C5q3b+5X46lV31asWLE5ItorLZtZ+KvnIOkCoB14a6XlEbEOWAfQ3t4eHR0dtSsuh+7ubqZbTVOhWfsFzds396vxTIe+FRkEO4ElZdOLs3kjSHo78L+At0bE/gLrMTOzCoo8R7ARWCrpOEmHAecBXeUNJL0Z+AqwKiKeKbAWMzMbQ2FBEBGDwFrgDuAR4JaIeFjSVZJWZc3+HJgH3Cppi6SuMVZnZmYFKfQcQURsADaMmndF2eO3F/n6jayvfz+9u/exeMEcFs6bNWVtG0mz9qsoPbv2svulAXp27eX4tvn1LscayLQ4WWwjfXvLTi5bv5XWlhYGhoa45pxlrDr56ENu20iatV9FueL2B7nxnh184qRBLr32h1y4/BiuWn1SvcuyBuEhJqaZvv79XLZ+Ky8PDLF3/yAvDwzxyfVb6et/5Xn0ato2kmbtV1F6du3lxnt2jJh344930LNrb50qskbjIJhmenfvo7Vl5GZpbWmhd/e+Q2rbSJq1X0XZ8uTzVc03G81BMM0sXjCHgaGhEfMGhoZYvGDOIbVtJM3ar6KcvOTIquabjeYgmGYWzpvFNecsY3ZrC/NnzWR2awvXnLOs4snSato2kmbtV1GOb5vPhcuPGTHvwuXH+ISx5eaTxdPQqpOP5qzjF+W6Yqaato2kWftVlKtWn8SFZx7Lg5vv4c5Lz3QIWFUcBNPUwnmzcv/yq6ZtI2nWfhXl+Lb59M5tdQhY1XxoyMwscQ4CM7PEOQjMzBLnIDAzS5yDwMwscQ4CM7PEOQjMzBLnIDAzS5yDwMwscQ4CM7PEOQjMzBLnIDAzS5yDwMwscQ4CM7PEOQjMzBLnIDAzS5yDwMwscQ4CM7PEOQjMzBLnIDAzS5yDwMwscQ4CM7PEOQjMzBLnIDAzS5yDwMwscQ4CM7PEOQjMzBLnIDAzS1yhQSBppaTtknokXV5h+SxJN2fL75V0bJH1mFWrr38/Dzz5PH39+8dtt+nxPj7/3e1serxvytZZbdueXXvZ/dIAPbv2Tti2GkXVW83r7xs4mPs9uG3Tk033HhS5XoCZU77GjKQZwHXAO4BeYKOkrojYVtbsw8DuiDhe0nnA1cD7iqrJrBrf3rKTy9ZvpbWlhYGhIa45ZxmrTj76Fe0uuP4eftRTCoAv3NXDW45fyE0XnXlI66y27RW3P8iN9+zgEycNcum1P+TC5cdw1eqTJtnz4uut9vX/8I0DXHr1Xbneg2HN8h4Uud5hRe4RnA70RMRjEXEA6ARWj2qzGvha9vg24G2SVGBNZrn09e/nsvVbeXlgiL37B3l5YIhPrt/6ir/GNj3e98sQGHZ3T1/FPYO866y2bc+uvSN+AQLc+OMdh/xXcVH1Tub1D0Yk+R4Uud5yiogpW9mIFUvnAisj4qJs+v3AGRGxtqzNQ1mb3mz637M2z41a1xpgTTZ5ArC9kKInbxHw3IStGk+z9gsm6Jta58ydueCoN6ilZcbwvBgaOji4++lHY2DfS8PzZsxf9GszDj/yqNHPP/ji808f3PvcU5NZZ7VtW+YesXDmq15zLMDBl/YwY+4RAAy+8OwTQy/tmfhY1SG+B9W2nczrD/crz3tQrkHegyn5WczhdRHxmkoLCjs0NJUiYh2wrt51jEXSpohor3cdU61Z+wXN2zdJmwb3PON+NZDp8LNY5KGhncCSsunF2byKbSTNBI4AJp3eZmZWvSKDYCOwVNJxkg4DzgO6RrXpAj6QPT4XuCuKOlZlZmYVFXZoKCIGJa0F7gBmADdExMOSrgI2RUQX8DfATZJ6gF9QCotGNG0PWx2iZu0XNG/f3K/GU/e+FXay2MzMGoO/WWxmljgHgZlZ4hwEVZL0hKQHJW2RtCmbd6Wkndm8LZLeVe86J0PSkZJuk/Rvkh6RtFzSqyV9T9JPs/8X1LvOao3Rr4bfZpJOKKt/i6QXJH2s0bfZOP1qhm12qaSHJT0k6ZuSZmcX1NybDbVzc3ZxTW3r8jmC6kh6Amgv/9KbpCuB/oj4i3rVNRUkfQ24OyKuz34Y5wKfAn4REX+WjRe1ICIuq2uhVRqjXx+jCbbZsGxIl53AGcAlNPg2GzaqXx+igbeZpKOBHwEnRsQ+SbcAG4B3Ad+KiE5JXwYeiIgv1bI27xEYAJKOAH6L0pVcRMSBiHiekcOAfA34vfpUODnj9KvZvA3494j4GQ2+zUYp71czmAnMyb43NRd4GvhtSkPsQJ22l4OgegF8V9LmbOiLYWslbZV0Q6PtimeOA54F/lbS/ZKul3Q40BYRT2dtfg601a3CyRmrX9D426zcecA3s8eNvs3KlfcLGnibRcRO4C+AHZQCYA+wGXg+IgazZr3A1I0ml5ODoHq/GRGnAGcDl0j6LeBLwK8DJ1PawJ+rY32TNRM4BfhSRLwZeBEYMXR49mW/RjuWOFa/mmGbAZAd7loF3Dp6WYNuM6Bivxp6m2XBtZrSHye/BhwOrKxrURkHQZWyVCcingH+Hjg9InZFxMGIGAL+mtLIq42mF+iNiHuz6dso/QLdJekogOz/Z+pU32RV7FeTbLNhZwM/iYhd2XSjb7NhI/rVBNvs7cDjEfFsRAwA3wLOAo7MDhVB5aF4CucgqIKkwyXNH34MvBN4aPhDl3kP8FA96jsUEfFz4ElJJ2Sz3gZsY+QwIB8Avl2H8iZtrH41wzYrcz4jD5809DYrM6JfTbDNdgBnSporSfzqM/Z9SkPsQJ22l68aqoKk11PaC4DSIYe/i4jPSrqJ0u5qAE8AHyk7RtswJJ0MXA8cBjxG6SqNFuAW4BjgZ8B7I+IXdStyEsbo1xdojm12OKVfMK+PiD3ZvIU0/jar1K+G/5xJ+t+Ubr41CNwPXETpnEAn8Ops3gURMfW3IRuvLgeBmVnafGjIzCxxDgIzs8Q5CMzMEucgMDNLnIPAzCxxDXHzerO8sksn/yWbfC1wkNIQE1D68t+BuhRWgaQO4EBE/L9612JpcxBYU4mIPkrXmk+LUWElzSwbR2a0DqAfyB0EE6zPbFJ8aMianqRTJf0gGyjwjrLhF7olXStpU3afgtMkfSsbx/8zWZtjs/sYfCNrc5ukuTnW+5cq3a/io5LenY03f7+kOyW1SToW+G/ApdnY+m+R9FVJ55bV3Z/93yHpbkldlL4VPUPSn0vamA3A9pFavp/WfBwE1uwEfBE4NyJOBW4APlu2/EBEtANfpvTV/kuANwEfzA4zAZwA/FVEvBF4AfgfklonWO9hEdEeEZ+jNAb9mdmgd53AJyPiiew1r42IkyPi7gn6cQrw0Yh4A/BhYE9EnAacBlws6bjq3xqzEh8asmY3i9Iv9u+VhndhBqWRK4d1Zf8/CDw8PGSBpMeAJcDzwJMR8a9Zu68Dfwj88wTrvbns8WLg5myP4TDg8Un0476IGH7eO4FlZXsPRwBLJ7leMweBNT1R+gW/fIzlw2O6DJU9Hp4e/nyMHoclcqz3xbLHXwQ+HxFd2QniK8d4ziDZXrqkFkqhUWl9Av5nRNwxxnrMquJDQ9bs9gOvkbQcQFKrpN+och3HDD8f+ANKh3q2V7HeI/jV0MIfKJu/F5hfNv0EcGr2eBXQOsb67gD+e3Z4CklvKLvZjlnVHATW7IYoDfF7taQHgC3Af65yHdsp3YToEWABpZvcHKhivVcCt0raDDxXNv8fgPcMnyymNMb+W7P1LWfkXkC56ykNX/wTSQ8BX8F793YIPPqo2Tiyq3v+MSLeVOdSzArjPQIzs8R5j8DMLHHeIzAzS5yDwMwscQ4CM7PEOQjMzBLnIDAzS9z/Bxqb4z828mHUAAAAAElFTkSuQmCC\n", + "text/plain": [ + "
" + ] + }, + "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\",ylim=[0,1])\n", + "plt.grid(True)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Logistic regression\n", + "\n", + "Let's assume O-rings independently fail with the same probability which solely depends on temperature. A logistic regression should allow us to estimate the influence of temperature." + ] + }, + { + "cell_type": "code", + "execution_count": 5, + "metadata": {}, + "outputs": [ + { + "name": "stderr", + "output_type": "stream", + "text": [ + "/home/aschmide/miniconda3/lib/python3.7/site-packages/ipykernel_launcher.py:7: DeprecationWarning: Calling Family(..) with a link class as argument is deprecated.\n", + "Use an instance of a link class instead.\n", + " import sys\n" + ] + }, + { + "data": { + "text/html": [ + "\n", + "\n", + "\n", + " \n", + "\n", + "\n", + " \n", + "\n", + "\n", + " \n", + "\n", + "\n", + " \n", + "\n", + "\n", + " \n", + "\n", + "\n", + " \n", + "\n", + "\n", + " \n", + "\n", + "\n", + " \n", + "\n", + "\n", + " \n", + "\n", + "
Generalized Linear Model Regression Results
Dep. Variable: Frequency No. Observations: 23
Model: GLM Df Residuals: 21
Model Family: Binomial Df Model: 1
Link Function: logit Scale: 1.0000
Method: IRLS Log-Likelihood: -3.9210
Date: Tue, 14 Apr 2020 Deviance: 3.0144
Time: 10:27:51 Pearson chi2: 5.00
No. Iterations: 6
Covariance Type: nonrobust
\n", + "\n", + "\n", + " \n", + "\n", + "\n", + " \n", + "\n", + "\n", + " \n", + "\n", + "
coef std err z P>|z| [0.025 0.975]
Intercept 5.0850 7.477 0.680 0.496 -9.570 19.740
Temperature -0.1156 0.115 -1.004 0.316 -0.341 0.110
" + ], + "text/plain": [ + "\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: Tue, 14 Apr 2020 Deviance: 3.0144\n", + "Time: 10:27:51 Pearson chi2: 5.00\n", + "No. Iterations: 6 \n", + "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": 5, + "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']], \n", + " family=sm.families.Binomial(sm.families.links.logit)).fit()\n", + "\n", + "logmodel.summary()" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "The maximum likelyhood estimator of the intercept and of Temperature are thus $\\hat{\\alpha}=5.0849$ and $\\hat{\\beta}=-0.1156$. This **corresponds** to the values from the article of Dalal *et al.* The standard errors are $s_{\\hat{\\alpha}} = 7.477$ and $s_{\\hat{\\beta}} = 0.115$, which is **different** from the $3.052$ and $0.04702$ reported by Dallal *et al.* The deviance is $3.01444$ with 21 degrees of freedom. I cannot find any value similar to the Goodness of fit ($G^2=18.086$) reported by Dalal *et al.* There seems to be something wrong. Oh I know, I haven't indicated that my observations are actually the result of 6 observations for each rocket launch. Let's indicate these weights (since the weights are always the same throughout all experiments, it does not change the estimates of the fit but it does influence the variance estimates)." + ] + }, + { + "cell_type": "code", + "execution_count": 6, + "metadata": {}, + "outputs": [ + { + "name": "stderr", + "output_type": "stream", + "text": [ + "/home/aschmide/miniconda3/lib/python3.7/site-packages/ipykernel_launcher.py:2: DeprecationWarning: Calling Family(..) with a link class as argument is deprecated.\n", + "Use an instance of a link class instead.\n", + " \n" + ] + }, + { + "data": { + "text/html": [ + "\n", + "\n", + "\n", + " \n", + "\n", + "\n", + " \n", + "\n", + "\n", + " \n", + "\n", + "\n", + " \n", + "\n", + "\n", + " \n", + "\n", + "\n", + " \n", + "\n", + "\n", + " \n", + "\n", + "\n", + " \n", + "\n", + "\n", + " \n", + "\n", + "
Generalized Linear Model Regression Results
Dep. Variable: Frequency No. Observations: 23
Model: GLM Df Residuals: 21
Model Family: Binomial Df Model: 1
Link Function: logit Scale: 1.0000
Method: IRLS Log-Likelihood: -23.526
Date: Tue, 14 Apr 2020 Deviance: 18.086
Time: 10:27:51 Pearson chi2: 30.0
No. Iterations: 6
Covariance Type: nonrobust
\n", + "\n", + "\n", + " \n", + "\n", + "\n", + " \n", + "\n", + "\n", + " \n", + "\n", + "
coef std err z P>|z| [0.025 0.975]
Intercept 5.0850 3.052 1.666 0.096 -0.898 11.068
Temperature -0.1156 0.047 -2.458 0.014 -0.208 -0.023
" + ], + "text/plain": [ + "\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: -23.526\n", + "Date: Tue, 14 Apr 2020 Deviance: 18.086\n", + "Time: 10:27:51 Pearson chi2: 30.0\n", + "No. Iterations: 6 \n", + "Covariance Type: nonrobust \n", + "===============================================================================\n", + " coef std err z P>|z| [0.025 0.975]\n", + "-------------------------------------------------------------------------------\n", + "Intercept 5.0850 3.052 1.666 0.096 -0.898 11.068\n", + "Temperature -0.1156 0.047 -2.458 0.014 -0.208 -0.023\n", + "===============================================================================\n", + "\"\"\"" + ] + }, + "execution_count": 6, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "logmodel=sm.GLM(data['Frequency'], data[['Intercept','Temperature']], \n", + " family=sm.families.Binomial(sm.families.links.logit),\n", + " var_weights=data['Count']).fit()\n", + "\n", + "logmodel.summary()" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Good, now I have recovered the asymptotic standard errors $s_{\\hat{\\alpha}}=3.052$ and $s_{\\hat{\\beta}}=0.047$.\n", + "The Goodness of fit (Deviance) indicated for this model is $G^2=18.086$ with 21 degrees of freedom (Df Residuals).\n", + "\n", + "**I have therefore managed to fully replicate the results of the Dalal *et al.* article**." + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Predicting failure probability\n", + "The temperature when launching the shuttle was 31°F. Let's try to estimate the failure probability for such temperature using our model.:" + ] + }, + { + "cell_type": "code", + "execution_count": 7, + "metadata": {}, + "outputs": [ + { + "data": { + "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXQAAAEKCAYAAAACS67iAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4yLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+j8jraAAAZqUlEQVR4nO3df3TV9Z3n8efbECT8ECwqIwaF2SKOg/Ij4VdxbGirYNtB3GFE6tDaU0p3t7Q6jsyRM93qOHrOduOOzjqOI6usnXU1IMdJsYfT0Dpk7bqrBgqIwAbQppJoB8XyIzZICO/94/u96SUkuTc39+be++H1OCcn9/u9n+/3+3nnS15887nf+7nm7oiISPE7L98dEBGR7FCgi4gEQoEuIhIIBbqISCAU6CIigVCgi4gEImWgm9laMztkZm/18LyZ2X81swNm9qaZTc9+N0VEJJV0rtCfARb08vxNwMT4awXwRP+7JSIifZUy0N39FeCjXprcDPyTR14DRpnZpdnqoIiIpGdQFvZxGXAwabk5Xvd+14ZmtoLoKp6ysrKKcePG9flgH51wPulwLLO+FhwH1VJgQqkDVEuhKj0PRpdl9hLmvn37PnT3i7t7LhuBnjZ3XwOsAaisrPStW7dmtJ/6+nqqqqqy2LP8US2FJ5Q6QLUUqv7UYma/6um5bNzl0gIkX2qXx+tERGQAZSPQNwJfje92mQ0cdfezhltERCS3Ug65mNnzQBVwkZk1A/cBpQDu/o/AJuCLwAHgt8DXc9VZERHpWcpAd/elKZ534NtZ65GIFI329naam5s5ceJEzo81cuRI9u7dm/PjDIR0ahkyZAjl5eWUlpamvd8BfVFURMLS3NzMiBEjGD9+PGa5vQfl+PHjjBgxIqfHGCipanF3Dh8+THNzMxMmTEh7v3rrv4hk7MSJE4wePTrnYX6uMTNGjx7d5798FOgi0i8K89zI5OeqQBcRCYTG0EWkqJWUlHDNNdd0LtfW1jJ+/Pj8dSiPFOgiUtTKysrYsWNHt8+5O+7OeeedG4MR50aVInLOaGpqYtKkSXz1q19l8uTJHDx4kOrqambMmMG1117Lfffd19n2oYce4sorr+S6665j6dKlPPzwwwBUVVWRmJrkww8/7Lzi7+joYNWqVZ37evLJJ4HfvZV/8eLFXHXVVdx+++1Ed3RDQ0MDn/nMZ5gyZQozZ87k+PHjLFiw4Iz/hK677jp27tzZ79p1hS4iWfHXL+1mz3vHsrrPq8dewH1//Ie9tmlra2Pq1KkATJgwgUceeYT9+/fzwx/+kNmzZ7N582b279/PG2+8gbuzcOFCXnnlFYYNG0ZNTQ07duzg1KlTTJ8+nYqKil6P9fTTTzNy5EgaGhr45JNPmDt3LjfeeCMA27dvZ/fu3YwdO5a5c+fy6quvMnPmTJYsWcK6deuYMWMGx44do6ysjGXLlvHMM8/w6KOPsm/fPk6cOMGUKVP6/fNSoItIUes65NLU1MQVV1zB7NmzAdi8eTObN29m2rRpALS2trJ//36OHz/OLbfcwtChQwFYuHBhymNt3ryZN998kw0bNgBw9OhR9u/fz+DBg5k5cybl5eUATJ06laamJkaOHMmll17KjBkzALjgggsAuOWWW5g7dy7V1dWsXbuWO+64Iys/CwW6iGRFqivpgTRs2LDOx+7O6tWr+da3vnVGm0cffbTH7QcNGsTp06cBzrgX3N157LHHmD9//hnt6+vrOf/88zuXS0pKOHXqVI/7Hzp0KDfccAM/+tGPWL9+Pdu2bUuvsBQ0hi4iQZs/fz5r166ltbUVgJaWFg4dOsT1119PbW0tbW1tHD9+nJdeeqlzm/Hjx3eGbOJqPLGvJ554gvb2dgD27dvHxx9/3OOxJ02axPvvv09DQwMQvUM0EfTLly/nu9/9LjNmzODCCy/MSq26QheRoN14443s3buXOXPmADB8+HCeffZZpk+fzpIlS5gyZQqXXHJJ57AIwD333MOtt97KmjVr+NKXvtS5fvny5TQ1NTF9+nTcnYsvvpja2toejz148GDWrVvHd77zHdra2igrK+NnP/sZABUVFVxwwQV8/etZnM8wcVvPQH9VVFR4prZs2ZLxtoVGtRSeUOpwz30te/bsyen+kx07diyn+7/vvvu8uro6p8dIOHbsmLe0tPjEiRO9o6Ojx3bd/XyBrd5DrmrIRURkgD333HPMmjWLhx56KKv3yGvIRUQEuP/++wfsWF/5ylfOepE2G3SFLiL94vEbaCS7Mvm5KtBFJGNDhgzh8OHDCvUs83g+9CFDhvRpOw25iEjGysvLaW5u5oMPPsj5sU6cONHngCtU6dSS+MSivlCgi0jGSktL+/SJOv1RX1/f+W7PYperWjTkIiISCAW6iEggFOgiIoFQoIuIBEKBLiISCAW6iEggFOgiIoFQoIuIBEKBLiISCAW6iEggFOgiIoFQoIuIBEKBLiISCAW6iEggFOgiIoFQoIuIBCKtQDezBWbWaGYHzOzebp6/3My2mNl2M3vTzL6Y/a6KiEhvUga6mZUAjwM3AVcDS83s6i7Nvgesd/dpwG3AP2S7oyIi0rt0rtBnAgfc/R13PwnUADd3aePABfHjkcB72euiiIikw1J9WreZLQYWuPvyeHkZMMvdVya1uRTYDFwIDAO+4O7butnXCmAFwJgxYypqamoy6nRrayvDhw/PaNtCo1oKTyh1gGopVP2pZd68edvcvbLbJ9291y9gMfBU0vIy4O+7tLkb+Iv48RxgD3Beb/utqKjwTG3ZsiXjbQuNaik8odThrloKVX9qAbZ6D7mazpBLCzAuabk8XpfsG8D6+D+I/wsMAS5KY98iIpIl6QR6AzDRzCaY2WCiFz03dmnzLvB5ADP7A6JA/yCbHRURkd6lDHR3PwWsBOqAvUR3s+w2swfMbGHc7C+Ab5rZTuB54I74TwMRERkgg9Jp5O6bgE1d1n0/6fEeYG52uyYiIn2hd4qKiARCgS4iEggFuohIIBToIiKBUKCLiARCgS4iEggFuohIIBToIiKBUKCLiARCgS4iEggFuohIIBToIiKBUKCLiARCgS4iEggFuohIIBToIiKBUKCLiARCgS4iEggFuohIIBToIiKBUKCLiARCgS4iEggFuohIIBToIiKBUKCLiARCgS4iEggFuohIIBToIiKBUKCLiARCgS4iEggFuohIIBToIiKBUKCLiARCgS4iEoi0At3MFphZo5kdMLN7e2hzq5ntMbPdZvZcdrspIiKpDErVwMxKgMeBG4BmoMHMNrr7nqQ2E4HVwFx3/42ZXZKrDouISPfSuUKfCRxw93fc/SRQA9zcpc03gcfd/TcA7n4ou90UEZFUzN17b2C2GFjg7svj5WXALHdfmdSmFtgHzAVKgPvd/Sfd7GsFsAJgzJgxFTU1NRl1urW1leHDh2e0baFRLYUnlDpAtRSq/tQyb968be5e2d1zKYdc0jQImAhUAeXAK2Z2jbsfSW7k7muANQCVlZVeVVWV0cHq6+vJdNtCo1oKTyh1gGopVLmqJZ0hlxZgXNJyebwuWTOw0d3b3f2XRFfrE7PTRRERSUc6gd4ATDSzCWY2GLgN2NilTS3R1TlmdhFwJfBOFvspIiIppAx0dz8FrATqgL3AenffbWYPmNnCuFkdcNjM9gBbgFXufjhXnRYRkbOlNYbu7puATV3WfT/psQN3x18iIpIHeqeoiEggFOgiIoFQoIuIBEKBLiISCAW6iEggFOgiIoFQoIuIBEKBLiISCAW6iEggFOgiIoFQoIuIBEKBLiISCAW6iEggFOgiIoFQoIuIBEKBLiISCAW6iEggFOgiIoFQoIuIBEKBLiISCAW6iEggFOgiIoFQoIuIBEKBLiISCAW6iEggFOgiIoFQoIuIBEKBLiISCAW6iEggFOgiIoFQoIuIBEKBLiISCAW6iEggFOgiIoFQoIuIBCKtQDezBWbWaGYHzOzeXtr9iZm5mVVmr4siIpKOlIFuZiXA48BNwNXAUjO7upt2I4A7gdez3UkREUktnSv0mcABd3/H3U8CNcDN3bT7G+AHwIks9k9ERNJk7t57A7PFwAJ3Xx4vLwNmufvKpDbTgb9y9z8xs3rgHnff2s2+VgArAMaMGVNRU1OTUadbW1sZPnx4RtsWGtVSeEKpA1RLoepPLfPmzdvm7t0Oaw/qV68AMzsP+FvgjlRt3X0NsAagsrLSq6qqMjpmfX09mW5baFRL4QmlDlAthSpXtaQz5NICjEtaLo/XJYwAJgP1ZtYEzAY26oVREZGBlU6gNwATzWyCmQ0GbgM2Jp5096PufpG7j3f38cBrwMLuhlxERCR3Uga6u58CVgJ1wF5gvbvvNrMHzGxhrjsoIiLpSWsM3d03AZu6rPt+D22r+t8tERHpK71TVEQkEAp0EZFAKNBFRAKhQBcRCYQCXUQkEP1+p6hIf9Rub6G6rpH3jrQxdlQZq+ZPYtG0y/LdLUmTzl9hUaBL3tRub2H1i7toa+8AoOVIG6tf3AWgUCgCOn+FR0MukjfVdY2dYZDQ1t5BdV1jnnokfaHzV3gU6JI37x1p69N6KSw6f4VHgS55M3ZUWZ/WS2HR+Ss8CnTJm1XzJ1FWWnLGurLSElbNn5SnHklf6PwVHr0oKnmTeOFMd0kUJ52/wqNAl7xaNO0yBUAR0/krLBpyEREJhAJdRCQQCnQRkUAo0EVEAqFAFxEJhAJdRCQQCnQRkUAo0EVEAqFAFxEJhAJdRCQQCnQRkUAo0EVEAqFAFxEJhAJdRCQQCnQRkUAo0EVEAqFAFxEJhAJdRCQQCnQRkUAo0EVEAqFAFxEJRFqBbmYLzKzRzA6Y2b3dPH+3me0xszfN7GUzuyL7XRURkd6kDHQzKwEeB24CrgaWmtnVXZptByrd/VpgA/Cfs91RERHpXTpX6DOBA+7+jrufBGqAm5MbuPsWd/9tvPgaUJ7dboqISCrm7r03MFsMLHD35fHyMmCWu6/sof3fA7929we7eW4FsAJgzJgxFTU1NRl1urW1leHDh2e0baFRLYUnlDpAtRSq/tQyb968be5e2d1zg/rVqy7M7M+ASuCz3T3v7muANQCVlZVeVVWV0XHq6+vJdNtCo1oKTyh1gGopVLmqJZ1AbwHGJS2Xx+vOYGZfAP4K+Ky7f5Kd7omISLrSGUNvACaa2QQzGwzcBmxMbmBm04AngYXufij73RQRkVRSBrq7nwJWAnXAXmC9u+82swfMbGHcrBoYDrxgZjvMbGMPuxMRkRxJawzd3TcBm7qs+37S4y9kuV8iGand3kJ1XSPvHWlj7KgyVs2fBHDWukXTLhuQY+fiOOn4Xu0unn/9IHdNbucbqzexdNY4Hlx0TV76IgMnqy+KiuRT7fYWVr+4i7b2DgBajrSx6oWdYNDe4Z3rVr+4CyCrYdvdsXNxnHR8r3YXz772budyh3vnskI9bHrrvwSjuq6xM1AT2k97Z5gntLV3UF3XmPNj5+I46Xj+9YN9Wi/hUKBLMN470paTtv3ZX7aPk46OHt5b0tN6CYcCXYIxdlRZTtr2Z3/ZPk46Ssz6tF7CoUCXYKyaP4my0pIz1pWeZ5SWnBlkZaUlnS+W5vLYuThOOpbOGten9RIOvSgqwUi8+JiPu1x6OnY+7nJJvPCZGDMvMdNdLucIBboEZdG0y7oN0YEI1p6OnQ8PLrqGBxddQ319PW/fXpXv7sgA0ZCLiEggFOgiIoFQoIuIBEKBLiISCAW6iEggFOgiIoFQoIuIBEKBLiISCAW6iEggFOgiIoFQoIuIBEKBLiISCAW6iEggFOgiIoFQoIuIBEKBLiISCAW6iEggFOgiIoFQoIuIBEKBLiISCAW6iEggFOgiIoFQoIuIBEKBLiISCAW6iEggFOgiIoFQoIuIBEKBLiISiEHpNDKzBcDfASXAU+7+n7o8fz7wT0AFcBhY4u5N2e2qSLhqt7dQXdfIe0faGDuqjFXzJ/HC1nd59e2POtvM/Tef4k8rLz+rHXDWuq2/+ojnXz/IXZPb+cbqTSydNY4HF12T1nG729+iaZel3e/EsTvcKTHLybG727anPp5LUga6mZUAjwM3AM1Ag5ltdPc9Sc2+AfzG3T9tZrcBPwCW5KLDIqGp3d7C6hd30dbeAUDLkTbuWrfjrHavvv3RGQHfcqSNVRt2gkP7ae9cd/e6HZxO2q7DnWdfexfgjGDt7rirXtgJBu0dv9vf6hd3AZwVmN1tPxDH7m7bnvp4rklnyGUmcMDd33H3k0ANcHOXNjcDP4wfbwA+b2aWvW6KhKu6rrEznPqqvcM7wzzhdA9tn3/9YMrjtp/2zkBNaGvvoLqu8az9dbf9QBy7u2176uO5xty99wZmi4EF7r48Xl4GzHL3lUlt3orbNMfLb8dtPuyyrxXAinhxEpDpGbgI+DBlq+KgWgrPgNYx+Pc+XZGrfXf89iglQ0d2Lp/89YFtmR43edv+bp/hthcBH/a2bdc+FrD+/Bu7wt0v7u6JtMbQs8Xd1wBr+rsfM9vq7pVZ6FLeqZbCE0odENVy6uihYGoJ6bzkopZ0hlxagHFJy+Xxum7bmNkgYCTRi6MiIjJA0gn0BmCimU0ws8HAbcDGLm02Al+LHy8G/sVTjeWIiEhWpRxycfdTZrYSqCO6bXGtu+82sweAre6+EXga+B9mdgD4iCj0c6nfwzYFRLUUnlDqANVSqHJSS8oXRUVEpDjonaIiIoFQoIuIBKLgA93MhpjZG2a208x2m9lfx+snmNnrZnbAzNbFL9gWPDMrMbPtZvbjeLlY62gys11mtsPMtsbrPmVmPzWz/fH3C/Pdz3SY2Sgz22Bm/8/M9prZnGKsxcwmxecj8XXMzO4q0lr+PP59f8vMno9zoFh/V+6M69htZnfF63JyTgo+0IFPgM+5+xRgKrDAzGYTTS/wiLt/GvgN0fQDxeBOYG/ScrHWATDP3acm3U97L/Cyu08EXo6Xi8HfAT9x96uAKUTnp+hqcffG+HxMJZpX6bfAP1NktZjZZcB3gUp3n0x0M0ZiSpGi+l0xs8nAN4necT8F+LKZfZpcnRN3L5ovYCjwC2AW0busBsXr5wB1+e5fGv0vj0/e54AfA1aMdcR9bQIu6rKuEbg0fnwp0JjvfqZRx0jgl8Q3CBRzLV36fyPwajHWAlwGHAQ+RXQn3o+B+cX4uwL8KfB00vJ/BP4yV+ekGK7QE8MUO4BDwE+Bt4Ej7n4qbtJM9I+g0D1KdDITU16MpjjrAHBgs5lti6d0ABjj7u/Hj38NjMlP1/pkAvAB8N/jobCnzGwYxVlLstuA5+PHRVWLu7cADwPvAu8DR4FtFOfvylvAH5nZaDMbCnyR6E2YOTknRRHo7t7h0Z+R5UR/ulyV5y71mZl9GTjk7sUy10Qq17n7dOAm4Ntmdn3ykx5dehTDPbGDgOnAE+4+DfiYLn/+FlEtAMRjywuBF7o+Vwy1xOPJNxP9ZzsWGAYsyGunMuTue4mGijYDPwF2AB1d2mTtnBRFoCe4+xFgC9GfW6PiaQag++kICs1cYKGZNRHNWPk5orHbYqsD6LyKwt0PEY3TzgT+1cwuBYi/H8pfD9PWDDS7++vx8gaigC/GWhJuAn7h7v8aLxdbLV8AfunuH7h7O/Ai0e9Psf6uPO3uFe5+PdHY/z5ydE4KPtDN7GIzGxU/LiOal30vUbAvjpt9DfhRfnqYHndf7e7l7j6e6M/hf3H32ymyOgDMbJiZjUg8JhqvfYszp4Aoilrc/dfAQTObFK/6PLCHIqwlyVJ+N9wCxVfLu8BsMxsaT8OdOCdF97sCYGaXxN8vB/4t8Bw5OicF/05RM7uWaK71EqL/gNa7+wNm9vtEV7qfArYDf+bun+Svp+kzsyrgHnf/cjHWEff5n+PFQcBz7v6QmY0G1gOXA78CbnX3j3rYTcEws6nAU8Bg4B3g68T/1ii+WoYRBeLvu/vReF3RnZf49uQlwCmi34vlRGPmRfW7AmBmPyd6vawduNvdX87VOSn4QBcRkfQU/JCLiIikR4EuIhIIBbqISCAU6CIigVCgi4gEYkA/JFokXfFtXS/Hi79H9O66D+Llme5+Mi8d60Z8G+pJd/8/+e6LnNsU6FKQ3P0w0eyamNn9QKu7P5yv/pjZoKR5RLqqAlqBtAM9xf5EMqIhFykaZlZhZv8rnhCsLumt0/Vm9oiZbY3nM59hZi/Gc00/GLcZH893/j/jNhviyZJS7fdRi+Z7v9PM/jiej3u7mf3MzMaY2Xjg3wF/Hs9B/kdm9oyZLU7qd2v8vcrMfm5mG4E98aRz1WbWYGZvmtm3BvLnKeFRoEuxMOAxYLG7VwBrgYeSnj/p0bzs/0j0NupvA5OBO+LhG4BJwD+4+x8Ax4D/YGalKfY72N0r3f2/AP8bmB1P4lUD/KW7N8XHfMSjuch/nqKO6cCd7n4l0XzeR919BjAD+KaZTej7j0YkoiEXKRbnEwX0T6PpPSghmlo1YWP8fRewOzE1qZm9QzRd6RHgoLu/Grd7luhDFH6SYr/rkh6XA+viK/jBRPOo99Ub7p7Y7kbg2qSr+ZHAxAz3K6JAl6JhREE9p4fnE3N6nE56nFhO/DvvOs+Fp7Hfj5MePwb8rbtvjF8Ivb+HbU4R//VrZucRhX93+zPgO+5e18N+RPpEQy5SLD4BLjazOQBmVmpmf9jHfVye2B74CtEQSmMf9juS303Z+rWk9ceBEUnLTUQfAQfRvOSlPeyvDvj38bAPZnZlPLmWSEYU6FIsThNNnfoDM9tJ9EEBn+njPhqJPoxjL3Ah0YdanOzDfu8HXjCzbUQfh5bwEnBL4kVR4L8Bn433N4czr8qTPUU0LewvzOwt4En0V7P0g2ZblHNCfDfKjz360GGRIOkKXUQkELpCFxEJhK7QRUQCoUAXEQmEAl1EJBAKdBGRQCjQRUQC8f8Bc2QgqTmZRacAAAAASUVORK5CYII=\n", + "text/plain": [ + "
" + ] + }, + "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)\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": {}, + "source": [ + "The predictiction is not correct. Let's see what went wrong" + ] + }, + { + "cell_type": "code", + "execution_count": 8, + "metadata": {}, + "outputs": [ + { + "data": { + "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXQAAAEKCAYAAAACS67iAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4yLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+j8jraAAAapUlEQVR4nO3df3xV9Z3n8dfHABJ+W9QMGDTMFnEckB/hZ3FsaKtg20HcYUTq0OpDSne3tDqOzEMe2xXr6OMx3bijs661ssraWVcCskyKfTANrUPWrh0VKCC/NoA2lUQ7KDZAbIAQPvvHOTe9xCT3JrmX3Pvl/Xw88uCec773nO8nh7xz8r3nfq+5OyIikv8u6u0OiIhIZijQRUQCoUAXEQmEAl1EJBAKdBGRQCjQRUQCkTLQzWy1mR0xsz0dbL/DzN4ys91m9gszm5D5boqISCrpXKE/D8ztZPuvgM+6+3jgb4BVGeiXiIh0UZ9UDdz9VTMr6WT7L5IWXweKe94tERHpqpSB3kV3A//U0UYzWwosBSgsLCwdNWpUtw5y9uxZLroojOF/1ZKbQqkllDpAtSQcOHDgQ3e/rN2N7p7yCygB9qRoMxvYDwxPZ5+lpaXeXVu2bOn2c3ONaslNodQSSh3uqiUB2OYd5GpGrtDN7DrgWeBmdz+aiX2KiEjX9PjvFzO7EtgALHb3Az3vkoiIdEfKK3QzWwOUAZeaWR2wEugL4O4/AB4EhgPfNzOAM+4+JVsdFhGR9qVzl8uiFNuXAEsy1iMRyRvNzc3U1dVx8uTJrB9r6NCh7N+/P+vHOR/SqaV///4UFxfTt2/ftPeb6btcROQCUldXx+DBgykpKSH+Cz1rTpw4weDBg7N6jPMlVS3uztGjR6mrq2P06NFp7zeMe4BEpFecPHmS4cOHZz3MLzRmxvDhw7v8l48CXUR6RGGeHd35virQRUQCoTF0EclrBQUFjB8/vnW5srKSkpKS3utQL1Kgi0heKywsZOfOne1uS7yDMpQpA1K5MKoUkQtGbW0tY8eO5atf/Srjxo3j8OHDlJeXM3XqVK677jpWrlzZ2vbRRx/l6quv5vrrr2fRokU89thjAJSVlbFt2zYAPvzww9Yr/paWFpYvX966r2eeeQaA6upqysrKWLBgAddccw133HFHYkoUtm7dymc+8xkmTJjAtGnTOHHiBHPnzj3nl9D111/Prl27ely7rtBFJCO++/Je9r13PKP7vHbkEFb+6R932qapqYmJEycCMHr0aB5//HEOHjzID3/4Q2bMmMHmzZs5ePAgb775Ju7OvHnzePXVVxk4cCAVFRXs3LmTM2fOMHnyZEpLSzs91nPPPcfQoUPZunUrp06dYtasWdx0000A7Nixg7179zJy5EhmzZrFa6+9xrRp01i4cCFr165l6tSpHD9+nMLCQhYvXszzzz/PE088wYEDBzh58iQTJvT8oyQU6CKS19oOudTW1nLVVVcxY8YMADZv3szmzZuZNGkSAI2NjRw8eJATJ05w6623MmDAAADmzZuX8libN2/mrbfeYv369QAcO3aMgwcP0q9fP6ZNm0ZxcTR7+MSJE6mtrWXo0KGMGDGCqVOnAjBkyBAAbr31VmbNmkV5eTmrV6/mzjvvzMj3QoEuIhmR6kr6fBo4cGDrY3dnxYoVfOMb3zinzRNPPNHh8/v06cPZs2cBzrkX3N158sknmTNnzjntq6urufjii1uXCwoKOHPmTIf7HzBgADfeeCM/+tGPWLduHdu3b0+vsBQ0hi4iQZszZw6rV6+msbERgPr6eo4cOcINN9xAZWUlTU1NnDhxgpdffrn1OSUlJa0hm7gaT+zr6aefprm5GYADBw7w8ccfd3jssWPH8v7777N161YgeodoIuiXLFnCt7/9baZOncoll1ySkVp1hS4iQbvpppvYv38/M2fOBGDQoEG88MILTJ48mYULFzJhwgQuv/zy1mERgPvvv5/bbruNVatW8aUvfal1/ZIlS6itrWXy5Mm4O5dddhmVlZUdHrtfv36sXbuWb33rWzQ1NVFYWMjPfvYzAEpLSxkyZAh33XVX5ortaKL0bH/pAy4iqiU3hVJLtuvYt29fVvef7Pjx41nd/8qVK728vDyrx0g4fvy419fX+5gxY7ylpaXDdu19f+nkAy405CIicp69+OKLTJ8+nUcffTSj98hryEVEBHjooYfO27G+8pWvfOJF2kzQFbqI9IjHb6CRzOrO91WBLiLd1r9/f44ePapQzzCP50Pv379/l56nIRcR6bbi4mLq6ur44IMPsn6skydPdjngclU6tSQ+sagrFOgi0m19+/bt0ifq9ER1dXXruz3zXbZq0ZCLiEggFOgiIoFQoIuIBEKBLiISCAW6iEggFOgiIoFQoIuIBEKBLiISiJSBbmarzeyIme3pYLuZ2X81s0Nm9paZTc58N0VEJJV0rtCfB+Z2sv1mYEz8tRR4uufdEhGRrkoZ6O7+KvBRJ01uAf4hnnv9dWCYmY3IVAdFRCQ9mZjL5QrgcNJyXbzu/Qzs+xO++/JefrGviadr/iUbuz/vGhpUSy4KpZZQ6oCwahly9hRlZZnf73mdnMvMlhINy1BUVER1dXWX91FXd4qWlhYaGhoy3LveoVpyUyi1hFIHhFVLYWFLt/IvpY4+my75CygB9nSw7RlgUdJyDTAi1T71maIR1ZKbQqkllDrcVUsCWf5M0Y3AV+O7XWYAx9w9K8MtIiLSsZRDLma2BigDLjWzOmAl0BfA3X8AbAK+CBwCfgfcla3OiohIx1IGursvSrHdgW9mrEciItIteqeoiEggFOgiIoFQoIuIBEKBLiISCAW6iEggFOgiIoFQoIuIBEKBLiISCAW6iEggFOgiIoFQoIuIBEKBLiISCAW6iEggFOgiIoFQoIuIBEKBLiISCAW6iEggFOgiIoFQoIuIBEKBLiISCAW6iEggFOgiIoFQoIuIBEKBLiISCAW6iEggFOgiIoFQoIuIBEKBLiISiLQC3czmmlmNmR0yswfa2X6lmW0xsx1m9paZfTHzXRURkc6kDHQzKwCeAm4GrgUWmdm1bZp9B1jn7pOA24HvZ7qjIiLSuXSu0KcBh9z9HXc/DVQAt7Rp48CQ+PFQ4L3MdVFERNJh7t55A7MFwFx3XxIvLwamu/uypDYjgM3AJcBA4Avuvr2dfS0FlgIUFRWVVlRUdKvTjY2NDBo0qFvPzTWqJTeFUksodYBqSZg9e/Z2d5/S7kZ37/QLWAA8m7S8GPhvbdrcB/xV/HgmsA+4qLP9lpaWendt2bKl28/NNaolN4VSSyh1uKuWBGCbd5Cr6Qy51AOjkpaL43XJ7gbWxb8g/gXoD1yaxr5FRCRD0gn0rcAYMxttZv2IXvTc2KbNu8DnAczsj4gC/YNMdlRERDqXMtDd/QywDKgC9hPdzbLXzB42s3lxs78Cvm5mu4A1wJ3xnwYiInKe9EmnkbtvAja1Wfdg0uN9wKzMdk1ERLpC7xQVEQmEAl1EJBAKdBGRQCjQRUQCoUAXEQmEAl1EJBAKdBGRQCjQRUQCoUAXEQmEAl1EJBAKdBGRQCjQRUQCoUAXEQmEAl1EJBAKdBGRQCjQRUQCoUAXEQmEAl1EJBAKdBGRQCjQRUQCoUAXEQmEAl1EJBAKdBGRQCjQRUQCoUAXEQmEAl1EJBAKdBGRQKQV6GY218xqzOyQmT3QQZvbzGyfme01sxcz200REUmlT6oGZlYAPAXcCNQBW81so7vvS2ozBlgBzHL335rZ5dnqsIiItC+dK/RpwCF3f8fdTwMVwC1t2nwdeMrdfwvg7kcy200REUnF3L3zBmYLgLnuviReXgxMd/dlSW0qgQPALKAAeMjdf9LOvpYCSwGKiopKKyoqutXpxsZGBg0a1K3n5hrVkptCqSWUOkC1JMyePXu7u09pb1vKIZc09QHGAGVAMfCqmY1394bkRu6+ClgFMGXKFC8rK+vWwaqrq+nuc3ONaslNodQSSh2gWtKRzpBLPTAqabk4XpesDtjo7s3u/iuiq/UxmemiiIikI51A3wqMMbPRZtYPuB3Y2KZNJdHVOWZ2KXA18E4G+ykiIimkDHR3PwMsA6qA/cA6d99rZg+b2by4WRVw1Mz2AVuA5e5+NFudFhGRT0prDN3dNwGb2qx7MOmxA/fFXyIi0gv0TlERkUAo0EVEAqFAFxEJhAJdRCQQCnQRkUAo0EVEAqFAFxEJhAJdRCQQCnQRkUAo0EVEAqFAFxEJhAJdRCQQCnQRkUAo0EVEAqFAFxEJhAJdRCQQCnQRkUAo0EVEAqFAFxEJhAJdRCQQCnQRkUAo0EVEAqFAFxEJhAJdRCQQCnQRkUAo0EVEAqFAFxEJRFqBbmZzzazGzA6Z2QOdtPszM3Mzm5K5LoqISDpSBrqZFQBPATcD1wKLzOzadtoNBu4B3sh0J0VEJLV0rtCnAYfc/R13Pw1UALe00+5vgO8BJzPYPxERSZO5e+cNzBYAc919Sby8GJju7suS2kwG/qO7/5mZVQP3u/u2dva1FFgKUFRUVFpRUdGtTjc2NjJo0KBuPTfXqJbcFEotodQBqiVh9uzZ29293WHtPj3qFWBmFwF/B9yZqq27rwJWAUyZMsXLysq6dczq6mq6+9xco1pyUyi1hFIHqJZ0pDPkUg+MSloujtclDAbGAdVmVgvMADbqhVERkfMrnUDfCowxs9Fm1g+4HdiY2Ojux9z9UncvcfcS4HVgXntDLiIikj0pA93dzwDLgCpgP7DO3fea2cNmNi/bHRQRkfSkNYbu7puATW3WPdhB27Ked0tERLpK7xQVEQmEAl1EJBAKdBGRQPT4PnSRnqjcUU95VQ3vNTQxclghy+eMZf6kK3q7W5Imnb/cokCXXlO5o54VG3bT1NwCQH1DEys27AZQKOQBnb/coyEX6TXlVTWtYZDQ1NxCeVVNL/VIukLnL/co0KXXvNfQ1KX1klt0/nKPAl16zchhhV1aL7lF5y/3KNCl1yyfM5bCvgXnrCvsW8DyOWN7qUfSFTp/uUcvikqvSbxwprsk8pPOX+5RoEuvmj/pCgVAHtP5yy0achERCYQCXUQkEAp0EZFAKNBFRAKhQBcRCYQCXUQkEAp0EZFAKNBFRAKhQBcRCYQCXUQkEAp0EZFAKNBFRAKhQBcRCYQCXUQkEAp0EZFAKNBFRAKRVqCb2VwzqzGzQ2b2QDvb7zOzfWb2lpm9YmZXZb6rIiLSmZSBbmYFwFPAzcC1wCIzu7ZNsx3AFHe/DlgP/OdMd1RERDqXzhX6NOCQu7/j7qeBCuCW5AbuvsXdfxcvvg4UZ7abIiKSSjqBfgVwOGm5Ll7XkbuBf+pJp0REpOvM3TtvYLYAmOvuS+LlxcB0d1/WTtu/AJYBn3X3U+1sXwosBSgqKiqtqKjoVqcbGxsZNGhQt56ba1RLbgqlllDqANWSMHv27O3uPqXdje7e6RcwE6hKWl4BrGin3ReA/cDlqfbp7pSWlnp3bdmypdvPzTWqJTeFUksodbirlgRgm3eQq+kMuWwFxpjZaDPrB9wObExuYGaTgGeAee5+pDu/dUREpGf6pGrg7mfMbBlQBRQAq919r5k9TPSbYiNQDgwCXjIzgHfdfV4W+y3Srsod9ZRX1fBeQxMjhxWyfM5YgE+smz+ps5eBMnfsbBwnHd+p3M2aNw5z77hm7l6xiUXTR/HI/PG90hc5f1IGOoC7bwI2tVn3YNLjL2S4XyJdVrmjnhUbdtPU3AJAfUMTy1/aBQbNLd66bsWG3QAZDdv2jp2N46TjO5W7eeH1d1uXW9xblxXqYdM7RSUY5VU1rYGa0HzWW8M8oam5hfKqmqwfOxvHSceaNw53ab2EQ4EuwXivoSkrbXuyv0wfJx0tHdy51tF6CYcCXYIxclhhVtr2ZH+ZPk46CqLXsdJeL+FQoEswls8ZS2HfgnPW9b3I6FtwbpAV9i1ofbE0m8fOxnHSsWj6qC6tl3Ck9aKoSD5IvPjYG3e5dHTs3rjLJfHCZ2LMvMBMd7lcIBToEpT5k65oN0TPR7B2dOze8Mj88TwyfzzV1dW8fUdZb3dHzhMNuYiIBEKBLiISCAW6iEggFOgiIoFQoIuIBEKBLiISCAW6iEggFOgiIoFQoIuIBEKBLiISCAW6iEggFOgiIoFQoIuIBEKBLiISCAW6iEggFOgiIoFQoIuIBEKBLiISCAW6iEggFOgiIoFQoIuIBEKBLiISiD7pNDKzucDfAwXAs+7+t222Xwz8A1AKHAUWunttZrsqEq7KHfWUV9XwXkMTI4cVsnzOWF7a9i6vvf1Ra5tZ/+ZT/PmUKz/RDvjEum2//og1bxzm3nHN3L1iE4umj+KR+ePTOu78SVd0uD6d5yeO3eJOgVmXjt1eLeket712F5qUgW5mBcBTwI1AHbDVzDa6+76kZncDv3X3T5vZ7cD3gIXZ6LBIaCp31LNiw26amlsAqG9o4t61Oz/R7rW3Pzon4Osbmli+fhc4NJ/11nX3rd3J2aTntbjzwuvvApwTrO0dd8WG3Wz79Uf87+31n1gPnBOa7T2/J8de/tIuMGhu+X0t6R63vXYXonSGXKYBh9z9HXc/DVQAt7Rpcwvww/jxeuDzZmaZ66ZIuMqralrDqauaW7w1zBPOdtB2zRuHUx63qbmFNW8cbnd9eVVNyuf35NjNZ701zLt63PbaXYjM3TtvYLYAmOvuS+LlxcB0d1+W1GZP3KYuXn47bvNhm30tBZbGi2OB7p6BS4EPU7bKD6olN523Wvr9wadLs7Xvlt8do2DA0Nbl0785tL0nx+3J8zPw3EuBDzt7bvIxclxP/n9d5e6XtbchrTH0THH3VcCqnu7HzLa5+5QMdKnXqZbcFEotZrbtzLEjeV8HhHNOIHu1pDPkUg+MSloujte128bM+gBDiV4cFRGR8ySdQN8KjDGz0WbWD7gd2NimzUbga/HjBcA/e6qxHBERyaiUQy7ufsbMlgFVRLctrnb3vWb2MLDN3TcCzwH/08wOAR8RhX429XjYJoeoltwUSi2h1AGqJaWUL4qKiEh+0DtFRUQCoUAXEQlEzge6mfU3szfNbJeZ7TWz78brR5vZG2Z2yMzWxi/Y5jwzKzCzHWb243g5X+uoNbPdZrbTzLbF6z5lZj81s4Pxv5f0dj/TYWbDzGy9mf0/M9tvZjPzsRYzGxufj8TXcTO7N09r+cv4532Pma2JcyBff1buievYa2b3xuuyck5yPtCBU8Dn3H0CMBGYa2YziKYXeNzdPw38lmj6gXxwD7A/aTlf6wCY7e4Tk+6nfQB4xd3HAK/Ey/ng74GfuPs1wASi85N3tbh7TXw+JhLNq/Q74B/Js1rM7Arg28AUdx9HdDNGYkqRvPpZMbNxwNeJ3nE/AfiymX2abJ0Td8+bL2AA8EtgOtG7rPrE62cCVb3dvzT6XxyfvM8BPwYsH+uI+1oLXNpmXQ0wIn48Aqjp7X6mUcdQ4FfENwjkcy1t+n8T8Fo+1gJcARwGPkV0J96PgTn5+LMC/DnwXNLyfwL+OlvnJB+u0BPDFDuBI8BPgbeBBnc/EzepI/pPkOueIDqZiSkvhpOfdQA4sNnMtsdTOgAUufv78ePfAEW907UuGQ18APyPeCjsWTMbSH7Wkux2YE38OK9qcfd64DHgXeB94Biwnfz8WdkD/ImZDTezAcAXid6EmZVzkheB7u4tHv0ZWUz0p8s1vdylLjOzLwNH3D1f5ppI5Xp3nwzcDHzTzG5I3ujRpUc+3BPbB5gMPO3uk4CPafPnbx7VAkA8tjwPeKnttnyoJR5PvoXol+1IYCAwt1c71U3uvp9oqGgz8BNgJ9DSpk3GzkleBHqCuzcAW4j+3BoWTzMA7U9HkGtmAfPMrJZoxsrPEY3d5lsdQOtVFO5+hGicdhrwr2Y2AiD+90jv9TBtdUCdu78RL68nCvh8rCXhZuCX7v6v8XK+1fIF4Ffu/oG7NwMbiH5+8vVn5Tl3L3X3G4jG/g+QpXOS84FuZpeZ2bD4cSHRvOz7iYJ9Qdzsa8CPeqeH6XH3Fe5e7O4lRH8O/7O730Ge1QFgZgPNbHDiMdF47R7OnQIiL2px998Ah81sbLzq88A+8rCWJIv4/XAL5F8t7wIzzGxAPA134pzk3c8KgJldHv97JfBvgRfJ0jnJ+XeKmtl1RHOtFxD9Alrn7g+b2R8SXel+CtgB/IW7n+q9nqbPzMqA+939y/lYR9znf4wX+wAvuvujZjYcWAdcCfwauM3dP+pgNznDzCYCzwL9gHeAu4j/r5F/tQwkCsQ/dPdj8bq8Oy/x7ckLgTNEPxdLiMbM8+pnBcDMfk70elkzcJ+7v5Ktc5LzgS4iIunJ+SEXERFJjwJdRCQQCnQRkUAo0EVEAqFAFxEJxHn9kGiRdMW3db0SL/4B0bvrPoiXp7n76V7pWDvi21BPu/sversvcmFToEtOcvejRLNrYmYPAY3u/lhv9cfM+iTNI9JWGdAIpB3oKfYn0i0acpG8YWalZvZ/4gnBqpLeOl1tZo+b2bZ4PvOpZrYhnmv6kbhNSTzf+f+K26yPJ0tKtd8nLJrv/R4z+9N4Pu4dZvYzMysysxLg3wF/Gc9B/idm9ryZLUjqd2P8b5mZ/dzMNgL74knnys1sq5m9ZWbfOJ/fTwmPAl3yhQFPAgvcvRRYDTyatP20R/Oy/4DobdTfBMYBd8bDNwBjge+7+x8Bx4H/YGZ9U+y3n7tPcff/AvxfYEY8iVcF8NfuXhsf83GP5iL/eYo6JgP3uPvVRPN5H3P3qcBU4OtmNrrr3xqRiIZcJF9cTBTQP42m96CAaGrVhI3xv7uBvYmpSc3sHaLpShuAw+7+WtzuBaIPUfhJiv2uTXpcDKyNr+D7Ec2j3lVvunvieTcB1yVdzQ8FxnRzvyIKdMkbRhTUMzvYnpjT42zS48Ry4v9523kuPI39fpz0+Eng79x9Y/xC6EMdPOcM8V+/ZnYRUfi3tz8DvuXuVR3sR6RLNOQi+eIUcJmZzQQws75m9sdd3MeViecDXyEaQqnpwn6H8vspW7+WtP4EMDhpuZboI+Agmpe8bwf7qwL+fTzsg5ldHU+uJdItCnTJF2eJpk79npntIvqggM90cR81RB/GsR+4hOhDLU53Yb8PAS+Z2Xaij0NLeBm4NfGiKPDfgc/G+5vJuVflyZ4lmhb2l2a2B3gG/dUsPaDZFuWCEN+N8mOPPnRYJEi6QhcRCYSu0EVEAqErdBGRQCjQRUQCoUAXEQmEAl1EJBAKdBGRQPx/tRIgu/jEeqIAAAAASUVORK5CYII=\n", + "text/plain": [ + "
" + ] + }, + "metadata": { + "needs_background": "light" + }, + "output_type": "display_data" + } + ], + "source": [ + "data_pred.plot(x=\"Temperature\",y=\"Frequency\",kind=\"line\",ylim=[0,1.2])\n", + "plt.scatter(x=data[\"Temperature\"],y=data[\"Frequency\"])\n", + "plt.grid(True)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "There were warnings during the construction of the log model. let's try and change it" + ] + }, + { + "cell_type": "code", + "execution_count": 12, + "metadata": {}, + "outputs": [ + { + "ename": "DistributionNotFound", + "evalue": "The 'statsmodel==0.9.0' distribution was not found and is required by the application", + "output_type": "error", + "traceback": [ + "\u001b[0;31m---------------------------------------------------------------------------\u001b[0m", + "\u001b[0;31mDistributionNotFound\u001b[0m Traceback (most recent call last)", + "\u001b[0;32m\u001b[0m in \u001b[0;36m\u001b[0;34m\u001b[0m\n\u001b[1;32m 1\u001b[0m \u001b[0;32mimport\u001b[0m \u001b[0mpkg_resources\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m----> 2\u001b[0;31m \u001b[0mpkg_resources\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mrequire\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m\"statsmodel==0.9.0\"\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m 3\u001b[0m \u001b[0;32mimport\u001b[0m \u001b[0mstatsmodel\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n", + "\u001b[0;32m~/miniconda3/lib/python3.7/site-packages/pkg_resources/__init__.py\u001b[0m in \u001b[0;36mrequire\u001b[0;34m(self, *requirements)\u001b[0m\n\u001b[1;32m 899\u001b[0m \u001b[0mincluded\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0meven\u001b[0m \u001b[0;32mif\u001b[0m \u001b[0mthey\u001b[0m \u001b[0mwere\u001b[0m \u001b[0malready\u001b[0m \u001b[0mactivated\u001b[0m \u001b[0;32min\u001b[0m \u001b[0mthis\u001b[0m \u001b[0mworking\u001b[0m \u001b[0mset\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 900\u001b[0m \"\"\"\n\u001b[0;32m--> 901\u001b[0;31m \u001b[0mneeded\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mresolve\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mparse_requirements\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mrequirements\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m 902\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 903\u001b[0m \u001b[0;32mfor\u001b[0m \u001b[0mdist\u001b[0m \u001b[0;32min\u001b[0m \u001b[0mneeded\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n", + "\u001b[0;32m~/miniconda3/lib/python3.7/site-packages/pkg_resources/__init__.py\u001b[0m in \u001b[0;36mresolve\u001b[0;34m(self, requirements, env, installer, replace_conflicting, extras)\u001b[0m\n\u001b[1;32m 785\u001b[0m \u001b[0;32mif\u001b[0m \u001b[0mdist\u001b[0m \u001b[0;32mis\u001b[0m \u001b[0;32mNone\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 786\u001b[0m \u001b[0mrequirers\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mrequired_by\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mget\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mreq\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0;32mNone\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m--> 787\u001b[0;31m \u001b[0;32mraise\u001b[0m \u001b[0mDistributionNotFound\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mreq\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mrequirers\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m 788\u001b[0m \u001b[0mto_activate\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mappend\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mdist\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 789\u001b[0m \u001b[0;32mif\u001b[0m \u001b[0mdist\u001b[0m \u001b[0;32mnot\u001b[0m \u001b[0;32min\u001b[0m \u001b[0mreq\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n", + "\u001b[0;31mDistributionNotFound\u001b[0m: The 'statsmodel==0.9.0' distribution was not found and is required by the application" + ] + } + ], + "source": [ + "import pkg_resources\n", + "pkg_resources.require(\"statsmodel==0.9.0\")\n", + "import statsmodel" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "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']], \n", + " family=sm.families.Binomial(sm.families.links.logit)).fit()\n", + "\n", + "logmodel.summary()\n" + ] + }, + { + "cell_type": "markdown", + "metadata": { + "hideCode": false, + "hidePrompt": false, + "scrolled": true + }, + "source": [ + "This figure is very similar to the Figure 4 of Dalal *et al.* **I have managed to replicate the Figure 4 of the Dalal *et al.* article.**" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Computing and plotting uncertainty" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Following the documentation of [Seaborn](https://seaborn.pydata.org/generated/seaborn.regplot.html), I use regplot." + ] + }, + { + "cell_type": "code", + "execution_count": 9, + "metadata": {}, + "outputs": [ + { + "data": { + "image/png": "iVBORw0KGgoAAAANSUhEUgAAAZAAAAEQCAYAAACeDyIUAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4yLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+j8jraAAAgAElEQVR4nO3deXwUZZ4G8Keq+sp9n9ygQJRLATOjKGNAkoGAqIs4uLrrgbOKos44HxEdOVQcdFdFvEadGceF1V3GC6ICooMCIociAcKhIZAAnatzd/qsevePJm3CIZ0i6aQ7z/fzQehOVfXvtTt5Um+99b6SEEKAiIioneSuLoCIiEITA4SIiHRhgBARkS4MECIi0oUBQkREujBAiIhIl6AEyNKlS5GTk4MhQ4bg0KFDZ9xGVVUsWrQIEydOxDXXXINVq1YFozQiItIpKAEyYcIErFy5Er169TrrNmvWrEFpaSnWr1+P//3f/8Xy5ctx7NixYJRHREQ6BCVAxowZg4yMjJ/d5pNPPsGMGTMgyzISExMxceJErF27NhjlERGRDt3mGojVakVmZqb/cUZGBsrLy7uwIiIi+jndJkCIiCi0GLq6gBYZGRk4ceIERowYAeD0M5JA1dbaoWnhOb1XUlI0bLamri6j04Rz+8K5bQDbF8pkWUJCQpSufbtNgOTl5WHVqlWYNGkS6urqsGHDBqxcubLdx9E0EbYBAiCs2waEd/vCuW0A29cTBaUL68knn8RVV12F8vJy3HbbbZgyZQoAYPbs2dizZw8A4Nprr0Xv3r0xadIk3HjjjZgzZw769OkTjPKIiEgHKdymc7fZmsL2N4WUlBhUVTV2dRmdJpzbF85tA9i+UCbLEpKSovXt28G1EBFRD8EAISIiXRggRESkCwOEiIh0YYAQEZEuDBAiItKFAUJERLowQIiISBcGCBER6cIAISIiXRggRESkCwOEiIh0YYAQEZEuDBAiItKFAUJERLowQIiISBcGCBER6cIAISIiXRggRESkCwOEiIh0YYAQEZEuDBAiItKFAUJERLowQIiISBcGCBER6cIAISIiXRggRESkCwOEiIh0YYAQEZEuDBAiItKFAUJERLowQIiISBcGCBER6cIAISIiXRggRESkiyFYL1RSUoJ58+ahrq4O8fHxWLp0Kfr3799mG5vNhkceeQRWqxVerxfZ2dl47LHHYDAErUwiIgpQ0M5AFixYgFmzZmHdunWYNWsWHn/88dO2ee211zBo0CCsWbMGq1evxr59+7B+/fpglUhERO0QlACx2WwoKipCfn4+ACA/Px9FRUWoqalps50kSbDb7dA0DW63Gx6PB2lpacEokYiI2ikofUNWqxVpaWlQFAUAoCgKUlNTYbVakZiY6N/unnvuwX333Ydx48bB4XDg5ptvxujRo9v1WklJ0R1ae3eTkhLT1SV0qnBuXzi3DWD7eqJudXFh7dq1GDJkCP7+97/Dbrdj9uzZWLt2LfLy8gI+hs3WBE0TnVhl10lJiUFVVWNXl9Fpwrl94dw2gO0LZbIs6f7FOyhdWBkZGaioqICqqgAAVVVRWVmJjIyMNtutWLEC06ZNgyzLiImJQU5ODrZt2xaMEomIqJ2CEiBJSUnIyspCQUEBAKCgoABZWVltuq8AoHfv3vjqq68AAG63G1u3bsWFF14YjBKJiKidgjYKa+HChVixYgVyc3OxYsUKLFq0CAAwe/Zs7NmzBwAwf/58fPvtt5g6dSqmT5+O/v3748YbbwxWiURE1A6SECKsLhjwGkjoCuf2hXPbALYvlHX7ayBERBR+GCBERKQLA4SIiHRhgBARkS4MECIi0oUBQkREujBAiIhIFwYIERHpwgAhIiJdGCBERKQLA4SIiHRhgBARkS4MECIi0oUBQkREujBAiIhIFwYIERHpwgAhIiJdGCBERKQLA4SIiHRhgBARkS4MECIi0oUBQkREujBAiIhIFwYIERHpwgAhIiJdGCBERKQLA4SIiHRhgBARkS4MECIi0oUBQkREujBAiIhIFwYIERHpEnCAbNiwAV6vtzNrISKiEBJwgLz44osYN24cFi9ejN27d3dmTUREFAICDpDVq1fjrbfegtlsxn333Yfc3Fy88sorOHbsWED7l5SUYObMmcjNzcXMmTNx5MiRM273ySefYOrUqcjPz8fUqVNRXV0daIlERBREkhBCtHcnIQS2bt2KP/3pT/jhhx9w6aWXYubMmcjPz4csnzmTbr31Vtxwww249tpr8dFHH+G9997D22+/3WabPXv24OGHH8bf//53pKSkoLGxESaTCWazOeDabLYmaFq7mxQSUlJiUFXV2NVldJpwbl84tw1g+0KZLEtISorWt297dygtLcXLL7+MhQsXwuVyYe7cuZgxYwZWrlyJuXPnnnEfm82GoqIi5OfnAwDy8/NRVFSEmpqaNtu99dZbuP3225GSkgIAiImJaVd4EBFR8BgC3XDlypX46KOPcPToUfz617/GM888g1GjRvm/npubi8svv/yM+1qtVqSlpUFRFACAoihITU2F1WpFYmKif7vi4mL07t0bN998M5qbm3HNNdfg7rvvhiRJettHRESdJOAA+eqrr3DbbbdhwoQJMJlMp309IiICy5cvP69iVFXFwYMH8be//Q1utxt33nknMjMzMX369ICPofdULFSkpMR0dQmdKpzbF85tA9i+nijgAHnxxRchyzKMRqP/OY/HAyGEP1DGjRt3xn0zMjJQUVEBVVWhKApUVUVlZSUyMjLabJeZmYm8vDyYTCaYTCZMmDABhYWF7QoQXgMJXeHcvnBuG8D2hbKgXAO5/fbbsW/fvjbP7du3D3fcccc5901KSkJWVhYKCgoAAAUFBcjKymrTfQX4ro1s3rwZQgh4PB588803GDp0aKAlEhFREAUcIAcPHsTIkSPbPDdixAgcOHAgoP0XLlyIFStWIDc3FytWrMCiRYsAALNnz8aePXsAAFOmTEFSUhImT56M6dOn44ILLsC//Mu/BFoiEREFUcBdWLGxsaiurvaPkAKA6upqREREBLT/oEGDsGrVqtOef+ONN/z/lmUZjzzyCB555JFAyyIioi4S8BnIpEmT8Pvf/x6HDh2Cw+HAwYMH8fDDD+PXv/51Z9ZHRETdVMAB8uCDD2LQoEGYMWOG/8bBAQMG4He/+11n1kdERN1Uu+9EF0KgtrYWCQkJ3fL+DI7CCl3h3L5wbhvA9oWy8xmFFfA1EABobGxESUkJ7HZ7m+d/+ctf6npxIiIKXQEHyPvvv4/FixcjMjISFovF/7wkSfj88887pTgiIuq+Ag6Q559/HsuWLcP48eM7sx4iIgoRAV9EV1X1rHeaExFRzxNwgMyePRuvvvoqNE3rzHqIiChEBNyF9dZbb6G6uhpvvvkm4uPj23xt48aNHV0XERF1cwEHyLPPPtuZdRARUYgJOEAuu+yyzqyDiIhCTMDXQNxuN55//nlMmDABo0ePBgBs3rwZK1as6LTiiIio+wo4QJYsWYJDhw7hP//zP/13oF944YV45513Oq04IiLqvgLuwtqwYQPWr1+PyMhIyLIvd9LS0lBRUdFpxRERUfcV8BmI0WiEqqptnqupqTltRBYREfUMAQdIXl4eHn74YZSVlQEAKisrsXjxYkyZMqXTiiMiou6rXdO59+7dG9OmTUNDQwNyc3ORmpqKOXPmdGZ9RETUTbV7OnfA13XF6dyDL5ynlAbCu33h3DaA7QtlQZnOvaXrqkXrKd379Omj68WJiCh0BRwg11xzDSRJQusTlpYzkP3793d8ZURE1K0FHCAHDhxo87iqqgovvfQSxowZ0+FFERFR9xfwRfRTpaSk4NFHH8Vzzz3XkfUQEVGI0B0gAHD48GE4HI6OqoWIiEJIwF1Ys2bNajPqyuFw4Mcff+QwXiKiHirgAJkxY0abxxERERg6dCj69+/f0TUREVEICDhArrvuus6sg4iIQkzAAbJs2bKAtrv//vt1F0NERKEj4AA5evQo1q9fj2HDhqFXr144ceIE9uzZg0mTJsFsNndmjURE1A0FHCBCCPzXf/0XcnNz/c+tX78ea9euxdNPP90pxRERUfcV8DDer776ChMnTmzzXE5ODr788ssOL4qIiLq/gAOkX79+WLlyZZvn3nnnHfTt27fDiyIiou4v4C6sJ598Evfeey/efPNN/0qEBoMBy5cv78z6iIiomwo4QC666CKsW7cOu3fvRmVlJVJSUjBq1CgYjcbOrI+IiLop3VOZjB07Fh6PB83NzR1ZDxERhYiAz0AOHjyIu+++GyaTCRUVFZg8eTJ27NiBDz74AC+88EJn1khERN1QwGcgCxcuxNy5c7F27VoYDL7cGTt2LL799tuA9i8pKcHMmTORm5uLmTNn4siRI2fd9vDhwxg5ciSWLl0aaHlERBRkAQfIjz/+iGuvvRbATwtJRUZGwuVyBbT/ggULMGvWLKxbtw6zZs3C448/fsbtVFXFggULThsyTHQuQgBa+1doJiKdAg6QXr16Ye/evW2eKywsDGgYr81mQ1FREfLz8wEA+fn5KCoqQk1NzWnbvv766/jVr37FSRpJF5dH7eoSiHqMgK+B3H///fjtb3+Lm266CR6PB3/+85/x7rvv4oknnjjnvlarFWlpaVAUBQCgKApSU1NhtVqRmJjo3+7AgQPYvHkz3n77bbzyyis6mgPdi8OHipSUmK4uoVOdT/tUTcBW50ByQkSbpQe6C753oS3c26dHwAFy9dVX480338T//d//YezYsTh+/DiWL1+OYcOGdUghHo8Hf/zjH/H000/7g0YPm60Jmhae3RgpKTGoqmrs6jI6zfm2TwjA1uCE1+2BUTmvtdI6HN+70BbO7ZNlSfcv3gEFiKqqyM3NxSeffIKFCxe2+0UyMjJQUVEBVVWhKApUVUVlZSUyMjL821RVVaG0tBR33XUXAKChoQFCCDQ1NQV0lkME+K6BOFxeGCNNXV0KUdgLKEAURYGiKHC5XDCZ2v+NmZSUhKysLBQUFODaa69FQUEBsrKy2nRfZWZmYtu2bf7Hy5cvR3NzMx5++OF2vx71bE6XiiiLBkXuXmchROEm4O+wW2+9FQ888AC2b9+O0tJSlJWV+f8EYuHChVixYgVyc3OxYsUKLFq0CAAwe/Zs7NmzR1/1RGegCQG709vVZRCFPUmInx/3WFVVhZSUFAwdOtS3gySh9S6SJGH//v2dW2U78BpI6OqIayBV9Q5omoAsSUiMs8Agd4+L6XzvQls4t+98roGc8wykZf2PAwcO4MCBA8jJyfH/+8CBA90qPIhaaEKg0e5GeP4qQdQ9nDNATj1B2bFjR6cVQ9SRXB4VDhe7sog6yzkD5NTx9Ofo8SLqVpqaPfCGaZcmUVc75ygsVVXxzTff+IPj1McA8Mtf/rLzKiQ6Dy1dWfExZnSPqyFE4eOcAZKUlIT58+f7H8fHx7d5LEkSPv/8886pjqgDuDwqml1eRJkDvm+WiAJwzu+oL774Ihh1EHUqe7MHZqPSbUZlEYUD/kpGIa+wuBprt5XC6VYRZTFgzNBUDOmb0GYbTQg02N1IjDEB7MwKWS3vdXW9E8lxFuRl98WIQcldXVaPxVt1KaQVFldj5WeHUGd3I8KioMHhweotJThYWnvath6PCrdH64IqqSO0fq8jLQbU2d1Y+dkhFBZXd3VpPRYDhELa2m2lUBQZZqMCSZJgMihQFBmbdp84bVsBoNHhCX6R1CFOfa/NRt97vXZbaVeX1mMxQCikVdc7YTK0/RgbFRm1jWde6Mzj1bhmSIg603ttMsiornd2UUXEAKGQlhxngdvbtlvKo2pIiDGfdR+704tuuFwIncOZ3mu3V0NynKWLKiIGCIW0vOy+UFXfWYUQAm6vClXVcOXIzLPu4/GqcPFaSMg59b12eXzvdV72uVdFpc7BUVgU0lpG4KzdVgqHU0VshBE5l/Q6bRRWa0IATrcXJgPXDAklrd9rjsLqHhggFPJGDErGiEHJbWbjPReHy4tIsxEGhX1ZoaTlvabugV1Y1CMJATQ2u7u6DKKQxgChHsvlUdHs5my9RHoxQKhHa2r2QNV4QZ1IDwYI9WiaJtBg93DhKSIdGCDU43HhKSJ9GCBEAJoc7Moiai8GCBF8XVmN7MoiahcGCNFJzpMLTxFRYBggRK00NbvhUdmVRRQIBgiFBVu9Ext2lqHpPKdrFwJosLsh2JlFdE6cyoTCwuotJdhUaEVCjBn/ljcUqQkRuo/l8WpocngRG2mEYI4QnRXPQCgsXDE8Ayajbx2Q1z7aix+P1Z/X8ZqdHjjdXDeE6OcwQCgsDO4Tj4dnXYqYSCOcbhVvfbof24oqdB9PCKDR7oYawMSMRD0VA4TCRv/0WMy5fjgykiKhCeCjzSVY8/UR3SHg1YTveggzhOiMGCAUVuKjzbhr2sXI6udbD2Tr3nL8/dMDuu80d3lUNDrcXMGQ6AwYIBR2zEYFN18zGFedXJXwx+P1eOWDvaioadZ1PIfTy/tDiM6AAUJhSZYl5GX3xY1XXwCDIsHW4MSrH+3F3pKadh9LAGhs9sDL6yFEbTBAKKyNujAZd027GHFRJrg9Gv7ns0NYt700oFULW9M0gfomFzReECHyY4BQ2OudEo051w/HgIxYAMCX35/A3z7d3+6bDj1eDQ3Nbt5iSHQSA4R6hOgII26fkoVxwzMAAMXHG/DS+3twtLyxXcdxutSTd6oTUdACpKSkBDNnzkRubi5mzpyJI0eOnLbNyy+/jClTpmDq1Km4/vrrsWnTpmCVRz2AIkuY/Mt++M3EC2Eyymiwu/HGmiJ8tftEu7qmHC4vh/cSIYgBsmDBAsyaNQvr1q3DrFmz8Pjjj5+2zYgRI/CPf/wDa9aswZIlS/Dggw/C6XQGq0TqIYYPTMKc64YjLSECmhBYu60U/732YLu6tBwuL+qbXZ1YJVH3F5QAsdlsKCoqQn5+PgAgPz8fRUVFqKlpOyLmyiuvRESEbw6jIUOGQAiBurq6YJRIPUxKfATuvm4YxgxJAQAcLKvDS+8VovhE4FOgOF0q6u3uziqRqNsLSoBYrVakpaVBURQAgKIoSE1NhdVqPes+H374Ifr27Yv09PRglEg9kMmg4Prxg3BjzgW+Lq1mD/5asB/rtpfCG+CU7g6XF40OD8CrItQDdcvZeLdv345ly5bhr3/9a7v3TUqK7oSKuo+UlJiuLqFTnU/7VE3AK0vtvjaRc1kUhl2Ygr+s3oej1gZ8+f0JlFgbcfu0i5GeFBXQMQxmA+JjzJB+5pZ1vnehLdzbp0dQAiQjIwMVFRVQVRWKokBVVVRWViIjI+O0bXft2oU//OEPeOWVVzBw4MB2v5bN1tTuMf6hIiUlBlVV7Rs1FErOt31CALX1Dl3vvwHAnVOGYsPOY/jq+xMorWjEk3/dhrzsvvjFxemQzzGXSQ0AW40BsVFGSDh9W753oS2c2yfLku5fvIPShZWUlISsrCwUFBQAAAoKCpCVlYXExMQ22xUWFuLBBx/Eiy++iIsvvjgYpRH5KbKM3Mv64s6pFyE+2gSvKlDw9VH89eP9qG0892AOh8uL2kY33F4N7NKinkASIjiDEYuLizFv3jw0NDQgNjYWS5cuxcCBAzF79mzMnTsXw4cPxw033IDjx48jLS3Nv98zzzyDIUOGBPw6PAMJXR1xBlKl8wzkVE63Fx9/fRTfHqoCAJiMMvKy++KyrLRzno1IEmA0KIi2GGE2yRCC712oC+f2nc8ZSNACJFgYIKGrOwVIiw07y/Dl9yf8U8KnJUTg8mHp2P1jNWobXUiIMePKkZkY0jfhtH0lCTAbFERGGNErI0532wqLq7F2Wymq651IjrMgL7svRgxKPq92dZTVmw9j/Y5jcHpUWIwKJo3tjWnj2t/13N2F8/det+/CIgpFB0trseuHKsRFmxBh9o0grKh14INNJSivdcBsUtDg8GD1lhIcLK09bX8hAKdHRW2jE7Z6p655tAqLq7Hys0Oos7sRaTGgzu7Gys8OobC4+rzbd75Wbz6M1V8fgcujwiD7pr5f/fURrN58uKtLoyBhgBCdxabdJ6AoMiwmAxJiLEiMNfu/1uz0orreCQhAUWRs2n3irMcRwtclVlPvhMPtbdfVkbXbSqEoMsxGBZIkwWxUoCgy1m4rPY+WdYz1O45BggRFliBJsu9vSFi/41hXl0ZBwgAhOovaRheMyk/fIhaToc34Kq8qUF3vhN3hga3h3BfZvZpAfZMbNQ0OOD2BrbdeXe+EydD229RkkH3h1cWcbi/kUy4HyZLveeoZGCBEZ5EQY4bnlBsKDYoEoyIhOc7iDxeHS0VjswffFJUHdP3F4xWoa3ShptEFt0f92TOS5DjLyVFdP3F7NSTHWdrdno5mMRlwanM14XueegYGCNFZXDkyE6qqwe1VIYSA26vCZFRgMhkACUiKMyMqwndWoglg9eYjePmDPTh8oiGg47s9KmoaXbA1OGB3eqBq2mlL5+Zl94WqanB5fDW4PCpUVUNedt+Ob3A7TRrbGwICqiYghOb7GwKTxvbu6tIoSJSFCxcu7OoiOpLDEb6zpEZFmdHcHL5zL3VE+5pd3g57/5PjIpAcZ0GFrRlNzR7ER5nw61/0w0X9E3zPOTxIjrPgmrF9YDEZYD353HeHqlBua0av5ChEWny/jUdEmOA4y2SNmga4PRqcLhVeVUCWJfhObiSkJUYiLSECxyqbUG93IzHGjOuvGtgtRmEN6ZsACIGj5U1wqwIWo4LJv+gblqOwwvl7T5IkREaa9O3LYbyhI5yHEgLdcxhvexyrbELB1iMorWgC4Js+PvuiNFx9aS/0yYxHTY09oOO03EcSaTHAbJTPeGd7d8PPZug6n2G87Kwk6iC9U6Px22kXo7DYhnXbS1HX5MbXe8vx7cEq5P6iHy65IAlmo3LO4wjh695ye1TIsgSzQYHZrMBkkM95EyNRMDFAiDqQJEkYeUEyLuqfiK37yrFx13E43SpWbzqMz3eW4VejMnFZVhqMhsAuP2qagMPthcPthSxLiDAbEGFWYFTksO2qpdDBACHqBEaDjKtGZmLMkFR8tfs4tu6rgN3hwcdbj2JToRW/GpWJMUNTYVACH8eiaQJ2hwfNTg8MBt+9ISaDctowX6Jg4TWQEBLO/bBA6F8D+TmSUcEHX/yAbw9W+adFiY0y4aqRGRgzNBUmw7m7ts54XAlQFAkRJgPMRqVdgdSR+NkMXbwGQtTNJcRYMP3KgRg/KhP/3HUC3x2sQoPdjYKvj+Kfu07gimHpyL4oDRHm9n1LCgF4vQKNXg+aJA8UWYLZZIDZoMBklCFJYFcXdRoGCFEQJcRYcP1VA3H1JZn48vsT+PZgFewOD9bv8E3aOHZoKi4fno74aPO5D3YKIXx3x3sdHtjhgSxLMBlkmEwKTIoMRZFDYDwXhRIGCFEXaDkjybm0N7bssWLb/gq4PCo277Hi671WDBuYhCuGp6NPqv5V8DRNwOlW4XSrkCRAliQYDTIMBhkGRYZB8c1jJUvtX8WRCGCAEHWp2JM3J/7qkl7YVlSBrXvL0ejwoLDYhsJiG3qnROGXw9IxfGDSeV3fEAJQhYDqVgG3bx4uSfKNGpNlwCDLvrvsDb4zFZldXxQAXkQPIeF8IQ8I74voiYlRAd1I6FU1FBbbsGWPFVZbs//5SIsBY4akYGxWGpJiO28eLEkCFEmCQZGhGE6epUgyJNl3Y6RyljXn+dkMXbyIThQmDIqMSwen4JILk3GkvBFb95WjqKQGzU4vvtptxVe7rRjUKxZjhqTiov6JAd9PEighAK8Q8Goq0GrG4JazFYMiw2yUYVQUKCe7wKjnYoAQdUOSJGFARiwGZMSiwe7GjgOV2HmgEvV2N4qPN6D4eAMsJgUjL0jGpYOT0TslGlIn3qUuBHwTSmq+O+QBD+ST3V8mgwJLpG8teEX2rS3P0V89AwOEqJuLjTJhwujeuPqSXvjhWB12HKjEgaN1cLpVbCuqwLaiCiTHWTDygmSMvCAJyXERQalLEwKaCnhVLxqa3ahpcP50sf5kF1jLmYsiSzDIEmRZ8l13YcCEBQYIUYiQZQlD+iZgSN8ENDa7sftHm2/m35pmVNc78fm3x/D5t8fQKzkKwwclYfjARCTEBHfdEP/F+lO6wICfusEAX1ed0eD7I58MlJagaX0i1TpkWp7XBBhA3QQDhCgExUSaMG5EBsaNyIDVZsf3P1Rjd7ENDXY3jlfbcbzajrXbStErOQoXD0jERf0TkZoQnDOTs2npBgPQqivMR5Lgm3VYQpswkX1fAIRvRUdN0yAEIJ0cOaYosj90WgLFNxBAhiz7QleWJHapdRIGCFGIy0iKQkZSFHKz++JoeSMKi23YV1KDJofHHybrd5QhOc6CrH4JGNovAX3TYrrVBXAhAAEBCMC3/qLAmVdPOUkDvFABnH1p4JYzHgm+desNsi+M/OEkS5DR8hz8C4MJCEjCdwBFAiRZgqaJ0xb78gWSAHrw7ZkcxhtCwnkoIcBhvB1J0wSOlDdib4kN+4/Uot7edjGkCLOCC3vHY3CfeFzQOw6xOhcUahHs9nUkyf+fU7rMTv5HgoSExEjU1zWffOyLDQgJibHmkJ9in8N4iagNWZYwMDMWAzNjMfXy/jhebcf+I7U4UFoLq60ZDpfqv1kRANITI3FB7zgMyoxF/4zYgNYtCRf+E4mzPC8g/NPEtNaNTuC6DAOEKMxJkoTeKdHonRKNa8b2QV2TCwdL6/DDsTr8eLwebo+G8ppmlNc0Y3OhFbIkoXdqFAZmxmFARgz6psbAbOo5gUKBY4AQ9TDx0WZkX5SG7IvS4FU1lFU24Ydj9Sg+Xo9jVU3QhEBpRRNKK5qwcZfvN+2M5Cj0S4tBv/QY9E2LQVzU+XV5UXhggBD1YAZF9t+wiLF94HR7UWJtxOET9SixNsJqs0MTwPEqO45X2fH13nIAvntT+qREo09qNHqlRiEiqv2zB1PoY4BQ+JCAaIvhpyE1OHlRVIgzdXG33u3kP6SfHkttv3byMBBCnHE4qKYKeDQVQvNtp4Xo2BSLyYCsfgnI6pcAAHC6vTha3uj7U9GIY5V2eFQNDXY39ro0cfYAABGfSURBVNlrsO9Izck99yM5zoLM5ChkJkUhIzkSGUlRiI4wdl1jqNMxQChsSAAiLV33A6vlXgNVE/CoGjRNQBMCqle0zrSzjvrsjpljMRn8Ny8CgKppKLc1o7SyCccqm1BW2QRbvRMCQHW9E9X1Tv+FeQCIiTAiPSkSaYmRSEuIQFpCJFISInrURfpwxgAh6iAtAeCbtbbtD8jEhEiIkzfOydJPd1W37CMACA1QhQa3W4PLq571bKcrKbKMXinR6JUSDVzse87p9qLRpeFASTWs1c04Xm1Hdb0DQgCNDg8aj9Xjh2P1bY4TH21CSnwEUuIjkBxvQXJcBJLjLIiNMoX8sNiehAFCFASKIsN0rplzFQCQEWFq6QbT4FEFHE4v3F6124VJC4vJgMz0KKTE/HRh3ePVUFHbjHKbb3RXRW0zymscsDt8twfWNblR1+Q+LViMiozEWDMSYy1IirUgIdaMxBgzEmItSIg2d/jsw3R+GCBE3VDLdByKDESYFHhVDU63CofTC7UbnpmcymiQ/UOHW7M7PaisdaCy1oHqOgeq6h2oqnOirtEFAcCjaqiodaCi1nHG40ZHGBEfbUJ8jBnxUWbERZsQF21GXJQJcVEmREcYIfMGjaBhgBB1c0L4uo6iLDIiLUZoqua7zqIJeNwq3Cevt4SCKIsRAzKMvlFfrXi8GmoanbDVO2Fr8P1d2+hCTYMLdU0uqCfb1+TwoMnhwbGqM9/1Lkm+6y4xUSbERJgQE2k8+ccXLjGRRkRHGBEVYYTJIHfqFPg9AQOEKIS0zOukKIAJgGQxQNUEvF4Bj+YLEqEJaJqAV9N8czuFyBlLWkIk0hIiT/uapgk0NLtR2+hCbaMvUOqa3KhvcqHe7kZ9kxuuk9eXhAAamj1oaPYA+PmpVYyKjEiLAVERRkRZDIiyGBFhMSDSbEBkq78jzAaosgy3ywuzSeE1mlYYIEQhTAjf+hsmowQTfro+0PIzTtWEbwoOIXwrDXo0eFQVqhYawQL4pmWJjzYjPtqMARln3sbp9qLB7kGD3Y2GZjcam91osHvQ6HCjsdmDpmbfv90ezb+PR9V8AXTKPGE/RwJgNimwmBREWYy4bvxAXHphynm2MHQFLUBKSkowb9481NXVIT4+HkuXLkX//v3bbKOqKp588kls2rQJkiThrrvuwowZM4JVIlHYaAkGWZJgMkgoLK7Guu2lqG9yIy0xEhNG98bgPvH4x8YfsedwDTyqBlUVGNonDhcNSMLmwhOw1TsRG2XClSMzMaRvAg6W1mLT7hOobXQhIcaMK0dm4nhVEzYXlsPlVWE2KBg3Ih05o/ucsaYz7X+247YMGw70GF/vsfrq8KgwG311TL2i/2n7u70q9vxow9Z95ai3uxFpNqBvWjScHhVHyxvhcquQZQlmowKvJuB0edG6d1AAcLpVON0q6prc+NvH+7EhtQx52X0xYlDyeb5roSdos/HeeuutuOGGG3Dttdfio48+wnvvvYe33367zTYffvgh1qxZgzfeeAN1dXWYPn06/ud//ge9e/cO+HU4G2/oCuf2dWXbCoursfKzQ/6RYG6vBlXVkBRjxoGyesjyyRUElZNrnhskJMVHwGI0wKsJeFUVw/onYteP1RAAFEmCy6uiockFh9t356QE35mOKoAJl/Y6LUQOltZi9ZYSKIoMoyKfDCwNowen4NtDVac9P+2KAaeFyNmO0S81GrsP10CCb4i0b0r2wOtwOL2AJCHCrJxWw+A+8XB5VJgsJpRXNsLh8qL4eD22H6iAyaAgNSECDc0eqKqGm68ZHJIh0u1n47XZbCgqKsLf/vY3AEB+fj6eeOIJ1NTUIDEx0b/dJ598ghkzZkCWZSQmJmLixIlYu3Yt7rzzzoBfK9xHYLB9oaur2ratqAKpiZEwGX66N8XtVVFV60BqQkSb+xpbfvWKMPtuyDSeXKxp1482xESZYDEqvmVpAVQaDZAl3/WL1m2raXIjIdbi6yLTfLPZllY2YVDveBhk2Te7LQCPR0NppR1902NhVGT/7LduVcXBsjoMH5TUprCio7VIT4qC0aDg5DS5cHtVWGscSEuIaFODpgkcLKtHbna/Nv8v9pXUnPb/wlbvG/GV1GopYLdXxb6SGgwbmASTUUF8fCTMJ3sId/9YjT5pMYiOMCHCrCAqwjfMeltRBUaFYHfW+XwugxIgVqsVaWlpUBTfm6YoClJTU2G1WtsEiNVqRWZmpv9xRkYGysvL2/VaCQlRHVN0N6X3N4VQEc7t66q2zb/9F13yuq397uYx532Mxwed/w/nBedxjPSkqPM+RrjhXTlERKRLUAIkIyMDFRUVUFXfUDtVVVFZWYmMjIzTtjtx4oT/sdVqRXp6ejBKJCKidgpKgCQlJSErKwsFBQUAgIKCAmRlZbXpvgKAvLw8rFq1CpqmoaamBhs2bEBubm4wSiQionYK2iis4uJizJs3Dw0NDYiNjcXSpUsxcOBAzJ49G3PnzsXw4cOhqioWL16MLVu2AABmz56NmTNnBqM8IiJqp6AFCBERhRdeRCciIl0YIEREpAsDhIiIdGGAEBGRLiE7G+8999yDY8eOQZZlREZG4o9//COysrICmrQxVLz00ktYvnw51qxZg8GDB+P777/H448/DpfLhV69euHZZ59FUlLSuQ/UDeXk5MBkMsFsNgMAHnroIVx55ZVh0UaXy4UlS5Zg69atMJvNGDVqFJ544omw+GweO3YMc+bM8T9ubGxEU1MTtm/fHhbtA4B//vOfWLZs2cnZigXuvfdeTJo0KWzat3HjRixbtgxerxdxcXF4+umn0adPH33tEyGqoaHB/+/PPvtMTJ8+XQghxC233CI+/PBDIYQQH374objlllu6pL7ztXfvXnHHHXeIq6++Whw8eFCoqiomTpwoduzYIYQQ4uWXXxbz5s3r4ir1a2lXa+HSxieeeEI89dRTQtM0IYQQVVVVQojw+Wy29uSTT4pFixYJIcKjfZqmiTFjxvg/m/v37xejRo0SqqqGRfvq6urEZZddJg4fPiyE8LXj9ttvF0Loe/9CNkBa++CDD8R1110nqqurxejRo4XX6xVCCOH1esXo0aOFzWbr4grbx+VyiRtvvFGUlZX5f9Du3r1bTJkyxb+NzWYTo0aN6sIqz8+ZAiQc2tjU1CRGjx4tmpqa2jwfLp/N1lwul8jOzhZ79+4Nm/ZpmiYuu+wysXPnTiGEENu3bxeTJk0Km/bt3r1bTJ482f+4trZWDB48WHf7QrYLCwAeffRRbNmyBUIIvPnmmwFP2tjdLVu2DNOmTWszjf2pE00mJiZC0zT/6WYoeuihhyCEwOjRo/G73/0uLNpYVlaG+Ph4vPTSS9i2bRuioqJw//33w2KxhMVns7UvvvgCaWlpuPjii7F3796waJ8kSXjhhRdwzz33IDIyEna7Ha+//nrY/GwZMGAAqqurUVhYiBEjRmDNmjUAAp/w9lQhfRH9qaeewsaNG/Hggw/imWee6epyOsSuXbuwd+9ezJo1q6tL6VQrV67E6tWr8d5770EIgcWLF3d1SR1CVVWUlZXhoosuwvvvv4+HHnoI9913H5qbm7u6tA733nvv4YYbbujqMjqU1+vFn//8Z7zyyiv45z//iVdffRUPPPBA2Lx/MTExeP755/H000/j+uuvh81mQ2xsrO72hXSAtJg+fTq2bduG9PT0gCZt7M527NiB4uJiTJgwATk5OSgvL8cdd9yBo0ePtplosqamBrIsh8xv5qdqeU9MJhNmzZqF77777rTJNEOxjRkZGTAYDMjPzwcAjBw5EgkJCbBYLCH/2WytoqICO3bswNSpUwEEPmFqd7d//35UVlZi9OjRAIDRo0cjIiICZrM5LNoHAJdffjneeecdvP/++/jXf/1XOJ1O9OrVS1f7QjJA7HY7rFar//EXX3yBuLi4gCdt7M7uuusubN68GV988QW++OILpKen4y9/+QvuvPNOOJ1O7Ny5EwDw7rvvIi8vr4ur1ae5uRmNjb7V+YQQ+OSTT5CVlYVhw4aFfBsTExORnZ3tn8+tpKQENpsN/fv3D/nPZmsffPABxo8fj4QE36qB4fC9BwDp6ekoLy/H4cOHAfjm8LPZbOjXr19YtA8AqqqqAACapuG5557DTTfdhF69eulqX0jOhVVdXY177rkHDocDsiwjLi4ODz/8MC6++OKzTtoYqnJycvDaa69h8ODB+O6777BgwYI2Q1yTk0NvCc2ysjLcd999UFUVmqZh0KBBeOyxx5CamhoWbSwrK8P8+fNRV1cHg8GABx54AOPHjw+rz2Zubi4effRRXHXVVf7nwqV9q1evxhtvvAFJ8q3UN3fuXEycODFs2vfoo4/iu+++g8fjwRVXXIH58+fDbDbral9IBggREXW9kOzCIiKirscAISIiXRggRESkCwOEiIh0YYAQEZEuDBAiItIlpOfCIjrVJZdc4v+3w+GAyWTyz++zaNEiTJs2ratK0y0nJwdPPvkkLr/88q4uhagNBgiFlV27dvn/HQo/eL1eLwyGzv02DMZrUM/ELizqETRNw+uvv46JEyciOzsb999/P+rq6gD4FkkaMmQI3nvvPYwfPx5jx47FO++8g8LCQkydOhVjxoxpM9nj+++/j5tuugmLFy/G6NGjkZeXh61bt/q/3tjYiPnz52PcuHG48sor8fzzz/vnGGrZd8mSJcjOzsby5ctRWlqKW2+9FdnZ2cjOzsbvf/97NDQ0AAD+8Ic/4MSJE/iP//gPXHLJJXjjjTewbdu2NneAA76w/PrrrwEAy5cvx9y5c/HQQw/h0ksvxQcffPCzNRHpxQChHuG///u/sWHDBqxYsQKbNm1CXFzcaTMA7969G+vXr8fzzz+PJUuW4LXXXsNbb72Fjz/+GJ9++im2b9/u37awsBB9+/bFN998g7lz5+Lee+/1B9K8efNgMBiwfv16fPjhh9iyZQtWrVrVZt8+ffpgy5YtuPvuuyGEwG9/+1ts2rQJn376KcrLy7F8+XIAwLPPPovMzEy89tpr2LVrF2bPnh1Qez///HPk5eVh586dmDp16jlrItKDAUI9wrvvvosHH3wQ6enpMJlMuPfee7Fu3Tp4vV7/NnPmzIHZbMa4ceMQGRmJ/Px8JCUlIS0tDWPGjEFRUZF/28TERPzbv/0bjEYjJk+ejAEDBmDjxo2orq7Gl19+ifnz5yMyMhJJSUn493//d3z88cf+fVNTU3HLLbfAYDDAYrGgX79+uOKKK2AymZCYmIjbbrsNO3bsOK/2jho1ChMnToQsy2hqajpnTUR6sGOUeoQTJ05gzpw5kOWffmeSZRk2m83/uPXa62az+bTHrddMSEtL80+2BwCZmZmorKzEiRMn4PV6MW7cOP/XNE1rMy12enp6m9qqq6vx1FNPYefOnbDb7RBCIDY29rza2/o1AqmJSA8GCPUI6enpWLJkiX+dh9aOHTvW7uNVVFRACOEPEavVipycHP8ZzjfffHPWC9etgwcAnnvuOUiShDVr1iA+Ph4bNmz42QW2IiIi4HQ6/Y9VVUVNTc1ZXyOQmoj0YBcW9Qi/+c1v8MILL+D48eMAfItVbdiwQffxampq8Pbbb8Pj8eDTTz9FcXExxo8fj9TUVFxxxRX405/+hKamJmiahtLS0jbXT05lt9sRGRmJmJgYVFRU4M0332zz9eTkZJSVlfkfDxgwAC6XCxs3boTH48Grr74Kt9t91uPrqYkoEAwQ6hFuvfVW5OTk4Pbbb8cll1yCG2+8EYWFhbqPN2LECBw9ehS/+MUv8MILL+DFF1/0L670zDPPwOPxYPLkyRg7dizmzp3rX8TnTO69914UFRVhzJgxuOuuuzBp0qQ2X7/rrrvw6quvYsyYMfjLX/6CmJgYLFiwAI899hiuuuoqREREnNYtdqr21kQUCK4HQtRO77//PlatWoV33nmnq0sh6lI8AyEiIl0YIEREpAu7sIiISBeegRARkS4MECIi0oUBQkREujBAiIhIFwYIERHpwgAhIiJd/h8ZPsfu8SN/bwAAAABJRU5ErkJggg==\n", + "text/plain": [ + "
" + ] + }, + "metadata": {}, + "output_type": "display_data" + } + ], + "source": [ + "sns.set(color_codes=True)\n", + "plt.xlim(30,90)\n", + "plt.ylim(0,1)\n", + "sns.regplot(x='Temperature', y='Frequency', data=data, logistic=True)\n", + "plt.show()" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "**I think I have managed to correctly compute and plot the uncertainty of my prediction.** Although the shaded area seems very similar to [the one obtained by with R](https://app-learninglab.inria.fr/moocrr/gitlab/moocrr-session3/moocrr-reproducibility-study/tree/master/challenger.pdf), I can spot a few differences (e.g., the blue point for temperature 63 is outside)... Could this be a numerical error ? Or a difference in the statistical method ? It is not clear which one is \"right\"." + ] + }, + { + "cell_type": "code", + "execution_count": 10, + "metadata": {}, + "outputs": [ + { + "data": { + "text/html": [ + "
\n", + "\n", + "\n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + "
TemperatureInterceptFrequency
030.011.0
130.511.0
231.011.0
331.511.0
432.011.0
............
11688.011.0
11788.511.0
11889.011.0
11989.511.0
12090.011.0
\n", + "

121 rows × 3 columns

\n", + "
" + ], + "text/plain": [ + " Temperature Intercept Frequency\n", + "0 30.0 1 1.0\n", + "1 30.5 1 1.0\n", + "2 31.0 1 1.0\n", + "3 31.5 1 1.0\n", + "4 32.0 1 1.0\n", + ".. ... ... ...\n", + "116 88.0 1 1.0\n", + "117 88.5 1 1.0\n", + "118 89.0 1 1.0\n", + "119 89.5 1 1.0\n", + "120 90.0 1 1.0\n", + "\n", + "[121 rows x 3 columns]" + ] + }, + "execution_count": 10, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "data_pred" + ] + }, + { + "cell_type": "code", + "execution_count": 11, + "metadata": {}, + "outputs": [ + { + "data": { + "text/plain": [ + "0 1.0\n", + "1 1.0\n", + "2 1.0\n", + "3 1.0\n", + "4 1.0\n", + " ... \n", + "116 1.0\n", + "117 1.0\n", + "118 1.0\n", + "119 1.0\n", + "120 1.0\n", + "Length: 121, dtype: float64" + ] + }, + "execution_count": 11, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "data_pred = pd.DataFrame({'Temperature': np.linspace(start=30, stop=90, num=121), 'Intercept': 1})\n", + "logmodel.predict(data_pred)" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [] + } + ], + "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.7.6" + }, + "toc": { + "base_numbering": 1, + "nav_menu": {}, + "number_sections": true, + "sideBar": true, + "skip_h1_title": false, + "title_cell": "Table of Contents", + "title_sidebar": "Contents", + "toc_cell": false, + "toc_position": {}, + "toc_section_display": true, + "toc_window_display": false + } + }, + "nbformat": 4, + "nbformat_minor": 2 +} diff --git a/module4/challenger.ipynb b/module4/challenger.ipynb index 8ce3c9d37ea640101486e2f776b953dcc1b370ea..1200c129f650f7c3fc443ec027f30c17b407c0a9 100644 --- a/module4/challenger.ipynb +++ b/module4/challenger.ipynb @@ -181,34 +181,14 @@ }, { "cell_type": "code", - "execution_count": 7, - "metadata": {}, - "outputs": [ - { - "data": { - "text/plain": [ - "'/home/aschmide/Documents/formation_MOOC_RR/mooc-rr/module4'" - ] - }, - "execution_count": 7, - "metadata": {}, - "output_type": "execute_result" - } - ], - "source": [ - "pwd" - ] - }, - { - "cell_type": "code", - "execution_count": 10, + "execution_count": 3, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ - "file downloaded on Sun Apr 12 17:22:23 2020\n" + "file downloaded on Tue Apr 14 10:06:35 2020\n" ] }, { @@ -455,15 +435,14 @@ "22 1/12/86 6 58 200 1" ] }, - "execution_count": 10, + "execution_count": 3, "metadata": {}, "output_type": "execute_result" } ], "source": [ "#data = pd.read_csv(\"https://app-learninglab.inria.fr/moocrr/gitlab/moocrr-session3/moocrr-reproducibility-study/blob/master/data/shuttle.csv\")\n", - "data_url = 'https://app-learninglab.inria.fr/moocrr/gitlab/moocrr-session3/moocrr-reproducibility-study/blob/master/data/shuttle.csv'\n", - "\n", + "data_url ='https://app-learninglab.inria.fr/moocrr/gitlab/moocrr-session3/moocrr-reproducibility-study/raw/7e01583b99527ad27cbeae0f9d2085fe8f2f1d15/data/shuttle.csv?inline=false'\n", "data_file='/home/aschmide/Documents/formation_MOOC_RR/mooc-rr/module4/shuttle.csv'\n", "import os\n", "import urllib.request\n", @@ -489,9 +468,22 @@ }, { "cell_type": "code", - "execution_count": null, + "execution_count": 4, "metadata": {}, - "outputs": [], + "outputs": [ + { + "data": { + "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYIAAAEKCAYAAAAfGVI8AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4yLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+j8jraAAAYAElEQVR4nO3dfZRddX3v8fdnkiEPJEJMvCMlQbBErlRyEYaHXGqd1IcGu0xq4Sp0IeoV4r2X3Cq6KtTrolyqaxVapVcvVVNKFbQOD7E4rWlRiqPYK5BEQoDQ4BQwTMAAYwgZCMlM5nv/OHv0zHBmZp/J7HPmnN/ntVZWzt77d/b5/s6eM5/ZD+e3FRGYmVm6WupdgJmZ1ZeDwMwscQ4CM7PEOQjMzBLnIDAzS5yDwMwscYUFgaQbJD0j6aExlkvSFyT1SNoq6ZSiajEzs7EVuUfwVWDlOMvPBpZm/9YAXyqwFjMzG0NhQRARPwR+MU6T1cCNUXIPcKSko4qqx8zMKptZx9c+GniybLo3m/f06IaS1lDaa2DOnDmnLlmypCYF5jU0NERLS/OdbmnWfkHz9s39ajy16tujjz76XES8ptKyegZBbhGxDlgH0N7eHps2bapzRSN1d3fT0dFR7zKmXLP2C5q3b+5X46lV3yT9bKxl9YzYnUD5n/aLs3lmZlZD9QyCLuDC7OqhM4E9EfGKw0JmZlaswg4NSfom0AEsktQL/AnQChARXwY2AO8CeoCXgA8VVYuZmY2tsCCIiPMnWB7AJUW9vpmZ5dOcp+HNzCw3B4GZWeIcBGZmiXMQmJklzkFgZpY4B4GZWeIcBGZmiXMQmJklzkFgZpY4B4GZWeIcBGZmiXMQmJklzkFgZpY4B4GZWeIcBGZmiXMQmJklzkFgZpY4B4GZWeIcBGZmiXMQmJklzkFgZpY4B4GZWeIcBGZmiXMQmJklzkFgZpY4B4GZWeIcBGZmiXMQmJklzkFgZpY4B4GZWeIcBGZmiXMQmJklzkFgZpY4B4GZWeIcBGZmiSs0CCStlLRdUo+kyyssP0bS9yXdL2mrpHcVWY+Zmb1SYUEgaQZwHXA2cCJwvqQTRzX7NHBLRLwZOA/4q6LqMTOzyorcIzgd6ImIxyLiANAJrB7VJoBXZY+PAJ4qsB4zM6tAEVHMiqVzgZURcVE2/X7gjIhYW9bmKOC7wALgcODtEbG5wrrWAGsA2traTu3s7Cyk5snq7+9n3rx59S5jyjVrv6B5++Z+NZ5a9W3FihWbI6K90rKZhb/6+M4HvhoRn5O0HLhJ0psiYqi8UUSsA9YBtLe3R0dHR+0rHUd3dzfTraap0Kz9gubtm/vVeKZD34o8NLQTWFI2vTibV+7DwC0AEfFjYDawqMCazMxslCKDYCOwVNJxkg6jdDK4a1SbHcDbACS9kVIQPFtgTWZmNkphQRARg8Ba4A7gEUpXBz0s6SpJq7JmnwAulvQA8E3gg1HUSQszM6uo0HMEEbEB2DBq3hVlj7cBZxVZg5mZjc/fLDYzS5yDwMwscQ4CM7PEOQjMzBLnIDAzS5yDwMwscQ4CM7PEOQjMzBLnIDAzS5yDwMwscQ4CM7PEOQjMzBLnIDAzS5yDwMwscQ4CM7PEOQjMzBLnIDAzS5yDwMwscQ4CM7PEOQjMzBLnIDAzS5yDwMwscQ4CM7PEOQjMzBLnIDAzS5yDwMwscQ4CM7PEOQjMzBLnIDAzS5yDwMwscQ4CM7PEOQjMzBLnIDAzS5yDwMwscYUGgaSVkrZL6pF0+Rht3itpm6SHJf1dkfWYmdkrzSxqxZJmANcB7wB6gY2SuiJiW1mbpcAfA2dFxG5J/6GoeszMrLJcewSSTprEuk8HeiLisYg4AHQCq0e1uRi4LiJ2A0TEM5N4HTMzOwSKiIkbSXcDs4CvAt+IiD05nnMusDIiLsqm3w+cERFry9rcDjwKnAXMAK6MiH+usK41wBqAtra2Uzs7OyfuWQ319/czb968epcx5Zq1X9C8fXO/Gk+t+rZixYrNEdFeaVmuQ0MR8ZbsMM5/BTZLug/424j43iHWNhNYCnQAi4EfSjopIp4f9frrgHUA7e3t0dHRcYgvO7W6u7uZbjVNhWbtFzRv39yvxjMd+pb7ZHFE/BT4NHAZ8FbgC5L+TdLvj/GUncCSsunF2bxyvUBXRAxExOOU9g6W5q3JzMwOXd5zBMskXQs8Avw28O6IeGP2+NoxnrYRWCrpOEmHAecBXaPa3E5pbwBJi4A3AI9V2wkzM5u8vFcNfRG4HvhUROwbnhkRT0n6dKUnRMSgpLXAHZSO/98QEQ9LugrYFBFd2bJ3StoGHAT+KCL6DqE/ZmZWpbxB8LvAvog4CCCpBZgdES9FxE1jPSkiNgAbRs27ouxxAB/P/pmZWR3kPUdwJzCnbHpuNs/MzBpc3iCYHRH9wxPZ47nFlGRmZrWUNwhelHTK8ISkU4F947Q3M7MGkfccwceAWyU9BQh4LfC+wqoyM7OayfuFso2S/iNwQjZre0QMFFeWmZnVSjWDzp0GHJs95xRJRMSNhVRlZmY1kysIJN0E/DqwhdL1/gABOAjMzBpc3j2CduDEyDNCnZmZNZS8Vw09ROkEsZmZNZm8ewSLgG3ZqKP7h2dGxKpCqjIzs5rJGwRXFlmEmZnVT97LR38g6XXA0oi4U9JcSgPJmZlZg8s7DPXFwG3AV7JZR1MaQtrMzBpc3pPFl1C6neQL8Mub1PhG82ZmTSBvEOzPbkAPgKSZlL5HYGZmDS5vEPxA0qeAOZLeAdwK/ENxZZmZWa3kDYLLgWeBB4GPULrZTMU7k5mZWWPJe9XQEPDX2T8zM2sieccaepwK5wQi4vVTXpGZmdVUNWMNDZsN/Bfg1VNfjpmZ1VqucwQR0Vf2b2dE/CWlG9qbmVmDy3to6JSyyRZKewjV3MvAzMymqby/zD9X9ngQeAJ475RXY2ZmNZf3qqEVRRdiZmb1kffQ0MfHWx4Rn5+acszMrNaquWroNKArm343cB/w0yKKMjOz2skbBIuBUyJiL4CkK4HvRMQFRRVmZma1kXeIiTbgQNn0gWyemZk1uLx7BDcC90n6+2z694CvFVOSmZnVUt6rhj4r6Z+At2SzPhQR9xdXlpmZ1UreQ0MAc4EXIuL/AL2SjiuoJjMzq6G8t6r8E+Ay4I+zWa3A14sqyszMaifvHsF7gFXAiwAR8RQwv6iizMysdvIGwYGICLKhqCUdXlxJZmZWS3mD4BZJXwGOlHQxcCe+SY2ZWVOYMAgkCbgZuA1YD5wAXBERX8zx3JWStkvqkXT5OO3OkRSS2sdqY2ZmxZjw8tGICEkbIuIk4Ht5VyxpBnAd8A6gF9goqSsito1qNx/4KHBvVZWbmdmUyHto6CeSTqty3acDPRHxWEQcADqB1RXa/SlwNfByles3M7MpkPebxWcAF0h6gtKVQ6K0s7BsnOccDTxZNt2breeXshveLImI70j6o7FWJGkNsAagra2N7u7unGXXRn9//7SraSo0a7+gefvmfjWe6dC3cYNA0jERsQP4nal+YUktwOeBD07UNiLWAesA2tvbo6OjY6rLOSTd3d1Mt5qmQrP2C5q3b+5X45kOfZtoj+B2SqOO/kzS+og4p4p17wSWlE0vzuYNmw+8CegunY/mtUCXpFURsamK1zEzs0Mw0TkClT1+fZXr3ggslXScpMOA8/jV/QyIiD0RsSgijo2IY4F7AIeAmVmNTRQEMcbjCUXEILAWuAN4BLglIh6WdJWkVdWVaWZmRZno0NB/kvQCpT2DOdlj+NXJ4leN9+SI2ABsGDXvijHaduSq2MzMptS4QRARM2pViJmZ1Uc1w1CbmVkTchCYmSXOQWBmljgHgZlZ4pIJgr7+/Tzw5PP09e+vdylmVoW+/v3sGzjoz26BkgiCb2/ZyVlX38UF19/LWVffRdeWnRM/yczqbviz+/izL/qzW6CmD4K+/v1ctn4rLw8MsXf/IC8PDPHJ9Vv914XZNFf+2T0Y4c9ugZo+CHp376O1ZWQ3W1ta6N29r04VmVke/uzWTtMHweIFcxgYGhoxb2BoiMUL5tSpIjPLw5/d2mn6IFg4bxbXnLOM2a0tzJ81k9mtLVxzzjIWzptV79LMbBzln90Zkj+7Bcp7Y5qGturkoznr+EX07t7H4gVz/INk1iCGP7v3/fhH/Ouq3/RntyBJBAGU/rrwD5FZ41k4bxZzWmf481ugpj80ZGZm43MQmJklzkFgZpY4B4GZWeIcBGZmiXMQmJklzkFgZpY4B4GZWeIcBGZmiXMQmJklzkFgZpY4B4GZWeIcBGZmiXMQmJklzkFgZpY4B4GZWeIcBGZmiXMQmJklzkFgZpY4B4GZWeIcBGZmiSs0CCStlLRdUo+kyyss/7ikbZK2SvoXSa8rsh4zM3ulwoJA0gzgOuBs4ETgfEknjmp2P9AeEcuA24BriqrHzMwqK3KP4HSgJyIei4gDQCewurxBRHw/Il7KJu8BFhdYj5mZVaCIKGbF0rnAyoi4KJt+P3BGRKwdo/3/BX4eEZ+psGwNsAagra3t1M7OzkJqnqz+/n7mzZtX7zKmXLP2C5q3b+5X46lV31asWLE5ItorLZtZ+KvnIOkCoB14a6XlEbEOWAfQ3t4eHR0dtSsuh+7ubqZbTVOhWfsFzds396vxTIe+FRkEO4ElZdOLs3kjSHo78L+At0bE/gLrMTOzCoo8R7ARWCrpOEmHAecBXeUNJL0Z+AqwKiKeKbAWMzMbQ2FBEBGDwFrgDuAR4JaIeFjSVZJWZc3+HJgH3Cppi6SuMVZnZmYFKfQcQURsADaMmndF2eO3F/n6jayvfz+9u/exeMEcFs6bNWVtG0mz9qsoPbv2svulAXp27eX4tvn1LscayLQ4WWwjfXvLTi5bv5XWlhYGhoa45pxlrDr56ENu20iatV9FueL2B7nxnh184qRBLr32h1y4/BiuWn1SvcuyBuEhJqaZvv79XLZ+Ky8PDLF3/yAvDwzxyfVb6et/5Xn0ato2kmbtV1F6du3lxnt2jJh344930LNrb50qskbjIJhmenfvo7Vl5GZpbWmhd/e+Q2rbSJq1X0XZ8uTzVc03G81BMM0sXjCHgaGhEfMGhoZYvGDOIbVtJM3ar6KcvOTIquabjeYgmGYWzpvFNecsY3ZrC/NnzWR2awvXnLOs4snSato2kmbtV1GOb5vPhcuPGTHvwuXH+ISx5eaTxdPQqpOP5qzjF+W6Yqaato2kWftVlKtWn8SFZx7Lg5vv4c5Lz3QIWFUcBNPUwnmzcv/yq6ZtI2nWfhXl+Lb59M5tdQhY1XxoyMwscQ4CM7PEOQjMzBLnIDAzS5yDwMwscQ4CM7PEOQjMzBLnIDAzS5yDwMwscQ4CM7PEOQjMzBLnIDAzS5yDwMwscQ4CM7PEOQjMzBLnIDAzS5yDwMwscQ4CM7PEOQjMzBLnIDAzS5yDwMwscQ4CM7PEOQjMzBLnIDAzS5yDwMwscQ4CM7PEOQjMzBLnIDAzS1yhQSBppaTtknokXV5h+SxJN2fL75V0bJH1mFWrr38/Dzz5PH39+8dtt+nxPj7/3e1serxvytZZbdueXXvZ/dIAPbv2Tti2GkXVW83r7xs4mPs9uG3Tk033HhS5XoCZU77GjKQZwHXAO4BeYKOkrojYVtbsw8DuiDhe0nnA1cD7iqrJrBrf3rKTy9ZvpbWlhYGhIa45ZxmrTj76Fe0uuP4eftRTCoAv3NXDW45fyE0XnXlI66y27RW3P8iN9+zgEycNcum1P+TC5cdw1eqTJtnz4uut9vX/8I0DXHr1Xbneg2HN8h4Uud5hRe4RnA70RMRjEXEA6ARWj2qzGvha9vg24G2SVGBNZrn09e/nsvVbeXlgiL37B3l5YIhPrt/6ir/GNj3e98sQGHZ3T1/FPYO866y2bc+uvSN+AQLc+OMdh/xXcVH1Tub1D0Yk+R4Uud5yiogpW9mIFUvnAisj4qJs+v3AGRGxtqzNQ1mb3mz637M2z41a1xpgTTZ5ArC9kKInbxHw3IStGk+z9gsm6Jta58ydueCoN6ilZcbwvBgaOji4++lHY2DfS8PzZsxf9GszDj/yqNHPP/ji808f3PvcU5NZZ7VtW+YesXDmq15zLMDBl/YwY+4RAAy+8OwTQy/tmfhY1SG+B9W2nczrD/crz3tQrkHegyn5WczhdRHxmkoLCjs0NJUiYh2wrt51jEXSpohor3cdU61Z+wXN2zdJmwb3PON+NZDp8LNY5KGhncCSsunF2byKbSTNBI4AJp3eZmZWvSKDYCOwVNJxkg4DzgO6RrXpAj6QPT4XuCuKOlZlZmYVFXZoKCIGJa0F7gBmADdExMOSrgI2RUQX8DfATZJ6gF9QCotGNG0PWx2iZu0XNG/f3K/GU/e+FXay2MzMGoO/WWxmljgHgZlZ4hwEVZL0hKQHJW2RtCmbd6Wkndm8LZLeVe86J0PSkZJuk/Rvkh6RtFzSqyV9T9JPs/8X1LvOao3Rr4bfZpJOKKt/i6QXJH2s0bfZOP1qhm12qaSHJT0k6ZuSZmcX1NybDbVzc3ZxTW3r8jmC6kh6Amgv/9KbpCuB/oj4i3rVNRUkfQ24OyKuz34Y5wKfAn4REX+WjRe1ICIuq2uhVRqjXx+jCbbZsGxIl53AGcAlNPg2GzaqXx+igbeZpKOBHwEnRsQ+SbcAG4B3Ad+KiE5JXwYeiIgv1bI27xEYAJKOAH6L0pVcRMSBiHiekcOAfA34vfpUODnj9KvZvA3494j4GQ2+zUYp71czmAnMyb43NRd4GvhtSkPsQJ22l4OgegF8V9LmbOiLYWslbZV0Q6PtimeOA54F/lbS/ZKul3Q40BYRT2dtfg601a3CyRmrX9D426zcecA3s8eNvs3KlfcLGnibRcRO4C+AHZQCYA+wGXg+IgazZr3A1I0ml5ODoHq/GRGnAGcDl0j6LeBLwK8DJ1PawJ+rY32TNRM4BfhSRLwZeBEYMXR49mW/RjuWOFa/mmGbAZAd7loF3Dp6WYNuM6Bivxp6m2XBtZrSHye/BhwOrKxrURkHQZWyVCcingH+Hjg9InZFxMGIGAL+mtLIq42mF+iNiHuz6dso/QLdJekogOz/Z+pU32RV7FeTbLNhZwM/iYhd2XSjb7NhI/rVBNvs7cDjEfFsRAwA3wLOAo7MDhVB5aF4CucgqIKkwyXNH34MvBN4aPhDl3kP8FA96jsUEfFz4ElJJ2Sz3gZsY+QwIB8Avl2H8iZtrH41wzYrcz4jD5809DYrM6JfTbDNdgBnSporSfzqM/Z9SkPsQJ22l68aqoKk11PaC4DSIYe/i4jPSrqJ0u5qAE8AHyk7RtswJJ0MXA8cBjxG6SqNFuAW4BjgZ8B7I+IXdStyEsbo1xdojm12OKVfMK+PiD3ZvIU0/jar1K+G/5xJ+t+Ubr41CNwPXETpnEAn8Ops3gURMfW3IRuvLgeBmVnafGjIzCxxDgIzs8Q5CMzMEucgMDNLnIPAzCxxDXHzerO8sksn/yWbfC1wkNIQE1D68t+BuhRWgaQO4EBE/L9612JpcxBYU4mIPkrXmk+LUWElzSwbR2a0DqAfyB0EE6zPbFJ8aMianqRTJf0gGyjwjrLhF7olXStpU3afgtMkfSsbx/8zWZtjs/sYfCNrc5ukuTnW+5cq3a/io5LenY03f7+kOyW1SToW+G/ApdnY+m+R9FVJ55bV3Z/93yHpbkldlL4VPUPSn0vamA3A9pFavp/WfBwE1uwEfBE4NyJOBW4APlu2/EBEtANfpvTV/kuANwEfzA4zAZwA/FVEvBF4AfgfklonWO9hEdEeEZ+jNAb9mdmgd53AJyPiiew1r42IkyPi7gn6cQrw0Yh4A/BhYE9EnAacBlws6bjq3xqzEh8asmY3i9Iv9u+VhndhBqWRK4d1Zf8/CDw8PGSBpMeAJcDzwJMR8a9Zu68Dfwj88wTrvbns8WLg5myP4TDg8Un0476IGH7eO4FlZXsPRwBLJ7leMweBNT1R+gW/fIzlw2O6DJU9Hp4e/nyMHoclcqz3xbLHXwQ+HxFd2QniK8d4ziDZXrqkFkqhUWl9Av5nRNwxxnrMquJDQ9bs9gOvkbQcQFKrpN+och3HDD8f+ANKh3q2V7HeI/jV0MIfKJu/F5hfNv0EcGr2eBXQOsb67gD+e3Z4CklvKLvZjlnVHATW7IYoDfF7taQHgC3Af65yHdsp3YToEWABpZvcHKhivVcCt0raDDxXNv8fgPcMnyymNMb+W7P1LWfkXkC56ykNX/wTSQ8BX8F793YIPPqo2Tiyq3v+MSLeVOdSzArjPQIzs8R5j8DMLHHeIzAzS5yDwMwscQ4CM7PEOQjMzBLnIDAzS9z/Bxqb4z828mHUAAAAAElFTkSuQmCC\n", + "text/plain": [ + "
" + ] + }, + "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", @@ -513,9 +505,91 @@ }, { "cell_type": "code", - "execution_count": null, + "execution_count": 5, "metadata": {}, - "outputs": [], + "outputs": [ + { + "name": "stderr", + "output_type": "stream", + "text": [ + "/home/aschmide/miniconda3/lib/python3.7/site-packages/ipykernel_launcher.py:7: DeprecationWarning: Calling Family(..) with a link class as argument is deprecated.\n", + "Use an instance of a link class instead.\n", + " import sys\n" + ] + }, + { + "data": { + "text/html": [ + "\n", + "\n", + "\n", + " \n", + "\n", + "\n", + " \n", + "\n", + "\n", + " \n", + "\n", + "\n", + " \n", + "\n", + "\n", + " \n", + "\n", + "\n", + " \n", + "\n", + "\n", + " \n", + "\n", + "\n", + " \n", + "\n", + "\n", + " \n", + "\n", + "
Generalized Linear Model Regression Results
Dep. Variable: Frequency No. Observations: 23
Model: GLM Df Residuals: 21
Model Family: Binomial Df Model: 1
Link Function: logit Scale: 1.0000
Method: IRLS Log-Likelihood: -3.9210
Date: Tue, 14 Apr 2020 Deviance: 3.0144
Time: 10:27:51 Pearson chi2: 5.00
No. Iterations: 6
Covariance Type: nonrobust
\n", + "\n", + "\n", + " \n", + "\n", + "\n", + " \n", + "\n", + "\n", + " \n", + "\n", + "
coef std err z P>|z| [0.025 0.975]
Intercept 5.0850 7.477 0.680 0.496 -9.570 19.740
Temperature -0.1156 0.115 -1.004 0.316 -0.341 0.110
" + ], + "text/plain": [ + "\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: Tue, 14 Apr 2020 Deviance: 3.0144\n", + "Time: 10:27:51 Pearson chi2: 5.00\n", + "No. Iterations: 6 \n", + "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": 5, + "metadata": {}, + "output_type": "execute_result" + } + ], "source": [ "import statsmodels.api as sm\n", "\n", @@ -537,9 +611,91 @@ }, { "cell_type": "code", - "execution_count": null, + "execution_count": 6, "metadata": {}, - "outputs": [], + "outputs": [ + { + "name": "stderr", + "output_type": "stream", + "text": [ + "/home/aschmide/miniconda3/lib/python3.7/site-packages/ipykernel_launcher.py:2: DeprecationWarning: Calling Family(..) with a link class as argument is deprecated.\n", + "Use an instance of a link class instead.\n", + " \n" + ] + }, + { + "data": { + "text/html": [ + "\n", + "\n", + "\n", + " \n", + "\n", + "\n", + " \n", + "\n", + "\n", + " \n", + "\n", + "\n", + " \n", + "\n", + "\n", + " \n", + "\n", + "\n", + " \n", + "\n", + "\n", + " \n", + "\n", + "\n", + " \n", + "\n", + "\n", + " \n", + "\n", + "
Generalized Linear Model Regression Results
Dep. Variable: Frequency No. Observations: 23
Model: GLM Df Residuals: 21
Model Family: Binomial Df Model: 1
Link Function: logit Scale: 1.0000
Method: IRLS Log-Likelihood: -23.526
Date: Tue, 14 Apr 2020 Deviance: 18.086
Time: 10:27:51 Pearson chi2: 30.0
No. Iterations: 6
Covariance Type: nonrobust
\n", + "\n", + "\n", + " \n", + "\n", + "\n", + " \n", + "\n", + "\n", + " \n", + "\n", + "
coef std err z P>|z| [0.025 0.975]
Intercept 5.0850 3.052 1.666 0.096 -0.898 11.068
Temperature -0.1156 0.047 -2.458 0.014 -0.208 -0.023
" + ], + "text/plain": [ + "\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: -23.526\n", + "Date: Tue, 14 Apr 2020 Deviance: 18.086\n", + "Time: 10:27:51 Pearson chi2: 30.0\n", + "No. Iterations: 6 \n", + "Covariance Type: nonrobust \n", + "===============================================================================\n", + " coef std err z P>|z| [0.025 0.975]\n", + "-------------------------------------------------------------------------------\n", + "Intercept 5.0850 3.052 1.666 0.096 -0.898 11.068\n", + "Temperature -0.1156 0.047 -2.458 0.014 -0.208 -0.023\n", + "===============================================================================\n", + "\"\"\"" + ] + }, + "execution_count": 6, + "metadata": {}, + "output_type": "execute_result" + } + ], "source": [ "logmodel=sm.GLM(data['Frequency'], data[['Intercept','Temperature']], \n", " family=sm.families.Binomial(sm.families.links.logit),\n", @@ -568,9 +724,22 @@ }, { "cell_type": "code", - "execution_count": null, + "execution_count": 7, "metadata": {}, - "outputs": [], + "outputs": [ + { + "data": { + "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXQAAAEKCAYAAAACS67iAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4yLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+j8jraAAAZqUlEQVR4nO3df3TV9Z3n8efbECT8ECwqIwaF2SKOg/Ij4VdxbGirYNtB3GFE6tDaU0p3t7Q6jsyRM93qOHrOduOOzjqOI6usnXU1IMdJsYfT0Dpk7bqrBgqIwAbQppJoB8XyIzZICO/94/u96SUkuTc39+be++H1OCcn9/u9n+/3+3nnS15887nf+7nm7oiISPE7L98dEBGR7FCgi4gEQoEuIhIIBbqISCAU6CIigVCgi4gEImWgm9laMztkZm/18LyZ2X81swNm9qaZTc9+N0VEJJV0rtCfARb08vxNwMT4awXwRP+7JSIifZUy0N39FeCjXprcDPyTR14DRpnZpdnqoIiIpGdQFvZxGXAwabk5Xvd+14ZmtoLoKp6ysrKKcePG9flgH51wPulwLLO+FhwH1VJgQqkDVEuhKj0PRpdl9hLmvn37PnT3i7t7LhuBnjZ3XwOsAaisrPStW7dmtJ/6+nqqqqqy2LP8US2FJ5Q6QLUUqv7UYma/6um5bNzl0gIkX2qXx+tERGQAZSPQNwJfje92mQ0cdfezhltERCS3Ug65mNnzQBVwkZk1A/cBpQDu/o/AJuCLwAHgt8DXc9VZERHpWcpAd/elKZ534NtZ65GIFI329naam5s5ceJEzo81cuRI9u7dm/PjDIR0ahkyZAjl5eWUlpamvd8BfVFURMLS3NzMiBEjGD9+PGa5vQfl+PHjjBgxIqfHGCipanF3Dh8+THNzMxMmTEh7v3rrv4hk7MSJE4wePTrnYX6uMTNGjx7d5798FOgi0i8K89zI5OeqQBcRCYTG0EWkqJWUlHDNNdd0LtfW1jJ+/Pj8dSiPFOgiUtTKysrYsWNHt8+5O+7OeeedG4MR50aVInLOaGpqYtKkSXz1q19l8uTJHDx4kOrqambMmMG1117Lfffd19n2oYce4sorr+S6665j6dKlPPzwwwBUVVWRmJrkww8/7Lzi7+joYNWqVZ37evLJJ4HfvZV/8eLFXHXVVdx+++1Ed3RDQ0MDn/nMZ5gyZQozZ87k+PHjLFiw4Iz/hK677jp27tzZ79p1hS4iWfHXL+1mz3vHsrrPq8dewH1//Ie9tmlra2Pq1KkATJgwgUceeYT9+/fzwx/+kNmzZ7N582b279/PG2+8gbuzcOFCXnnlFYYNG0ZNTQ07duzg1KlTTJ8+nYqKil6P9fTTTzNy5EgaGhr45JNPmDt3LjfeeCMA27dvZ/fu3YwdO5a5c+fy6quvMnPmTJYsWcK6deuYMWMGx44do6ysjGXLlvHMM8/w6KOPsm/fPk6cOMGUKVP6/fNSoItIUes65NLU1MQVV1zB7NmzAdi8eTObN29m2rRpALS2trJ//36OHz/OLbfcwtChQwFYuHBhymNt3ryZN998kw0bNgBw9OhR9u/fz+DBg5k5cybl5eUATJ06laamJkaOHMmll17KjBkzALjgggsAuOWWW5g7dy7V1dWsXbuWO+64Iys/CwW6iGRFqivpgTRs2LDOx+7O6tWr+da3vnVGm0cffbTH7QcNGsTp06cBzrgX3N157LHHmD9//hnt6+vrOf/88zuXS0pKOHXqVI/7Hzp0KDfccAM/+tGPWL9+Pdu2bUuvsBQ0hi4iQZs/fz5r166ltbUVgJaWFg4dOsT1119PbW0tbW1tHD9+nJdeeqlzm/Hjx3eGbOJqPLGvJ554gvb2dgD27dvHxx9/3OOxJ02axPvvv09DQwMQvUM0EfTLly/nu9/9LjNmzODCCy/MSq26QheRoN14443s3buXOXPmADB8+HCeffZZpk+fzpIlS5gyZQqXXHJJ57AIwD333MOtt97KmjVr+NKXvtS5fvny5TQ1NTF9+nTcnYsvvpja2toejz148GDWrVvHd77zHdra2igrK+NnP/sZABUVFVxwwQV8/etZnM8wcVvPQH9VVFR4prZs2ZLxtoVGtRSeUOpwz30te/bsyen+kx07diyn+7/vvvu8uro6p8dIOHbsmLe0tPjEiRO9o6Ojx3bd/XyBrd5DrmrIRURkgD333HPMmjWLhx56KKv3yGvIRUQEuP/++wfsWF/5ylfOepE2G3SFLiL94vEbaCS7Mvm5KtBFJGNDhgzh8OHDCvUs83g+9CFDhvRpOw25iEjGysvLaW5u5oMPPsj5sU6cONHngCtU6dSS+MSivlCgi0jGSktL+/SJOv1RX1/f+W7PYperWjTkIiISCAW6iEggFOgiIoFQoIuIBEKBLiISCAW6iEggFOgiIoFQoIuIBEKBLiISCAW6iEggFOgiIoFQoIuIBEKBLiISCAW6iEggFOgiIoFQoIuIBCKtQDezBWbWaGYHzOzebp6/3My2mNl2M3vTzL6Y/a6KiEhvUga6mZUAjwM3AVcDS83s6i7Nvgesd/dpwG3AP2S7oyIi0rt0rtBnAgfc/R13PwnUADd3aePABfHjkcB72euiiIikw1J9WreZLQYWuPvyeHkZMMvdVya1uRTYDFwIDAO+4O7butnXCmAFwJgxYypqamoy6nRrayvDhw/PaNtCo1oKTyh1gGopVP2pZd68edvcvbLbJ9291y9gMfBU0vIy4O+7tLkb+Iv48RxgD3Beb/utqKjwTG3ZsiXjbQuNaik8odThrloKVX9qAbZ6D7mazpBLCzAuabk8XpfsG8D6+D+I/wsMAS5KY98iIpIl6QR6AzDRzCaY2WCiFz03dmnzLvB5ADP7A6JA/yCbHRURkd6lDHR3PwWsBOqAvUR3s+w2swfMbGHc7C+Ab5rZTuB54I74TwMRERkgg9Jp5O6bgE1d1n0/6fEeYG52uyYiIn2hd4qKiARCgS4iEggFuohIIBToIiKBUKCLiARCgS4iEggFuohIIBToIiKBUKCLiARCgS4iEggFuohIIBToIiKBUKCLiARCgS4iEggFuohIIBToIiKBUKCLiARCgS4iEggFuohIIBToIiKBUKCLiARCgS4iEggFuohIIBToIiKBUKCLiARCgS4iEggFuohIIBToIiKBUKCLiARCgS4iEggFuohIIBToIiKBUKCLiARCgS4iEoi0At3MFphZo5kdMLN7e2hzq5ntMbPdZvZcdrspIiKpDErVwMxKgMeBG4BmoMHMNrr7nqQ2E4HVwFx3/42ZXZKrDouISPfSuUKfCRxw93fc/SRQA9zcpc03gcfd/TcA7n4ou90UEZFUzN17b2C2GFjg7svj5WXALHdfmdSmFtgHzAVKgPvd/Sfd7GsFsAJgzJgxFTU1NRl1urW1leHDh2e0baFRLYUnlDpAtRSq/tQyb968be5e2d1zKYdc0jQImAhUAeXAK2Z2jbsfSW7k7muANQCVlZVeVVWV0cHq6+vJdNtCo1oKTyh1gGopVLmqJZ0hlxZgXNJyebwuWTOw0d3b3f2XRFfrE7PTRRERSUc6gd4ATDSzCWY2GLgN2NilTS3R1TlmdhFwJfBOFvspIiIppAx0dz8FrATqgL3AenffbWYPmNnCuFkdcNjM9gBbgFXufjhXnRYRkbOlNYbu7puATV3WfT/psQN3x18iIpIHeqeoiEggFOgiIoFQoIuIBEKBLiISCAW6iEggFOgiIoFQoIuIBEKBLiISCAW6iEggFOgiIoFQoIuIBEKBLiISCAW6iEggFOgiIoFQoIuIBEKBLiISCAW6iEggFOgiIoFQoIuIBEKBLiISCAW6iEggFOgiIoFQoIuIBEKBLiISCAW6iEggFOgiIoFQoIuIBEKBLiISCAW6iEggFOgiIoFQoIuIBEKBLiISCAW6iEggFOgiIoFQoIuIBCKtQDezBWbWaGYHzOzeXtr9iZm5mVVmr4siIpKOlIFuZiXA48BNwNXAUjO7upt2I4A7gdez3UkREUktnSv0mcABd3/H3U8CNcDN3bT7G+AHwIks9k9ERNJk7t57A7PFwAJ3Xx4vLwNmufvKpDbTgb9y9z8xs3rgHnff2s2+VgArAMaMGVNRU1OTUadbW1sZPnx4RtsWGtVSeEKpA1RLoepPLfPmzdvm7t0Oaw/qV68AMzsP+FvgjlRt3X0NsAagsrLSq6qqMjpmfX09mW5baFRL4QmlDlAthSpXtaQz5NICjEtaLo/XJYwAJgP1ZtYEzAY26oVREZGBlU6gNwATzWyCmQ0GbgM2Jp5096PufpG7j3f38cBrwMLuhlxERCR3Uga6u58CVgJ1wF5gvbvvNrMHzGxhrjsoIiLpSWsM3d03AZu6rPt+D22r+t8tERHpK71TVEQkEAp0EZFAKNBFRAKhQBcRCYQCXUQkEP1+p6hIf9Rub6G6rpH3jrQxdlQZq+ZPYtG0y/LdLUmTzl9hUaBL3tRub2H1i7toa+8AoOVIG6tf3AWgUCgCOn+FR0MukjfVdY2dYZDQ1t5BdV1jnnokfaHzV3gU6JI37x1p69N6KSw6f4VHgS55M3ZUWZ/WS2HR+Ss8CnTJm1XzJ1FWWnLGurLSElbNn5SnHklf6PwVHr0oKnmTeOFMd0kUJ52/wqNAl7xaNO0yBUAR0/krLBpyEREJhAJdRCQQCnQRkUAo0EVEAqFAFxEJhAJdRCQQCnQRkUAo0EVEAqFAFxEJhAJdRCQQCnQRkUAo0EVEAqFAFxEJhAJdRCQQCnQRkUAo0EVEAqFAFxEJhAJdRCQQCnQRkUAo0EVEAqFAFxEJRFqBbmYLzKzRzA6Y2b3dPH+3me0xszfN7GUzuyL7XRURkd6kDHQzKwEeB24CrgaWmtnVXZptByrd/VpgA/Cfs91RERHpXTpX6DOBA+7+jrufBGqAm5MbuPsWd/9tvPgaUJ7dboqISCrm7r03MFsMLHD35fHyMmCWu6/sof3fA7929we7eW4FsAJgzJgxFTU1NRl1urW1leHDh2e0baFRLYUnlDpAtRSq/tQyb968be5e2d1zg/rVqy7M7M+ASuCz3T3v7muANQCVlZVeVVWV0XHq6+vJdNtCo1oKTyh1gGopVLmqJZ1AbwHGJS2Xx+vOYGZfAP4K+Ky7f5Kd7omISLrSGUNvACaa2QQzGwzcBmxMbmBm04AngYXufij73RQRkVRSBrq7nwJWAnXAXmC9u+82swfMbGHcrBoYDrxgZjvMbGMPuxMRkRxJawzd3TcBm7qs+37S4y9kuV8iGand3kJ1XSPvHWlj7KgyVs2fBHDWukXTLhuQY+fiOOn4Xu0unn/9IHdNbucbqzexdNY4Hlx0TV76IgMnqy+KiuRT7fYWVr+4i7b2DgBajrSx6oWdYNDe4Z3rVr+4CyCrYdvdsXNxnHR8r3YXz772budyh3vnskI9bHrrvwSjuq6xM1AT2k97Z5gntLV3UF3XmPNj5+I46Xj+9YN9Wi/hUKBLMN470paTtv3ZX7aPk46OHt5b0tN6CYcCXYIxdlRZTtr2Z3/ZPk46Ssz6tF7CoUCXYKyaP4my0pIz1pWeZ5SWnBlkZaUlnS+W5vLYuThOOpbOGten9RIOvSgqwUi8+JiPu1x6OnY+7nJJvPCZGDMvMdNdLucIBboEZdG0y7oN0YEI1p6OnQ8PLrqGBxddQ319PW/fXpXv7sgA0ZCLiEggFOgiIoFQoIuIBEKBLiISCAW6iEggFOgiIoFQoIuIBEKBLiISCAW6iEggFOgiIoFQoIuIBEKBLiISCAW6iEggFOgiIoFQoIuIBEKBLiISCAW6iEggFOgiIoFQoIuIBEKBLiISCAW6iEggFOgiIoFQoIuIBEKBLiISCAW6iEggFOgiIoFQoIuIBEKBLiISiEHpNDKzBcDfASXAU+7+n7o8fz7wT0AFcBhY4u5N2e2qSLhqt7dQXdfIe0faGDuqjFXzJ/HC1nd59e2POtvM/Tef4k8rLz+rHXDWuq2/+ojnXz/IXZPb+cbqTSydNY4HF12T1nG729+iaZel3e/EsTvcKTHLybG727anPp5LUga6mZUAjwM3AM1Ag5ltdPc9Sc2+AfzG3T9tZrcBPwCW5KLDIqGp3d7C6hd30dbeAUDLkTbuWrfjrHavvv3RGQHfcqSNVRt2gkP7ae9cd/e6HZxO2q7DnWdfexfgjGDt7rirXtgJBu0dv9vf6hd3AZwVmN1tPxDH7m7bnvp4rklnyGUmcMDd33H3k0ANcHOXNjcDP4wfbwA+b2aWvW6KhKu6rrEznPqqvcM7wzzhdA9tn3/9YMrjtp/2zkBNaGvvoLqu8az9dbf9QBy7u2176uO5xty99wZmi4EF7r48Xl4GzHL3lUlt3orbNMfLb8dtPuyyrxXAinhxEpDpGbgI+DBlq+KgWgrPgNYx+Pc+XZGrfXf89iglQ0d2Lp/89YFtmR43edv+bp/hthcBH/a2bdc+FrD+/Bu7wt0v7u6JtMbQs8Xd1wBr+rsfM9vq7pVZ6FLeqZbCE0odENVy6uihYGoJ6bzkopZ0hlxagHFJy+Xxum7bmNkgYCTRi6MiIjJA0gn0BmCimU0ws8HAbcDGLm02Al+LHy8G/sVTjeWIiEhWpRxycfdTZrYSqCO6bXGtu+82sweAre6+EXga+B9mdgD4iCj0c6nfwzYFRLUUnlDqANVSqHJSS8oXRUVEpDjonaIiIoFQoIuIBKLgA93MhpjZG2a208x2m9lfx+snmNnrZnbAzNbFL9gWPDMrMbPtZvbjeLlY62gys11mtsPMtsbrPmVmPzWz/fH3C/Pdz3SY2Sgz22Bm/8/M9prZnGKsxcwmxecj8XXMzO4q0lr+PP59f8vMno9zoFh/V+6M69htZnfF63JyTgo+0IFPgM+5+xRgKrDAzGYTTS/wiLt/GvgN0fQDxeBOYG/ScrHWATDP3acm3U97L/Cyu08EXo6Xi8HfAT9x96uAKUTnp+hqcffG+HxMJZpX6bfAP1NktZjZZcB3gUp3n0x0M0ZiSpGi+l0xs8nAN4necT8F+LKZfZpcnRN3L5ovYCjwC2AW0busBsXr5wB1+e5fGv0vj0/e54AfA1aMdcR9bQIu6rKuEbg0fnwp0JjvfqZRx0jgl8Q3CBRzLV36fyPwajHWAlwGHAQ+RXQn3o+B+cX4uwL8KfB00vJ/BP4yV+ekGK7QE8MUO4BDwE+Bt4Ej7n4qbtJM9I+g0D1KdDITU16MpjjrAHBgs5lti6d0ABjj7u/Hj38NjMlP1/pkAvAB8N/jobCnzGwYxVlLstuA5+PHRVWLu7cADwPvAu8DR4FtFOfvylvAH5nZaDMbCnyR6E2YOTknRRHo7t7h0Z+R5UR/ulyV5y71mZl9GTjk7sUy10Qq17n7dOAm4Ntmdn3ykx5dehTDPbGDgOnAE+4+DfiYLn/+FlEtAMRjywuBF7o+Vwy1xOPJNxP9ZzsWGAYsyGunMuTue4mGijYDPwF2AB1d2mTtnBRFoCe4+xFgC9GfW6PiaQag++kICs1cYKGZNRHNWPk5orHbYqsD6LyKwt0PEY3TzgT+1cwuBYi/H8pfD9PWDDS7++vx8gaigC/GWhJuAn7h7v8aLxdbLV8AfunuH7h7O/Ai0e9Psf6uPO3uFe5+PdHY/z5ydE4KPtDN7GIzGxU/LiOal30vUbAvjpt9DfhRfnqYHndf7e7l7j6e6M/hf3H32ymyOgDMbJiZjUg8JhqvfYszp4Aoilrc/dfAQTObFK/6PLCHIqwlyVJ+N9wCxVfLu8BsMxsaT8OdOCdF97sCYGaXxN8vB/4t8Bw5OicF/05RM7uWaK71EqL/gNa7+wNm9vtEV7qfArYDf+bun+Svp+kzsyrgHnf/cjHWEff5n+PFQcBz7v6QmY0G1gOXA78CbnX3j3rYTcEws6nAU8Bg4B3g68T/1ii+WoYRBeLvu/vReF3RnZf49uQlwCmi34vlRGPmRfW7AmBmPyd6vawduNvdX87VOSn4QBcRkfQU/JCLiIikR4EuIhIIBbqISCAU6CIigVCgi4gEYkA/JFokXfFtXS/Hi79H9O66D+Llme5+Mi8d60Z8G+pJd/8/+e6LnNsU6FKQ3P0w0eyamNn9QKu7P5yv/pjZoKR5RLqqAlqBtAM9xf5EMqIhFykaZlZhZv8rnhCsLumt0/Vm9oiZbY3nM59hZi/Gc00/GLcZH893/j/jNhviyZJS7fdRi+Z7v9PM/jiej3u7mf3MzMaY2Xjg3wF/Hs9B/kdm9oyZLU7qd2v8vcrMfm5mG4E98aRz1WbWYGZvmtm3BvLnKeFRoEuxMOAxYLG7VwBrgYeSnj/p0bzs/0j0NupvA5OBO+LhG4BJwD+4+x8Ax4D/YGalKfY72N0r3f2/AP8bmB1P4lUD/KW7N8XHfMSjuch/nqKO6cCd7n4l0XzeR919BjAD+KaZTej7j0YkoiEXKRbnEwX0T6PpPSghmlo1YWP8fRewOzE1qZm9QzRd6RHgoLu/Grd7luhDFH6SYr/rkh6XA+viK/jBRPOo99Ub7p7Y7kbg2qSr+ZHAxAz3K6JAl6JhREE9p4fnE3N6nE56nFhO/DvvOs+Fp7Hfj5MePwb8rbtvjF8Ivb+HbU4R//VrZucRhX93+zPgO+5e18N+RPpEQy5SLD4BLjazOQBmVmpmf9jHfVye2B74CtEQSmMf9juS303Z+rWk9ceBEUnLTUQfAQfRvOSlPeyvDvj38bAPZnZlPLmWSEYU6FIsThNNnfoDM9tJ9EEBn+njPhqJPoxjL3Ah0YdanOzDfu8HXjCzbUQfh5bwEnBL4kVR4L8Bn433N4czr8qTPUU0LewvzOwt4En0V7P0g2ZblHNCfDfKjz360GGRIOkKXUQkELpCFxEJhK7QRUQCoUAXEQmEAl1EJBAKdBGRQCjQRUQC8f8Bc2QgqTmZRacAAAAASUVORK5CYII=\n", + "text/plain": [ + "
" + ] + }, + "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", @@ -580,6 +749,86 @@ "plt.grid(True)" ] }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "The predictiction is not correct. Let's see what went wrong" + ] + }, + { + "cell_type": "code", + "execution_count": 8, + "metadata": {}, + "outputs": [ + { + "data": { + "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXQAAAEKCAYAAAACS67iAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4yLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+j8jraAAAapUlEQVR4nO3df3xV9Z3n8dfHABJ+W9QMGDTMFnEckB/hZ3FsaKtg20HcYUTq0OpDSne3tDqOzEMe2xXr6OMx3bijs661ssraWVcCskyKfTANrUPWrh0VKCC/NoA2lUQ7KDZAbIAQPvvHOTe9xCT3JrmX3Pvl/Xw88uCec773nO8nh7xz8r3nfq+5OyIikv8u6u0OiIhIZijQRUQCoUAXEQmEAl1EJBAKdBGRQCjQRUQCkTLQzWy1mR0xsz0dbL/DzN4ys91m9gszm5D5boqISCrpXKE/D8ztZPuvgM+6+3jgb4BVGeiXiIh0UZ9UDdz9VTMr6WT7L5IWXweKe94tERHpqpSB3kV3A//U0UYzWwosBSgsLCwdNWpUtw5y9uxZLroojOF/1ZKbQqkllDpAtSQcOHDgQ3e/rN2N7p7yCygB9qRoMxvYDwxPZ5+lpaXeXVu2bOn2c3ONaslNodQSSh3uqiUB2OYd5GpGrtDN7DrgWeBmdz+aiX2KiEjX9PjvFzO7EtgALHb3Az3vkoiIdEfKK3QzWwOUAZeaWR2wEugL4O4/AB4EhgPfNzOAM+4+JVsdFhGR9qVzl8uiFNuXAEsy1iMRyRvNzc3U1dVx8uTJrB9r6NCh7N+/P+vHOR/SqaV///4UFxfTt2/ftPeb6btcROQCUldXx+DBgykpKSH+Cz1rTpw4weDBg7N6jPMlVS3uztGjR6mrq2P06NFp7zeMe4BEpFecPHmS4cOHZz3MLzRmxvDhw7v8l48CXUR6RGGeHd35virQRUQCoTF0EclrBQUFjB8/vnW5srKSkpKS3utQL1Kgi0heKywsZOfOne1uS7yDMpQpA1K5MKoUkQtGbW0tY8eO5atf/Srjxo3j8OHDlJeXM3XqVK677jpWrlzZ2vbRRx/l6quv5vrrr2fRokU89thjAJSVlbFt2zYAPvzww9Yr/paWFpYvX966r2eeeQaA6upqysrKWLBgAddccw133HFHYkoUtm7dymc+8xkmTJjAtGnTOHHiBHPnzj3nl9D111/Prl27ely7rtBFJCO++/Je9r13PKP7vHbkEFb+6R932qapqYmJEycCMHr0aB5//HEOHjzID3/4Q2bMmMHmzZs5ePAgb775Ju7OvHnzePXVVxk4cCAVFRXs3LmTM2fOMHnyZEpLSzs91nPPPcfQoUPZunUrp06dYtasWdx0000A7Nixg7179zJy5EhmzZrFa6+9xrRp01i4cCFr165l6tSpHD9+nMLCQhYvXszzzz/PE088wYEDBzh58iQTJvT8oyQU6CKS19oOudTW1nLVVVcxY8YMADZv3szmzZuZNGkSAI2NjRw8eJATJ05w6623MmDAAADmzZuX8libN2/mrbfeYv369QAcO3aMgwcP0q9fP6ZNm0ZxcTR7+MSJE6mtrWXo0KGMGDGCqVOnAjBkyBAAbr31VmbNmkV5eTmrV6/mzjvvzMj3QoEuIhmR6kr6fBo4cGDrY3dnxYoVfOMb3zinzRNPPNHh8/v06cPZs2cBzrkX3N158sknmTNnzjntq6urufjii1uXCwoKOHPmTIf7HzBgADfeeCM/+tGPWLduHdu3b0+vsBQ0hi4iQZszZw6rV6+msbERgPr6eo4cOcINN9xAZWUlTU1NnDhxgpdffrn1OSUlJa0hm7gaT+zr6aefprm5GYADBw7w8ccfd3jssWPH8v7777N161YgeodoIuiXLFnCt7/9baZOncoll1ySkVp1hS4iQbvpppvYv38/M2fOBGDQoEG88MILTJ48mYULFzJhwgQuv/zy1mERgPvvv5/bbruNVatW8aUvfal1/ZIlS6itrWXy5Mm4O5dddhmVlZUdHrtfv36sXbuWb33rWzQ1NVFYWMjPfvYzAEpLSxkyZAh33XVX5ortaKL0bH/pAy4iqiU3hVJLtuvYt29fVvef7Pjx41nd/8qVK728vDyrx0g4fvy419fX+5gxY7ylpaXDdu19f+nkAy405CIicp69+OKLTJ8+nUcffTSj98hryEVEBHjooYfO27G+8pWvfOJF2kzQFbqI9IjHb6CRzOrO91WBLiLd1r9/f44ePapQzzCP50Pv379/l56nIRcR6bbi4mLq6ur44IMPsn6skydPdjngclU6tSQ+sagrFOgi0m19+/bt0ifq9ER1dXXruz3zXbZq0ZCLiEggFOgiIoFQoIuIBEKBLiISCAW6iEggFOgiIoFQoIuIBEKBLiISiJSBbmarzeyIme3pYLuZ2X81s0Nm9paZTc58N0VEJJV0rtCfB+Z2sv1mYEz8tRR4uufdEhGRrkoZ6O7+KvBRJ01uAf4hnnv9dWCYmY3IVAdFRCQ9mZjL5QrgcNJyXbzu/Qzs+xO++/JefrGviadr/iUbuz/vGhpUSy4KpZZQ6oCwahly9hRlZZnf73mdnMvMlhINy1BUVER1dXWX91FXd4qWlhYaGhoy3LveoVpyUyi1hFIHhFVLYWFLt/IvpY4+my75CygB9nSw7RlgUdJyDTAi1T71maIR1ZKbQqkllDrcVUsCWf5M0Y3AV+O7XWYAx9w9K8MtIiLSsZRDLma2BigDLjWzOmAl0BfA3X8AbAK+CBwCfgfcla3OiohIx1IGursvSrHdgW9mrEciItIteqeoiEggFOgiIoFQoIuIBEKBLiISCAW6iEggFOgiIoFQoIuIBEKBLiISCAW6iEggFOgiIoFQoIuIBEKBLiISCAW6iEggFOgiIoFQoIuIBEKBLiISCAW6iEggFOgiIoFQoIuIBEKBLiISCAW6iEggFOgiIoFQoIuIBEKBLiISCAW6iEggFOgiIoFQoIuIBEKBLiISiLQC3czmmlmNmR0yswfa2X6lmW0xsx1m9paZfTHzXRURkc6kDHQzKwCeAm4GrgUWmdm1bZp9B1jn7pOA24HvZ7qjIiLSuXSu0KcBh9z9HXc/DVQAt7Rp48CQ+PFQ4L3MdVFERNJh7t55A7MFwFx3XxIvLwamu/uypDYjgM3AJcBA4Avuvr2dfS0FlgIUFRWVVlRUdKvTjY2NDBo0qFvPzTWqJTeFUksodYBqSZg9e/Z2d5/S7kZ37/QLWAA8m7S8GPhvbdrcB/xV/HgmsA+4qLP9lpaWendt2bKl28/NNaolN4VSSyh1uKuWBGCbd5Cr6Qy51AOjkpaL43XJ7gbWxb8g/gXoD1yaxr5FRCRD0gn0rcAYMxttZv2IXvTc2KbNu8DnAczsj4gC/YNMdlRERDqXMtDd/QywDKgC9hPdzbLXzB42s3lxs78Cvm5mu4A1wJ3xnwYiInKe9EmnkbtvAja1Wfdg0uN9wKzMdk1ERLpC7xQVEQmEAl1EJBAKdBGRQCjQRUQCoUAXEQmEAl1EJBAKdBGRQCjQRUQCoUAXEQmEAl1EJBAKdBGRQCjQRUQCoUAXEQmEAl1EJBAKdBGRQCjQRUQCoUAXEQmEAl1EJBAKdBGRQCjQRUQCoUAXEQmEAl1EJBAKdBGRQCjQRUQCoUAXEQmEAl1EJBAKdBGRQKQV6GY218xqzOyQmT3QQZvbzGyfme01sxcz200REUmlT6oGZlYAPAXcCNQBW81so7vvS2ozBlgBzHL335rZ5dnqsIiItC+dK/RpwCF3f8fdTwMVwC1t2nwdeMrdfwvg7kcy200REUnF3L3zBmYLgLnuviReXgxMd/dlSW0qgQPALKAAeMjdf9LOvpYCSwGKiopKKyoqutXpxsZGBg0a1K3n5hrVkptCqSWUOkC1JMyePXu7u09pb1vKIZc09QHGAGVAMfCqmY1394bkRu6+ClgFMGXKFC8rK+vWwaqrq+nuc3ONaslNodQSSh2gWtKRzpBLPTAqabk4XpesDtjo7s3u/iuiq/UxmemiiIikI51A3wqMMbPRZtYPuB3Y2KZNJdHVOWZ2KXA18E4G+ykiIimkDHR3PwMsA6qA/cA6d99rZg+b2by4WRVw1Mz2AVuA5e5+NFudFhGRT0prDN3dNwGb2qx7MOmxA/fFXyIi0gv0TlERkUAo0EVEAqFAFxEJhAJdRCQQCnQRkUAo0EVEAqFAFxEJhAJdRCQQCnQRkUAo0EVEAqFAFxEJhAJdRCQQCnQRkUAo0EVEAqFAFxEJhAJdRCQQCnQRkUAo0EVEAqFAFxEJhAJdRCQQCnQRkUAo0EVEAqFAFxEJhAJdRCQQCnQRkUAo0EVEAqFAFxEJRFqBbmZzzazGzA6Z2QOdtPszM3Mzm5K5LoqISDpSBrqZFQBPATcD1wKLzOzadtoNBu4B3sh0J0VEJLV0rtCnAYfc/R13Pw1UALe00+5vgO8BJzPYPxERSZO5e+cNzBYAc919Sby8GJju7suS2kwG/qO7/5mZVQP3u/u2dva1FFgKUFRUVFpRUdGtTjc2NjJo0KBuPTfXqJbcFEotodQBqiVh9uzZ29293WHtPj3qFWBmFwF/B9yZqq27rwJWAUyZMsXLysq6dczq6mq6+9xco1pyUyi1hFIHqJZ0pDPkUg+MSloujtclDAbGAdVmVgvMADbqhVERkfMrnUDfCowxs9Fm1g+4HdiY2Ojux9z9UncvcfcS4HVgXntDLiIikj0pA93dzwDLgCpgP7DO3fea2cNmNi/bHRQRkfSkNYbu7puATW3WPdhB27Ked0tERLpK7xQVEQmEAl1EJBAKdBGRQPT4PnSRnqjcUU95VQ3vNTQxclghy+eMZf6kK3q7W5Imnb/cokCXXlO5o54VG3bT1NwCQH1DEys27AZQKOQBnb/coyEX6TXlVTWtYZDQ1NxCeVVNL/VIukLnL/co0KXXvNfQ1KX1klt0/nKPAl16zchhhV1aL7lF5y/3KNCl1yyfM5bCvgXnrCvsW8DyOWN7qUfSFTp/uUcvikqvSbxwprsk8pPOX+5RoEuvmj/pCgVAHtP5yy0achERCYQCXUQkEAp0EZFAKNBFRAKhQBcRCYQCXUQkEAp0EZFAKNBFRAKhQBcRCYQCXUQkEAp0EZFAKNBFRAKhQBcRCYQCXUQkEAp0EZFAKNBFRAKRVqCb2VwzqzGzQ2b2QDvb7zOzfWb2lpm9YmZXZb6rIiLSmZSBbmYFwFPAzcC1wCIzu7ZNsx3AFHe/DlgP/OdMd1RERDqXzhX6NOCQu7/j7qeBCuCW5AbuvsXdfxcvvg4UZ7abIiKSSjqBfgVwOGm5Ll7XkbuBf+pJp0REpOvM3TtvYLYAmOvuS+LlxcB0d1/WTtu/AJYBn3X3U+1sXwosBSgqKiqtqKjoVqcbGxsZNGhQt56ba1RLbgqlllDqANWSMHv27O3uPqXdje7e6RcwE6hKWl4BrGin3ReA/cDlqfbp7pSWlnp3bdmypdvPzTWqJTeFUksodbirlgRgm3eQq+kMuWwFxpjZaDPrB9wObExuYGaTgGeAee5+pDu/dUREpGf6pGrg7mfMbBlQBRQAq919r5k9TPSbYiNQDgwCXjIzgHfdfV4W+y3Srsod9ZRX1fBeQxMjhxWyfM5YgE+smz+ps5eBMnfsbBwnHd+p3M2aNw5z77hm7l6xiUXTR/HI/PG90hc5f1IGOoC7bwI2tVn3YNLjL2S4XyJdVrmjnhUbdtPU3AJAfUMTy1/aBQbNLd66bsWG3QAZDdv2jp2N46TjO5W7eeH1d1uXW9xblxXqYdM7RSUY5VU1rYGa0HzWW8M8oam5hfKqmqwfOxvHSceaNw53ab2EQ4EuwXivoSkrbXuyv0wfJx0tHdy51tF6CYcCXYIxclhhVtr2ZH+ZPk46CqLXsdJeL+FQoEswls8ZS2HfgnPW9b3I6FtwbpAV9i1ofbE0m8fOxnHSsWj6qC6tl3Ck9aKoSD5IvPjYG3e5dHTs3rjLJfHCZ2LMvMBMd7lcIBToEpT5k65oN0TPR7B2dOze8Mj88TwyfzzV1dW8fUdZb3dHzhMNuYiIBEKBLiISCAW6iEggFOgiIoFQoIuIBEKBLiISCAW6iEggFOgiIoFQoIuIBEKBLiISCAW6iEggFOgiIoFQoIuIBEKBLiISCAW6iEggFOgiIoFQoIuIBEKBLiISCAW6iEggFOgiIoFQoIuIBEKBLiISiD7pNDKzucDfAwXAs+7+t222Xwz8A1AKHAUWunttZrsqEq7KHfWUV9XwXkMTI4cVsnzOWF7a9i6vvf1Ra5tZ/+ZT/PmUKz/RDvjEum2//og1bxzm3nHN3L1iE4umj+KR+ePTOu78SVd0uD6d5yeO3eJOgVmXjt1eLeket712F5qUgW5mBcBTwI1AHbDVzDa6+76kZncDv3X3T5vZ7cD3gIXZ6LBIaCp31LNiw26amlsAqG9o4t61Oz/R7rW3Pzon4Osbmli+fhc4NJ/11nX3rd3J2aTntbjzwuvvApwTrO0dd8WG3Wz79Uf87+31n1gPnBOa7T2/J8de/tIuMGhu+X0t6R63vXYXonSGXKYBh9z9HXc/DVQAt7Rpcwvww/jxeuDzZmaZ66ZIuMqralrDqauaW7w1zBPOdtB2zRuHUx63qbmFNW8cbnd9eVVNyuf35NjNZ701zLt63PbaXYjM3TtvYLYAmOvuS+LlxcB0d1+W1GZP3KYuXn47bvNhm30tBZbGi2OB7p6BS4EPU7bKD6olN523Wvr9wadLs7Xvlt8do2DA0Nbl0785tL0nx+3J8zPw3EuBDzt7bvIxclxP/n9d5e6XtbchrTH0THH3VcCqnu7HzLa5+5QMdKnXqZbcFEotZrbtzLEjeV8HhHNOIHu1pDPkUg+MSloujte128bM+gBDiV4cFRGR8ySdQN8KjDGz0WbWD7gd2NimzUbga/HjBcA/e6qxHBERyaiUQy7ufsbMlgFVRLctrnb3vWb2MLDN3TcCzwH/08wOAR8RhX429XjYJoeoltwUSi2h1AGqJaWUL4qKiEh+0DtFRUQCoUAXEQlEzge6mfU3szfNbJeZ7TWz78brR5vZG2Z2yMzWxi/Y5jwzKzCzHWb243g5X+uoNbPdZrbTzLbF6z5lZj81s4Pxv5f0dj/TYWbDzGy9mf0/M9tvZjPzsRYzGxufj8TXcTO7N09r+cv4532Pma2JcyBff1buievYa2b3xuuyck5yPtCBU8Dn3H0CMBGYa2YziKYXeNzdPw38lmj6gXxwD7A/aTlf6wCY7e4Tk+6nfQB4xd3HAK/Ey/ng74GfuPs1wASi85N3tbh7TXw+JhLNq/Q74B/Js1rM7Arg28AUdx9HdDNGYkqRvPpZMbNxwNeJ3nE/AfiymX2abJ0Td8+bL2AA8EtgOtG7rPrE62cCVb3dvzT6XxyfvM8BPwYsH+uI+1oLXNpmXQ0wIn48Aqjp7X6mUcdQ4FfENwjkcy1t+n8T8Fo+1gJcARwGPkV0J96PgTn5+LMC/DnwXNLyfwL+OlvnJB+u0BPDFDuBI8BPgbeBBnc/EzepI/pPkOueIDqZiSkvhpOfdQA4sNnMtsdTOgAUufv78ePfAEW907UuGQ18APyPeCjsWTMbSH7Wkux2YE38OK9qcfd64DHgXeB94Biwnfz8WdkD/ImZDTezAcAXid6EmZVzkheB7u4tHv0ZWUz0p8s1vdylLjOzLwNH3D1f5ppI5Xp3nwzcDHzTzG5I3ujRpUc+3BPbB5gMPO3uk4CPafPnbx7VAkA8tjwPeKnttnyoJR5PvoXol+1IYCAwt1c71U3uvp9oqGgz8BNgJ9DSpk3GzkleBHqCuzcAW4j+3BoWTzMA7U9HkGtmAfPMrJZoxsrPEY3d5lsdQOtVFO5+hGicdhrwr2Y2AiD+90jv9TBtdUCdu78RL68nCvh8rCXhZuCX7v6v8XK+1fIF4Ffu/oG7NwMbiH5+8vVn5Tl3L3X3G4jG/g+QpXOS84FuZpeZ2bD4cSHRvOz7iYJ9Qdzsa8CPeqeH6XH3Fe5e7O4lRH8O/7O730Ge1QFgZgPNbHDiMdF47R7OnQIiL2px998Ah81sbLzq88A+8rCWJIv4/XAL5F8t7wIzzGxAPA134pzk3c8KgJldHv97JfBvgRfJ0jnJ+XeKmtl1RHOtFxD9Alrn7g+b2R8SXel+CtgB/IW7n+q9nqbPzMqA+939y/lYR9znf4wX+wAvuvujZjYcWAdcCfwauM3dP+pgNznDzCYCzwL9gHeAu4j/r5F/tQwkCsQ/dPdj8bq8Oy/x7ckLgTNEPxdLiMbM8+pnBcDMfk70elkzcJ+7v5Ktc5LzgS4iIunJ+SEXERFJjwJdRCQQCnQRkUAo0EVEAqFAFxEJxHn9kGiRdMW3db0SL/4B0bvrPoiXp7n76V7pWDvi21BPu/sversvcmFToEtOcvejRLNrYmYPAY3u/lhv9cfM+iTNI9JWGdAIpB3oKfYn0i0acpG8YWalZvZ/4gnBqpLeOl1tZo+b2bZ4PvOpZrYhnmv6kbhNSTzf+f+K26yPJ0tKtd8nLJrv/R4z+9N4Pu4dZvYzMysysxLg3wF/Gc9B/idm9ryZLUjqd2P8b5mZ/dzMNgL74knnys1sq5m9ZWbfOJ/fTwmPAl3yhQFPAgvcvRRYDTyatP20R/Oy/4DobdTfBMYBd8bDNwBjge+7+x8Bx4H/YGZ9U+y3n7tPcff/AvxfYEY8iVcF8NfuXhsf83GP5iL/eYo6JgP3uPvVRPN5H3P3qcBU4OtmNrrr3xqRiIZcJF9cTBTQP42m96CAaGrVhI3xv7uBvYmpSc3sHaLpShuAw+7+WtzuBaIPUfhJiv2uTXpcDKyNr+D7Ec2j3lVvunvieTcB1yVdzQ8FxnRzvyIKdMkbRhTUMzvYnpjT42zS48Ry4v9523kuPI39fpz0+Eng79x9Y/xC6EMdPOcM8V+/ZnYRUfi3tz8DvuXuVR3sR6RLNOQi+eIUcJmZzQQws75m9sdd3MeViecDXyEaQqnpwn6H8vspW7+WtP4EMDhpuZboI+Agmpe8bwf7qwL+fTzsg5ldHU+uJdItCnTJF2eJpk79npntIvqggM90cR81RB/GsR+4hOhDLU53Yb8PAS+Z2Xaij0NLeBm4NfGiKPDfgc/G+5vJuVflyZ4lmhb2l2a2B3gG/dUsPaDZFuWCEN+N8mOPPnRYJEi6QhcRCYSu0EVEAqErdBGRQCjQRUQCoUAXEQmEAl1EJBAKdBGRQPx/tRIgu/jEeqIAAAAASUVORK5CYII=\n", + "text/plain": [ + "
" + ] + }, + "metadata": { + "needs_background": "light" + }, + "output_type": "display_data" + } + ], + "source": [ + "data_pred.plot(x=\"Temperature\",y=\"Frequency\",kind=\"line\",ylim=[0,1.2])\n", + "plt.scatter(x=data[\"Temperature\"],y=data[\"Frequency\"])\n", + "plt.grid(True)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "There were warnings during the construction of the log model. let's try and change it" + ] + }, + { + "cell_type": "code", + "execution_count": 12, + "metadata": {}, + "outputs": [ + { + "ename": "DistributionNotFound", + "evalue": "The 'statsmodel==0.9.0' distribution was not found and is required by the application", + "output_type": "error", + "traceback": [ + "\u001b[0;31m---------------------------------------------------------------------------\u001b[0m", + "\u001b[0;31mDistributionNotFound\u001b[0m Traceback (most recent call last)", + "\u001b[0;32m\u001b[0m in \u001b[0;36m\u001b[0;34m\u001b[0m\n\u001b[1;32m 1\u001b[0m \u001b[0;32mimport\u001b[0m \u001b[0mpkg_resources\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m----> 2\u001b[0;31m \u001b[0mpkg_resources\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mrequire\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m\"statsmodel==0.9.0\"\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m 3\u001b[0m \u001b[0;32mimport\u001b[0m \u001b[0mstatsmodel\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n", + "\u001b[0;32m~/miniconda3/lib/python3.7/site-packages/pkg_resources/__init__.py\u001b[0m in \u001b[0;36mrequire\u001b[0;34m(self, *requirements)\u001b[0m\n\u001b[1;32m 899\u001b[0m \u001b[0mincluded\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0meven\u001b[0m \u001b[0;32mif\u001b[0m \u001b[0mthey\u001b[0m \u001b[0mwere\u001b[0m \u001b[0malready\u001b[0m \u001b[0mactivated\u001b[0m \u001b[0;32min\u001b[0m \u001b[0mthis\u001b[0m \u001b[0mworking\u001b[0m \u001b[0mset\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 900\u001b[0m \"\"\"\n\u001b[0;32m--> 901\u001b[0;31m \u001b[0mneeded\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mresolve\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mparse_requirements\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mrequirements\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m 902\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 903\u001b[0m \u001b[0;32mfor\u001b[0m \u001b[0mdist\u001b[0m \u001b[0;32min\u001b[0m \u001b[0mneeded\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n", + "\u001b[0;32m~/miniconda3/lib/python3.7/site-packages/pkg_resources/__init__.py\u001b[0m in \u001b[0;36mresolve\u001b[0;34m(self, requirements, env, installer, replace_conflicting, extras)\u001b[0m\n\u001b[1;32m 785\u001b[0m \u001b[0;32mif\u001b[0m \u001b[0mdist\u001b[0m \u001b[0;32mis\u001b[0m \u001b[0;32mNone\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 786\u001b[0m \u001b[0mrequirers\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mrequired_by\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mget\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mreq\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0;32mNone\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m--> 787\u001b[0;31m \u001b[0;32mraise\u001b[0m \u001b[0mDistributionNotFound\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mreq\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mrequirers\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m 788\u001b[0m \u001b[0mto_activate\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mappend\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mdist\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 789\u001b[0m \u001b[0;32mif\u001b[0m \u001b[0mdist\u001b[0m \u001b[0;32mnot\u001b[0m \u001b[0;32min\u001b[0m \u001b[0mreq\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n", + "\u001b[0;31mDistributionNotFound\u001b[0m: The 'statsmodel==0.9.0' distribution was not found and is required by the application" + ] + } + ], + "source": [ + "import pkg_resources\n", + "pkg_resources.require(\"statsmodel==0.9.0\")\n", + "import statsmodel" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "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']], \n", + " family=sm.families.Binomial(sm.families.links.logit)).fit()\n", + "\n", + "logmodel.summary()\n" + ] + }, { "cell_type": "markdown", "metadata": { @@ -607,9 +856,20 @@ }, { "cell_type": "code", - "execution_count": null, + "execution_count": 9, "metadata": {}, - "outputs": [], + "outputs": [ + { + "data": { + "image/png": "iVBORw0KGgoAAAANSUhEUgAAAZAAAAEQCAYAAACeDyIUAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4yLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+j8jraAAAgAElEQVR4nO3deXwUZZ4G8Keq+sp9n9ygQJRLATOjKGNAkoGAqIs4uLrrgbOKos44HxEdOVQcdFdFvEadGceF1V3GC6ICooMCIociAcKhIZAAnatzd/qsevePJm3CIZ0i6aQ7z/fzQehOVfXvtTt5Um+99b6SEEKAiIioneSuLoCIiEITA4SIiHRhgBARkS4MECIi0oUBQkREujBAiIhIl6AEyNKlS5GTk4MhQ4bg0KFDZ9xGVVUsWrQIEydOxDXXXINVq1YFozQiItIpKAEyYcIErFy5Er169TrrNmvWrEFpaSnWr1+P//3f/8Xy5ctx7NixYJRHREQ6BCVAxowZg4yMjJ/d5pNPPsGMGTMgyzISExMxceJErF27NhjlERGRDt3mGojVakVmZqb/cUZGBsrLy7uwIiIi+jndJkCIiCi0GLq6gBYZGRk4ceIERowYAeD0M5JA1dbaoWnhOb1XUlI0bLamri6j04Rz+8K5bQDbF8pkWUJCQpSufbtNgOTl5WHVqlWYNGkS6urqsGHDBqxcubLdx9E0EbYBAiCs2waEd/vCuW0A29cTBaUL68knn8RVV12F8vJy3HbbbZgyZQoAYPbs2dizZw8A4Nprr0Xv3r0xadIk3HjjjZgzZw769OkTjPKIiEgHKdymc7fZmsL2N4WUlBhUVTV2dRmdJpzbF85tA9i+UCbLEpKSovXt28G1EBFRD8EAISIiXRggRESkCwOEiIh0YYAQEZEuDBAiItKFAUJERLowQIiISBcGCBER6cIAISIiXRggRESkCwOEiIh0YYAQEZEuDBAiItKFAUJERLowQIiISBcGCBER6cIAISIiXRggRESkCwOEiIh0YYAQEZEuDBAiItKFAUJERLowQIiISBcGCBER6cIAISIiXRggRESkCwOEiIh0YYAQEZEuDBAiItKFAUJERLowQIiISBcGCBER6cIAISIiXRggRESkiyFYL1RSUoJ58+ahrq4O8fHxWLp0Kfr3799mG5vNhkceeQRWqxVerxfZ2dl47LHHYDAErUwiIgpQ0M5AFixYgFmzZmHdunWYNWsWHn/88dO2ee211zBo0CCsWbMGq1evxr59+7B+/fpglUhERO0QlACx2WwoKipCfn4+ACA/Px9FRUWoqalps50kSbDb7dA0DW63Gx6PB2lpacEokYiI2ikofUNWqxVpaWlQFAUAoCgKUlNTYbVakZiY6N/unnvuwX333Ydx48bB4XDg5ptvxujRo9v1WklJ0R1ae3eTkhLT1SV0qnBuXzi3DWD7eqJudXFh7dq1GDJkCP7+97/Dbrdj9uzZWLt2LfLy8gI+hs3WBE0TnVhl10lJiUFVVWNXl9Fpwrl94dw2gO0LZbIs6f7FOyhdWBkZGaioqICqqgAAVVVRWVmJjIyMNtutWLEC06ZNgyzLiImJQU5ODrZt2xaMEomIqJ2CEiBJSUnIyspCQUEBAKCgoABZWVltuq8AoHfv3vjqq68AAG63G1u3bsWFF14YjBKJiKidgjYKa+HChVixYgVyc3OxYsUKLFq0CAAwe/Zs7NmzBwAwf/58fPvtt5g6dSqmT5+O/v3748YbbwxWiURE1A6SECKsLhjwGkjoCuf2hXPbALYvlHX7ayBERBR+GCBERKQLA4SIiHRhgBARkS4MECIi0oUBQkREujBAiIhIFwYIERHpwgAhIiJdGCBERKQLA4SIiHRhgBARkS4MECIi0oUBQkREujBAiIhIFwYIERHpwgAhIiJdGCBERKQLA4SIiHRhgBARkS4MECIi0oUBQkREujBAiIhIFwYIERHpwgAhIiJdGCBERKQLA4SIiHRhgBARkS4MECIi0oUBQkREujBAiIhIFwYIERHpEnCAbNiwAV6vtzNrISKiEBJwgLz44osYN24cFi9ejN27d3dmTUREFAICDpDVq1fjrbfegtlsxn333Yfc3Fy88sorOHbsWED7l5SUYObMmcjNzcXMmTNx5MiRM273ySefYOrUqcjPz8fUqVNRXV0daIlERBREkhBCtHcnIQS2bt2KP/3pT/jhhx9w6aWXYubMmcjPz4csnzmTbr31Vtxwww249tpr8dFHH+G9997D22+/3WabPXv24OGHH8bf//53pKSkoLGxESaTCWazOeDabLYmaFq7mxQSUlJiUFXV2NVldJpwbl84tw1g+0KZLEtISorWt297dygtLcXLL7+MhQsXwuVyYe7cuZgxYwZWrlyJuXPnnnEfm82GoqIi5OfnAwDy8/NRVFSEmpqaNtu99dZbuP3225GSkgIAiImJaVd4EBFR8BgC3XDlypX46KOPcPToUfz617/GM888g1GjRvm/npubi8svv/yM+1qtVqSlpUFRFACAoihITU2F1WpFYmKif7vi4mL07t0bN998M5qbm3HNNdfg7rvvhiRJettHRESdJOAA+eqrr3DbbbdhwoQJMJlMp309IiICy5cvP69iVFXFwYMH8be//Q1utxt33nknMjMzMX369ICPofdULFSkpMR0dQmdKpzbF85tA9i+nijgAHnxxRchyzKMRqP/OY/HAyGEP1DGjRt3xn0zMjJQUVEBVVWhKApUVUVlZSUyMjLabJeZmYm8vDyYTCaYTCZMmDABhYWF7QoQXgMJXeHcvnBuG8D2hbKgXAO5/fbbsW/fvjbP7du3D3fcccc5901KSkJWVhYKCgoAAAUFBcjKymrTfQX4ro1s3rwZQgh4PB588803GDp0aKAlEhFREAUcIAcPHsTIkSPbPDdixAgcOHAgoP0XLlyIFStWIDc3FytWrMCiRYsAALNnz8aePXsAAFOmTEFSUhImT56M6dOn44ILLsC//Mu/BFoiEREFUcBdWLGxsaiurvaPkAKA6upqREREBLT/oEGDsGrVqtOef+ONN/z/lmUZjzzyCB555JFAyyIioi4S8BnIpEmT8Pvf/x6HDh2Cw+HAwYMH8fDDD+PXv/51Z9ZHRETdVMAB8uCDD2LQoEGYMWOG/8bBAQMG4He/+11n1kdERN1Uu+9EF0KgtrYWCQkJ3fL+DI7CCl3h3L5wbhvA9oWy8xmFFfA1EABobGxESUkJ7HZ7m+d/+ctf6npxIiIKXQEHyPvvv4/FixcjMjISFovF/7wkSfj88887pTgiIuq+Ag6Q559/HsuWLcP48eM7sx4iIgoRAV9EV1X1rHeaExFRzxNwgMyePRuvvvoqNE3rzHqIiChEBNyF9dZbb6G6uhpvvvkm4uPj23xt48aNHV0XERF1cwEHyLPPPtuZdRARUYgJOEAuu+yyzqyDiIhCTMDXQNxuN55//nlMmDABo0ePBgBs3rwZK1as6LTiiIio+wo4QJYsWYJDhw7hP//zP/13oF944YV45513Oq04IiLqvgLuwtqwYQPWr1+PyMhIyLIvd9LS0lBRUdFpxRERUfcV8BmI0WiEqqptnqupqTltRBYREfUMAQdIXl4eHn74YZSVlQEAKisrsXjxYkyZMqXTiiMiou6rXdO59+7dG9OmTUNDQwNyc3ORmpqKOXPmdGZ9RETUTbV7OnfA13XF6dyDL5ynlAbCu33h3DaA7QtlQZnOvaXrqkXrKd379Omj68WJiCh0BRwg11xzDSRJQusTlpYzkP3793d8ZURE1K0FHCAHDhxo87iqqgovvfQSxowZ0+FFERFR9xfwRfRTpaSk4NFHH8Vzzz3XkfUQEVGI0B0gAHD48GE4HI6OqoWIiEJIwF1Ys2bNajPqyuFw4Mcff+QwXiKiHirgAJkxY0abxxERERg6dCj69+/f0TUREVEICDhArrvuus6sg4iIQkzAAbJs2bKAtrv//vt1F0NERKEj4AA5evQo1q9fj2HDhqFXr144ceIE9uzZg0mTJsFsNndmjURE1A0FHCBCCPzXf/0XcnNz/c+tX78ea9euxdNPP90pxRERUfcV8DDer776ChMnTmzzXE5ODr788ssOL4qIiLq/gAOkX79+WLlyZZvn3nnnHfTt27fDiyIiou4v4C6sJ598Evfeey/efPNN/0qEBoMBy5cv78z6iIiomwo4QC666CKsW7cOu3fvRmVlJVJSUjBq1CgYjcbOrI+IiLop3VOZjB07Fh6PB83NzR1ZDxERhYiAz0AOHjyIu+++GyaTCRUVFZg8eTJ27NiBDz74AC+88EJn1khERN1QwGcgCxcuxNy5c7F27VoYDL7cGTt2LL799tuA9i8pKcHMmTORm5uLmTNn4siRI2fd9vDhwxg5ciSWLl0aaHlERBRkAQfIjz/+iGuvvRbATwtJRUZGwuVyBbT/ggULMGvWLKxbtw6zZs3C448/fsbtVFXFggULThsyTHQuQgBa+1doJiKdAg6QXr16Ye/evW2eKywsDGgYr81mQ1FREfLz8wEA+fn5KCoqQk1NzWnbvv766/jVr37FSRpJF5dH7eoSiHqMgK+B3H///fjtb3+Lm266CR6PB3/+85/x7rvv4oknnjjnvlarFWlpaVAUBQCgKApSU1NhtVqRmJjo3+7AgQPYvHkz3n77bbzyyis6mgPdi8OHipSUmK4uoVOdT/tUTcBW50ByQkSbpQe6C753oS3c26dHwAFy9dVX480338T//d//YezYsTh+/DiWL1+OYcOGdUghHo8Hf/zjH/H000/7g0YPm60Jmhae3RgpKTGoqmrs6jI6zfm2TwjA1uCE1+2BUTmvtdI6HN+70BbO7ZNlSfcv3gEFiKqqyM3NxSeffIKFCxe2+0UyMjJQUVEBVVWhKApUVUVlZSUyMjL821RVVaG0tBR33XUXAKChoQFCCDQ1NQV0lkME+K6BOFxeGCNNXV0KUdgLKEAURYGiKHC5XDCZ2v+NmZSUhKysLBQUFODaa69FQUEBsrKy2nRfZWZmYtu2bf7Hy5cvR3NzMx5++OF2vx71bE6XiiiLBkXuXmchROEm4O+wW2+9FQ888AC2b9+O0tJSlJWV+f8EYuHChVixYgVyc3OxYsUKLFq0CAAwe/Zs7NmzR1/1RGegCQG709vVZRCFPUmInx/3WFVVhZSUFAwdOtS3gySh9S6SJGH//v2dW2U78BpI6OqIayBV9Q5omoAsSUiMs8Agd4+L6XzvQls4t+98roGc8wykZf2PAwcO4MCBA8jJyfH/+8CBA90qPIhaaEKg0e5GeP4qQdQ9nDNATj1B2bFjR6cVQ9SRXB4VDhe7sog6yzkD5NTx9Ofo8SLqVpqaPfCGaZcmUVc75ygsVVXxzTff+IPj1McA8Mtf/rLzKiQ6Dy1dWfExZnSPqyFE4eOcAZKUlIT58+f7H8fHx7d5LEkSPv/8886pjqgDuDwqml1eRJkDvm+WiAJwzu+oL774Ihh1EHUqe7MHZqPSbUZlEYUD/kpGIa+wuBprt5XC6VYRZTFgzNBUDOmb0GYbTQg02N1IjDEB7MwKWS3vdXW9E8lxFuRl98WIQcldXVaPxVt1KaQVFldj5WeHUGd3I8KioMHhweotJThYWnvath6PCrdH64IqqSO0fq8jLQbU2d1Y+dkhFBZXd3VpPRYDhELa2m2lUBQZZqMCSZJgMihQFBmbdp84bVsBoNHhCX6R1CFOfa/NRt97vXZbaVeX1mMxQCikVdc7YTK0/RgbFRm1jWde6Mzj1bhmSIg603ttMsiornd2UUXEAKGQlhxngdvbtlvKo2pIiDGfdR+704tuuFwIncOZ3mu3V0NynKWLKiIGCIW0vOy+UFXfWYUQAm6vClXVcOXIzLPu4/GqcPFaSMg59b12eXzvdV72uVdFpc7BUVgU0lpG4KzdVgqHU0VshBE5l/Q6bRRWa0IATrcXJgPXDAklrd9rjsLqHhggFPJGDErGiEHJbWbjPReHy4tIsxEGhX1ZoaTlvabugV1Y1CMJATQ2u7u6DKKQxgChHsvlUdHs5my9RHoxQKhHa2r2QNV4QZ1IDwYI9WiaJtBg93DhKSIdGCDU43HhKSJ9GCBEAJoc7Moiai8GCBF8XVmN7MoiahcGCNFJzpMLTxFRYBggRK00NbvhUdmVRRQIBgiFBVu9Ext2lqHpPKdrFwJosLsh2JlFdE6cyoTCwuotJdhUaEVCjBn/ljcUqQkRuo/l8WpocngRG2mEYI4QnRXPQCgsXDE8Ayajbx2Q1z7aix+P1Z/X8ZqdHjjdXDeE6OcwQCgsDO4Tj4dnXYqYSCOcbhVvfbof24oqdB9PCKDR7oYawMSMRD0VA4TCRv/0WMy5fjgykiKhCeCjzSVY8/UR3SHg1YTveggzhOiMGCAUVuKjzbhr2sXI6udbD2Tr3nL8/dMDuu80d3lUNDrcXMGQ6AwYIBR2zEYFN18zGFedXJXwx+P1eOWDvaioadZ1PIfTy/tDiM6AAUJhSZYl5GX3xY1XXwCDIsHW4MSrH+3F3pKadh9LAGhs9sDL6yFEbTBAKKyNujAZd027GHFRJrg9Gv7ns0NYt700oFULW9M0gfomFzReECHyY4BQ2OudEo051w/HgIxYAMCX35/A3z7d3+6bDj1eDQ3Nbt5iSHQSA4R6hOgII26fkoVxwzMAAMXHG/DS+3twtLyxXcdxutSTd6oTUdACpKSkBDNnzkRubi5mzpyJI0eOnLbNyy+/jClTpmDq1Km4/vrrsWnTpmCVRz2AIkuY/Mt++M3EC2Eyymiwu/HGmiJ8tftEu7qmHC4vh/cSIYgBsmDBAsyaNQvr1q3DrFmz8Pjjj5+2zYgRI/CPf/wDa9aswZIlS/Dggw/C6XQGq0TqIYYPTMKc64YjLSECmhBYu60U/732YLu6tBwuL+qbXZ1YJVH3F5QAsdlsKCoqQn5+PgAgPz8fRUVFqKlpOyLmyiuvRESEbw6jIUOGQAiBurq6YJRIPUxKfATuvm4YxgxJAQAcLKvDS+8VovhE4FOgOF0q6u3uziqRqNsLSoBYrVakpaVBURQAgKIoSE1NhdVqPes+H374Ifr27Yv09PRglEg9kMmg4Prxg3BjzgW+Lq1mD/5asB/rtpfCG+CU7g6XF40OD8CrItQDdcvZeLdv345ly5bhr3/9a7v3TUqK7oSKuo+UlJiuLqFTnU/7VE3AK0vtvjaRc1kUhl2Ygr+s3oej1gZ8+f0JlFgbcfu0i5GeFBXQMQxmA+JjzJB+5pZ1vnehLdzbp0dQAiQjIwMVFRVQVRWKokBVVVRWViIjI+O0bXft2oU//OEPeOWVVzBw4MB2v5bN1tTuMf6hIiUlBlVV7Rs1FErOt31CALX1Dl3vvwHAnVOGYsPOY/jq+xMorWjEk3/dhrzsvvjFxemQzzGXSQ0AW40BsVFGSDh9W753oS2c2yfLku5fvIPShZWUlISsrCwUFBQAAAoKCpCVlYXExMQ22xUWFuLBBx/Eiy++iIsvvjgYpRH5KbKM3Mv64s6pFyE+2gSvKlDw9VH89eP9qG0892AOh8uL2kY33F4N7NKinkASIjiDEYuLizFv3jw0NDQgNjYWS5cuxcCBAzF79mzMnTsXw4cPxw033IDjx48jLS3Nv98zzzyDIUOGBPw6PAMJXR1xBlKl8wzkVE63Fx9/fRTfHqoCAJiMMvKy++KyrLRzno1IEmA0KIi2GGE2yRCC712oC+f2nc8ZSNACJFgYIKGrOwVIiw07y/Dl9yf8U8KnJUTg8mHp2P1jNWobXUiIMePKkZkY0jfhtH0lCTAbFERGGNErI0532wqLq7F2Wymq651IjrMgL7svRgxKPq92dZTVmw9j/Y5jcHpUWIwKJo3tjWnj2t/13N2F8/det+/CIgpFB0trseuHKsRFmxBh9o0grKh14INNJSivdcBsUtDg8GD1lhIcLK09bX8hAKdHRW2jE7Z6p655tAqLq7Hys0Oos7sRaTGgzu7Gys8OobC4+rzbd75Wbz6M1V8fgcujwiD7pr5f/fURrN58uKtLoyBhgBCdxabdJ6AoMiwmAxJiLEiMNfu/1uz0orreCQhAUWRs2n3irMcRwtclVlPvhMPtbdfVkbXbSqEoMsxGBZIkwWxUoCgy1m4rPY+WdYz1O45BggRFliBJsu9vSFi/41hXl0ZBwgAhOovaRheMyk/fIhaToc34Kq8qUF3vhN3hga3h3BfZvZpAfZMbNQ0OOD2BrbdeXe+EydD229RkkH3h1cWcbi/kUy4HyZLveeoZGCBEZ5EQY4bnlBsKDYoEoyIhOc7iDxeHS0VjswffFJUHdP3F4xWoa3ShptEFt0f92TOS5DjLyVFdP3F7NSTHWdrdno5mMRlwanM14XueegYGCNFZXDkyE6qqwe1VIYSA26vCZFRgMhkACUiKMyMqwndWoglg9eYjePmDPTh8oiGg47s9KmoaXbA1OGB3eqBq2mlL5+Zl94WqanB5fDW4PCpUVUNedt+Ob3A7TRrbGwICqiYghOb7GwKTxvbu6tIoSJSFCxcu7OoiOpLDEb6zpEZFmdHcHL5zL3VE+5pd3g57/5PjIpAcZ0GFrRlNzR7ER5nw61/0w0X9E3zPOTxIjrPgmrF9YDEZYD353HeHqlBua0av5ChEWny/jUdEmOA4y2SNmga4PRqcLhVeVUCWJfhObiSkJUYiLSECxyqbUG93IzHGjOuvGtgtRmEN6ZsACIGj5U1wqwIWo4LJv+gblqOwwvl7T5IkREaa9O3LYbyhI5yHEgLdcxhvexyrbELB1iMorWgC4Js+PvuiNFx9aS/0yYxHTY09oOO03EcSaTHAbJTPeGd7d8PPZug6n2G87Kwk6iC9U6Px22kXo7DYhnXbS1HX5MbXe8vx7cEq5P6iHy65IAlmo3LO4wjh695ye1TIsgSzQYHZrMBkkM95EyNRMDFAiDqQJEkYeUEyLuqfiK37yrFx13E43SpWbzqMz3eW4VejMnFZVhqMhsAuP2qagMPthcPthSxLiDAbEGFWYFTksO2qpdDBACHqBEaDjKtGZmLMkFR8tfs4tu6rgN3hwcdbj2JToRW/GpWJMUNTYVACH8eiaQJ2hwfNTg8MBt+9ISaDctowX6Jg4TWQEBLO/bBA6F8D+TmSUcEHX/yAbw9W+adFiY0y4aqRGRgzNBUmw7m7ts54XAlQFAkRJgPMRqVdgdSR+NkMXbwGQtTNJcRYMP3KgRg/KhP/3HUC3x2sQoPdjYKvj+Kfu07gimHpyL4oDRHm9n1LCgF4vQKNXg+aJA8UWYLZZIDZoMBklCFJYFcXdRoGCFEQJcRYcP1VA3H1JZn48vsT+PZgFewOD9bv8E3aOHZoKi4fno74aPO5D3YKIXx3x3sdHtjhgSxLMBlkmEwKTIoMRZFDYDwXhRIGCFEXaDkjybm0N7bssWLb/gq4PCo277Hi671WDBuYhCuGp6NPqv5V8DRNwOlW4XSrkCRAliQYDTIMBhkGRYZB8c1jJUvtX8WRCGCAEHWp2JM3J/7qkl7YVlSBrXvL0ejwoLDYhsJiG3qnROGXw9IxfGDSeV3fEAJQhYDqVgG3bx4uSfKNGpNlwCDLvrvsDb4zFZldXxQAXkQPIeF8IQ8I74voiYlRAd1I6FU1FBbbsGWPFVZbs//5SIsBY4akYGxWGpJiO28eLEkCFEmCQZGhGE6epUgyJNl3Y6RyljXn+dkMXbyIThQmDIqMSwen4JILk3GkvBFb95WjqKQGzU4vvtptxVe7rRjUKxZjhqTiov6JAd9PEighAK8Q8Goq0GrG4JazFYMiw2yUYVQUKCe7wKjnYoAQdUOSJGFARiwGZMSiwe7GjgOV2HmgEvV2N4qPN6D4eAMsJgUjL0jGpYOT0TslGlIn3qUuBHwTSmq+O+QBD+ST3V8mgwJLpG8teEX2rS3P0V89AwOEqJuLjTJhwujeuPqSXvjhWB12HKjEgaN1cLpVbCuqwLaiCiTHWTDygmSMvCAJyXERQalLEwKaCnhVLxqa3ahpcP50sf5kF1jLmYsiSzDIEmRZ8l13YcCEBQYIUYiQZQlD+iZgSN8ENDa7sftHm2/m35pmVNc78fm3x/D5t8fQKzkKwwclYfjARCTEBHfdEP/F+lO6wICfusEAX1ed0eD7I58MlJagaX0i1TpkWp7XBBhA3QQDhCgExUSaMG5EBsaNyIDVZsf3P1Rjd7ENDXY3jlfbcbzajrXbStErOQoXD0jERf0TkZoQnDOTs2npBgPQqivMR5Lgm3VYQpswkX1fAIRvRUdN0yAEIJ0cOaYosj90WgLFNxBAhiz7QleWJHapdRIGCFGIy0iKQkZSFHKz++JoeSMKi23YV1KDJofHHybrd5QhOc6CrH4JGNovAX3TYrrVBXAhAAEBCMC3/qLAmVdPOUkDvFABnH1p4JYzHgm+desNsi+M/OEkS5DR8hz8C4MJCEjCdwBFAiRZgqaJ0xb78gWSAHrw7ZkcxhtCwnkoIcBhvB1J0wSOlDdib4kN+4/Uot7edjGkCLOCC3vHY3CfeFzQOw6xOhcUahHs9nUkyf+fU7rMTv5HgoSExEjU1zWffOyLDQgJibHmkJ9in8N4iagNWZYwMDMWAzNjMfXy/jhebcf+I7U4UFoLq60ZDpfqv1kRANITI3FB7zgMyoxF/4zYgNYtCRf+E4mzPC8g/NPEtNaNTuC6DAOEKMxJkoTeKdHonRKNa8b2QV2TCwdL6/DDsTr8eLwebo+G8ppmlNc0Y3OhFbIkoXdqFAZmxmFARgz6psbAbOo5gUKBY4AQ9TDx0WZkX5SG7IvS4FU1lFU24Ydj9Sg+Xo9jVU3QhEBpRRNKK5qwcZfvN+2M5Cj0S4tBv/QY9E2LQVzU+XV5UXhggBD1YAZF9t+wiLF94HR7UWJtxOET9SixNsJqs0MTwPEqO45X2fH13nIAvntT+qREo09qNHqlRiEiqv2zB1PoY4BQ+JCAaIvhpyE1OHlRVIgzdXG33u3kP6SfHkttv3byMBBCnHE4qKYKeDQVQvNtp4Xo2BSLyYCsfgnI6pcAAHC6vTha3uj7U9GIY5V2eFQNDXY39ro0cfYAABGfSURBVNlrsO9Izck99yM5zoLM5ChkJkUhIzkSGUlRiI4wdl1jqNMxQChsSAAiLV33A6vlXgNVE/CoGjRNQBMCqle0zrSzjvrsjpljMRn8Ny8CgKppKLc1o7SyCccqm1BW2QRbvRMCQHW9E9X1Tv+FeQCIiTAiPSkSaYmRSEuIQFpCJFISInrURfpwxgAh6iAtAeCbtbbtD8jEhEiIkzfOydJPd1W37CMACA1QhQa3W4PLq571bKcrKbKMXinR6JUSDVzse87p9qLRpeFASTWs1c04Xm1Hdb0DQgCNDg8aj9Xjh2P1bY4TH21CSnwEUuIjkBxvQXJcBJLjLIiNMoX8sNiehAFCFASKIsN0rplzFQCQEWFq6QbT4FEFHE4v3F6124VJC4vJgMz0KKTE/HRh3ePVUFHbjHKbb3RXRW0zymscsDt8twfWNblR1+Q+LViMiozEWDMSYy1IirUgIdaMxBgzEmItSIg2d/jsw3R+GCBE3VDLdByKDESYFHhVDU63CofTC7UbnpmcymiQ/UOHW7M7PaisdaCy1oHqOgeq6h2oqnOirtEFAcCjaqiodaCi1nHG40ZHGBEfbUJ8jBnxUWbERZsQF21GXJQJcVEmREcYIfMGjaBhgBB1c0L4uo6iLDIiLUZoqua7zqIJeNwq3Cevt4SCKIsRAzKMvlFfrXi8GmoanbDVO2Fr8P1d2+hCTYMLdU0uqCfb1+TwoMnhwbGqM9/1Lkm+6y4xUSbERJgQE2k8+ccXLjGRRkRHGBEVYYTJIHfqFPg9AQOEKIS0zOukKIAJgGQxQNUEvF4Bj+YLEqEJaJqAV9N8czuFyBlLWkIk0hIiT/uapgk0NLtR2+hCbaMvUOqa3KhvcqHe7kZ9kxuuk9eXhAAamj1oaPYA+PmpVYyKjEiLAVERRkRZDIiyGBFhMSDSbEBkq78jzAaosgy3ywuzSeE1mlYYIEQhTAjf+hsmowQTfro+0PIzTtWEbwoOIXwrDXo0eFQVqhYawQL4pmWJjzYjPtqMARln3sbp9qLB7kGD3Y2GZjcam91osHvQ6HCjsdmDpmbfv90ezb+PR9V8AXTKPGE/RwJgNimwmBREWYy4bvxAXHphynm2MHQFLUBKSkowb9481NXVIT4+HkuXLkX//v3bbKOqKp588kls2rQJkiThrrvuwowZM4JVIlHYaAkGWZJgMkgoLK7Guu2lqG9yIy0xEhNG98bgPvH4x8YfsedwDTyqBlUVGNonDhcNSMLmwhOw1TsRG2XClSMzMaRvAg6W1mLT7hOobXQhIcaMK0dm4nhVEzYXlsPlVWE2KBg3Ih05o/ucsaYz7X+247YMGw70GF/vsfrq8KgwG311TL2i/2n7u70q9vxow9Z95ai3uxFpNqBvWjScHhVHyxvhcquQZQlmowKvJuB0edG6d1AAcLpVON0q6prc+NvH+7EhtQx52X0xYlDyeb5roSdos/HeeuutuOGGG3Dttdfio48+wnvvvYe33367zTYffvgh1qxZgzfeeAN1dXWYPn06/ud//ge9e/cO+HU4G2/oCuf2dWXbCoursfKzQ/6RYG6vBlXVkBRjxoGyesjyyRUElZNrnhskJMVHwGI0wKsJeFUVw/onYteP1RAAFEmCy6uiockFh9t356QE35mOKoAJl/Y6LUQOltZi9ZYSKIoMoyKfDCwNowen4NtDVac9P+2KAaeFyNmO0S81GrsP10CCb4i0b0r2wOtwOL2AJCHCrJxWw+A+8XB5VJgsJpRXNsLh8qL4eD22H6iAyaAgNSECDc0eqKqGm68ZHJIh0u1n47XZbCgqKsLf/vY3AEB+fj6eeOIJ1NTUIDEx0b/dJ598ghkzZkCWZSQmJmLixIlYu3Yt7rzzzoBfK9xHYLB9oaur2ratqAKpiZEwGX66N8XtVVFV60BqQkSb+xpbfvWKMPtuyDSeXKxp1482xESZYDEqvmVpAVQaDZAl3/WL1m2raXIjIdbi6yLTfLPZllY2YVDveBhk2Te7LQCPR0NppR1902NhVGT/7LduVcXBsjoMH5TUprCio7VIT4qC0aDg5DS5cHtVWGscSEuIaFODpgkcLKtHbna/Nv8v9pXUnPb/wlbvG/GV1GopYLdXxb6SGgwbmASTUUF8fCTMJ3sId/9YjT5pMYiOMCHCrCAqwjfMeltRBUaFYHfW+XwugxIgVqsVaWlpUBTfm6YoClJTU2G1WtsEiNVqRWZmpv9xRkYGysvL2/VaCQlRHVN0N6X3N4VQEc7t66q2zb/9F13yuq397uYx532Mxwed/w/nBedxjPSkqPM+RrjhXTlERKRLUAIkIyMDFRUVUFXfUDtVVVFZWYmMjIzTtjtx4oT/sdVqRXp6ejBKJCKidgpKgCQlJSErKwsFBQUAgIKCAmRlZbXpvgKAvLw8rFq1CpqmoaamBhs2bEBubm4wSiQionYK2iis4uJizJs3Dw0NDYiNjcXSpUsxcOBAzJ49G3PnzsXw4cOhqioWL16MLVu2AABmz56NmTNnBqM8IiJqp6AFCBERhRdeRCciIl0YIEREpAsDhIiIdGGAEBGRLiE7G+8999yDY8eOQZZlREZG4o9//COysrICmrQxVLz00ktYvnw51qxZg8GDB+P777/H448/DpfLhV69euHZZ59FUlLSuQ/UDeXk5MBkMsFsNgMAHnroIVx55ZVh0UaXy4UlS5Zg69atMJvNGDVqFJ544omw+GweO3YMc+bM8T9ubGxEU1MTtm/fHhbtA4B//vOfWLZs2cnZigXuvfdeTJo0KWzat3HjRixbtgxerxdxcXF4+umn0adPH33tEyGqoaHB/+/PPvtMTJ8+XQghxC233CI+/PBDIYQQH374objlllu6pL7ztXfvXnHHHXeIq6++Whw8eFCoqiomTpwoduzYIYQQ4uWXXxbz5s3r4ir1a2lXa+HSxieeeEI89dRTQtM0IYQQVVVVQojw+Wy29uSTT4pFixYJIcKjfZqmiTFjxvg/m/v37xejRo0SqqqGRfvq6urEZZddJg4fPiyE8LXj9ttvF0Loe/9CNkBa++CDD8R1110nqqurxejRo4XX6xVCCOH1esXo0aOFzWbr4grbx+VyiRtvvFGUlZX5f9Du3r1bTJkyxb+NzWYTo0aN6sIqz8+ZAiQc2tjU1CRGjx4tmpqa2jwfLp/N1lwul8jOzhZ79+4Nm/ZpmiYuu+wysXPnTiGEENu3bxeTJk0Km/bt3r1bTJ482f+4trZWDB48WHf7QrYLCwAeffRRbNmyBUIIvPnmmwFP2tjdLVu2DNOmTWszjf2pE00mJiZC0zT/6WYoeuihhyCEwOjRo/G73/0uLNpYVlaG+Ph4vPTSS9i2bRuioqJw//33w2KxhMVns7UvvvgCaWlpuPjii7F3796waJ8kSXjhhRdwzz33IDIyEna7Ha+//nrY/GwZMGAAqqurUVhYiBEjRmDNmjUAAp/w9lQhfRH9qaeewsaNG/Hggw/imWee6epyOsSuXbuwd+9ezJo1q6tL6VQrV67E6tWr8d5770EIgcWLF3d1SR1CVVWUlZXhoosuwvvvv4+HHnoI9913H5qbm7u6tA733nvv4YYbbujqMjqU1+vFn//8Z7zyyiv45z//iVdffRUPPPBA2Lx/MTExeP755/H000/j+uuvh81mQ2xsrO72hXSAtJg+fTq2bduG9PT0gCZt7M527NiB4uJiTJgwATk5OSgvL8cdd9yBo0ePtplosqamBrIsh8xv5qdqeU9MJhNmzZqF77777rTJNEOxjRkZGTAYDMjPzwcAjBw5EgkJCbBYLCH/2WytoqICO3bswNSpUwEEPmFqd7d//35UVlZi9OjRAIDRo0cjIiICZrM5LNoHAJdffjneeecdvP/++/jXf/1XOJ1O9OrVS1f7QjJA7HY7rFar//EXX3yBuLi4gCdt7M7uuusubN68GV988QW++OILpKen4y9/+QvuvPNOOJ1O7Ny5EwDw7rvvIi8vr4ur1ae5uRmNjb7V+YQQ+OSTT5CVlYVhw4aFfBsTExORnZ3tn8+tpKQENpsN/fv3D/nPZmsffPABxo8fj4QE36qB4fC9BwDp6ekoLy/H4cOHAfjm8LPZbOjXr19YtA8AqqqqAACapuG5557DTTfdhF69eulqX0jOhVVdXY177rkHDocDsiwjLi4ODz/8MC6++OKzTtoYqnJycvDaa69h8ODB+O6777BgwYI2Q1yTk0NvCc2ysjLcd999UFUVmqZh0KBBeOyxx5CamhoWbSwrK8P8+fNRV1cHg8GABx54AOPHjw+rz2Zubi4effRRXHXVVf7nwqV9q1evxhtvvAFJ8q3UN3fuXEycODFs2vfoo4/iu+++g8fjwRVXXIH58+fDbDbral9IBggREXW9kOzCIiKirscAISIiXRggRESkCwOEiIh0YYAQEZEuDBAiItIlpOfCIjrVJZdc4v+3w+GAyWTyz++zaNEiTJs2ratK0y0nJwdPPvkkLr/88q4uhagNBgiFlV27dvn/HQo/eL1eLwyGzv02DMZrUM/ELizqETRNw+uvv46JEyciOzsb999/P+rq6gD4FkkaMmQI3nvvPYwfPx5jx47FO++8g8LCQkydOhVjxoxpM9nj+++/j5tuugmLFy/G6NGjkZeXh61bt/q/3tjYiPnz52PcuHG48sor8fzzz/vnGGrZd8mSJcjOzsby5ctRWlqKW2+9FdnZ2cjOzsbvf/97NDQ0AAD+8Ic/4MSJE/iP//gPXHLJJXjjjTewbdu2NneAA76w/PrrrwEAy5cvx9y5c/HQQw/h0ksvxQcffPCzNRHpxQChHuG///u/sWHDBqxYsQKbNm1CXFzcaTMA7969G+vXr8fzzz+PJUuW4LXXXsNbb72Fjz/+GJ9++im2b9/u37awsBB9+/bFN998g7lz5+Lee+/1B9K8efNgMBiwfv16fPjhh9iyZQtWrVrVZt8+ffpgy5YtuPvuuyGEwG9/+1ts2rQJn376KcrLy7F8+XIAwLPPPovMzEy89tpr2LVrF2bPnh1Qez///HPk5eVh586dmDp16jlrItKDAUI9wrvvvosHH3wQ6enpMJlMuPfee7Fu3Tp4vV7/NnPmzIHZbMa4ceMQGRmJ/Px8JCUlIS0tDWPGjEFRUZF/28TERPzbv/0bjEYjJk+ejAEDBmDjxo2orq7Gl19+ifnz5yMyMhJJSUn493//d3z88cf+fVNTU3HLLbfAYDDAYrGgX79+uOKKK2AymZCYmIjbbrsNO3bsOK/2jho1ChMnToQsy2hqajpnTUR6sGOUeoQTJ05gzpw5kOWffmeSZRk2m83/uPXa62az+bTHrddMSEtL80+2BwCZmZmorKzEiRMn4PV6MW7cOP/XNE1rMy12enp6m9qqq6vx1FNPYefOnbDb7RBCIDY29rza2/o1AqmJSA8GCPUI6enpWLJkiX+dh9aOHTvW7uNVVFRACOEPEavVipycHP8ZzjfffHPWC9etgwcAnnvuOUiShDVr1iA+Ph4bNmz42QW2IiIi4HQ6/Y9VVUVNTc1ZXyOQmoj0YBcW9Qi/+c1v8MILL+D48eMAfItVbdiwQffxampq8Pbbb8Pj8eDTTz9FcXExxo8fj9TUVFxxxRX405/+hKamJmiahtLS0jbXT05lt9sRGRmJmJgYVFRU4M0332zz9eTkZJSVlfkfDxgwAC6XCxs3boTH48Grr74Kt9t91uPrqYkoEAwQ6hFuvfVW5OTk4Pbbb8cll1yCG2+8EYWFhbqPN2LECBw9ehS/+MUv8MILL+DFF1/0L670zDPPwOPxYPLkyRg7dizmzp3rX8TnTO69914UFRVhzJgxuOuuuzBp0qQ2X7/rrrvw6quvYsyYMfjLX/6CmJgYLFiwAI899hiuuuoqREREnNYtdqr21kQUCK4HQtRO77//PlatWoV33nmnq0sh6lI8AyEiIl0YIEREpAu7sIiISBeegRARkS4MECIi0oUBQkREujBAiIhIFwYIERHpwgAhIiJd/h8ZPsfu8SN/bwAAAABJRU5ErkJggg==\n", + "text/plain": [ + "
" + ] + }, + "metadata": {}, + "output_type": "display_data" + } + ], "source": [ "sns.set(color_codes=True)\n", "plt.xlim(30,90)\n", @@ -624,6 +884,174 @@ "source": [ "**I think I have managed to correctly compute and plot the uncertainty of my prediction.** Although the shaded area seems very similar to [the one obtained by with R](https://app-learninglab.inria.fr/moocrr/gitlab/moocrr-session3/moocrr-reproducibility-study/tree/master/challenger.pdf), I can spot a few differences (e.g., the blue point for temperature 63 is outside)... Could this be a numerical error ? Or a difference in the statistical method ? It is not clear which one is \"right\"." ] + }, + { + "cell_type": "code", + "execution_count": 10, + "metadata": {}, + "outputs": [ + { + "data": { + "text/html": [ + "
\n", + "\n", + "\n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + "
TemperatureInterceptFrequency
030.011.0
130.511.0
231.011.0
331.511.0
432.011.0
............
11688.011.0
11788.511.0
11889.011.0
11989.511.0
12090.011.0
\n", + "

121 rows × 3 columns

\n", + "
" + ], + "text/plain": [ + " Temperature Intercept Frequency\n", + "0 30.0 1 1.0\n", + "1 30.5 1 1.0\n", + "2 31.0 1 1.0\n", + "3 31.5 1 1.0\n", + "4 32.0 1 1.0\n", + ".. ... ... ...\n", + "116 88.0 1 1.0\n", + "117 88.5 1 1.0\n", + "118 89.0 1 1.0\n", + "119 89.5 1 1.0\n", + "120 90.0 1 1.0\n", + "\n", + "[121 rows x 3 columns]" + ] + }, + "execution_count": 10, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "data_pred" + ] + }, + { + "cell_type": "code", + "execution_count": 11, + "metadata": {}, + "outputs": [ + { + "data": { + "text/plain": [ + "0 1.0\n", + "1 1.0\n", + "2 1.0\n", + "3 1.0\n", + "4 1.0\n", + " ... \n", + "116 1.0\n", + "117 1.0\n", + "118 1.0\n", + "119 1.0\n", + "120 1.0\n", + "Length: 121, dtype: float64" + ] + }, + "execution_count": 11, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "data_pred = pd.DataFrame({'Temperature': np.linspace(start=30, stop=90, num=121), 'Intercept': 1})\n", + "logmodel.predict(data_pred)" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [] } ], "metadata": { diff --git a/module4/shuttle.csv b/module4/shuttle.csv new file mode 100644 index 0000000000000000000000000000000000000000..9412ab30682565afb0ca3b4bcc6b0f99f632bbf2 --- /dev/null +++ b/module4/shuttle.csv @@ -0,0 +1,24 @@ +Date,Count,Temperature,Pressure,Malfunction +4/12/81,6,66,50,0 +11/12/81,6,70,50,1 +3/22/82,6,69,50,0 +11/11/82,6,68,50,0 +4/04/83,6,67,50,0 +6/18/82,6,72,50,0 +8/30/83,6,73,100,0 +11/28/83,6,70,100,0 +2/03/84,6,57,200,1 +4/06/84,6,63,200,1 +8/30/84,6,70,200,1 +10/05/84,6,78,200,0 +11/08/84,6,67,200,0 +1/24/85,6,53,200,2 +4/12/85,6,67,200,0 +4/29/85,6,75,200,0 +6/17/85,6,70,200,0 +7/2903/85,6,81,200,0 +8/27/85,6,76,200,0 +10/03/85,6,79,200,0 +10/30/85,6,75,200,2 +11/26/85,6,76,200,0 +1/12/86,6,58,200,1