Homework 02 /exo5(Improved version).ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Analysis of the risk of failure of the Challenger shuttle O-rings"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## I. First attempt of Improvement"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The original paper shows commendable efforts to understand the risk of O-ring failure, but some errors in reasoning led to an underestimation of the risk. We will try to improve some aspects of the analysis to provide more reliable results."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Importing necessary libraries"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 22,
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import pandas as pd\n",
    "import statsmodels.api as sm\n",
    "import matplotlib.pyplot as plt"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Loading data\n",
    "Data from previous spaceflights, including temperatures, number of O-rings observed, and number of malfunctions, are loaded into a DataFrame data."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 23,
   "metadata": {},
   "outputs": [],
   "source": [
    "data = pd.read_csv(\"shuttle.csv\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Graphical data inspection"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Calculation of failure frequency\n",
    "The \"failure frequency\" of O-rings for each flight is calculated by dividing the number of observed malfunctions (Malfunction) by the total number of O-rings tested (Count). This gives an idea of ​​the proportion of seals that failed during each flight."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 24,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>Date</th>\n",
       "      <th>Count</th>\n",
       "      <th>Temperature</th>\n",
       "      <th>Pressure</th>\n",
       "      <th>Malfunction</th>\n",
       "      <th>Frequency</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>4/12/81</td>\n",
       "      <td>6</td>\n",
       "      <td>66</td>\n",
       "      <td>50</td>\n",
       "      <td>0</td>\n",
       "      <td>0.000000</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>11/12/81</td>\n",
       "      <td>6</td>\n",
       "      <td>70</td>\n",
       "      <td>50</td>\n",
       "      <td>1</td>\n",
       "      <td>0.166667</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>3/22/82</td>\n",
       "      <td>6</td>\n",
       "      <td>69</td>\n",
       "      <td>50</td>\n",
       "      <td>0</td>\n",
       "      <td>0.000000</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td>11/11/82</td>\n",
       "      <td>6</td>\n",
       "      <td>68</td>\n",
       "      <td>50</td>\n",
       "      <td>0</td>\n",
       "      <td>0.000000</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4</th>\n",
       "      <td>4/04/83</td>\n",
       "      <td>6</td>\n",
       "      <td>67</td>\n",
       "      <td>50</td>\n",
       "      <td>0</td>\n",
       "      <td>0.000000</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>5</th>\n",
       "      <td>6/18/82</td>\n",
       "      <td>6</td>\n",
       "      <td>72</td>\n",
       "      <td>50</td>\n",
       "      <td>0</td>\n",
       "      <td>0.000000</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>6</th>\n",
       "      <td>8/30/83</td>\n",
       "      <td>6</td>\n",
       "      <td>73</td>\n",
       "      <td>100</td>\n",
       "      <td>0</td>\n",
       "      <td>0.000000</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>7</th>\n",
       "      <td>11/28/83</td>\n",
       "      <td>6</td>\n",
       "      <td>70</td>\n",
       "      <td>100</td>\n",
       "      <td>0</td>\n",
       "      <td>0.000000</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>8</th>\n",
       "      <td>2/03/84</td>\n",
       "      <td>6</td>\n",
       "      <td>57</td>\n",
       "      <td>200</td>\n",
       "      <td>1</td>\n",
       "      <td>0.166667</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>9</th>\n",
       "      <td>4/06/84</td>\n",
       "      <td>6</td>\n",
       "      <td>63</td>\n",
       "      <td>200</td>\n",
       "      <td>1</td>\n",
       "      <td>0.166667</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>10</th>\n",
       "      <td>8/30/84</td>\n",
       "      <td>6</td>\n",
       "      <td>70</td>\n",
       "      <td>200</td>\n",
       "      <td>1</td>\n",
       "      <td>0.166667</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>11</th>\n",
       "      <td>10/05/84</td>\n",
       "      <td>6</td>\n",
       "      <td>78</td>\n",
       "      <td>200</td>\n",
       "      <td>0</td>\n",
       "      <td>0.000000</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>12</th>\n",
       "      <td>11/08/84</td>\n",
       "      <td>6</td>\n",
       "      <td>67</td>\n",
       "      <td>200</td>\n",
       "      <td>0</td>\n",
       "      <td>0.000000</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>13</th>\n",
       "      <td>1/24/85</td>\n",
       "      <td>6</td>\n",
       "      <td>53</td>\n",
       "      <td>200</td>\n",
       "      <td>2</td>\n",
       "      <td>0.333333</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>14</th>\n",
       "      <td>4/12/85</td>\n",
       "      <td>6</td>\n",
       "      <td>67</td>\n",
       "      <td>200</td>\n",
       "      <td>0</td>\n",
       "      <td>0.000000</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>15</th>\n",
       "      <td>4/29/85</td>\n",
       "      <td>6</td>\n",
       "      <td>75</td>\n",
       "      <td>200</td>\n",
       "      <td>0</td>\n",
       "      <td>0.000000</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>16</th>\n",
       "      <td>6/17/85</td>\n",
       "      <td>6</td>\n",
       "      <td>70</td>\n",
       "      <td>200</td>\n",
       "      <td>0</td>\n",
       "      <td>0.000000</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>17</th>\n",
       "      <td>7/29/85</td>\n",
       "      <td>6</td>\n",
       "      <td>81</td>\n",
       "      <td>200</td>\n",
       "      <td>0</td>\n",
       "      <td>0.000000</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>18</th>\n",
       "      <td>8/27/85</td>\n",
       "      <td>6</td>\n",
       "      <td>76</td>\n",
       "      <td>200</td>\n",
       "      <td>0</td>\n",
       "      <td>0.000000</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>19</th>\n",
       "      <td>10/03/85</td>\n",
       "      <td>6</td>\n",
       "      <td>79</td>\n",
       "      <td>200</td>\n",
       "      <td>0</td>\n",
       "      <td>0.000000</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>20</th>\n",
       "      <td>10/30/85</td>\n",
       "      <td>6</td>\n",
       "      <td>75</td>\n",
       "      <td>200</td>\n",
       "      <td>2</td>\n",
       "      <td>0.333333</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>21</th>\n",
       "      <td>11/26/85</td>\n",
       "      <td>6</td>\n",
       "      <td>76</td>\n",
       "      <td>200</td>\n",
       "      <td>0</td>\n",
       "      <td>0.000000</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>22</th>\n",
       "      <td>1/12/86</td>\n",
       "      <td>6</td>\n",
       "      <td>58</td>\n",
       "      <td>200</td>\n",
       "      <td>1</td>\n",
       "      <td>0.166667</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "        Date  Count  Temperature  Pressure  Malfunction  Frequency\n",
       "0    4/12/81      6           66        50            0   0.000000\n",
       "1   11/12/81      6           70        50            1   0.166667\n",
       "2    3/22/82      6           69        50            0   0.000000\n",
       "3   11/11/82      6           68        50            0   0.000000\n",
       "4    4/04/83      6           67        50            0   0.000000\n",
       "5    6/18/82      6           72        50            0   0.000000\n",
       "6    8/30/83      6           73       100            0   0.000000\n",
       "7   11/28/83      6           70       100            0   0.000000\n",
       "8    2/03/84      6           57       200            1   0.166667\n",
       "9    4/06/84      6           63       200            1   0.166667\n",
       "10   8/30/84      6           70       200            1   0.166667\n",
       "11  10/05/84      6           78       200            0   0.000000\n",
       "12  11/08/84      6           67       200            0   0.000000\n",
       "13   1/24/85      6           53       200            2   0.333333\n",
       "14   4/12/85      6           67       200            0   0.000000\n",
       "15   4/29/85      6           75       200            0   0.000000\n",
       "16   6/17/85      6           70       200            0   0.000000\n",
       "17   7/29/85      6           81       200            0   0.000000\n",
       "18   8/27/85      6           76       200            0   0.000000\n",
       "19  10/03/85      6           79       200            0   0.000000\n",
       "20  10/30/85      6           75       200            2   0.333333\n",
       "21  11/26/85      6           76       200            0   0.000000\n",
       "22   1/12/86      6           58       200            1   0.166667"
      ]
     },
     "execution_count": 24,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Do not exclude flights without failure\n",
    "data[\"Frequency\"] = data.Malfunction / data.Count\n",
    "data"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Including all flights, including those without failures, provides a more complete view of the operating conditions of the O-rings. This avoids biasing the analysis by focusing only on malfunction cases and provides a better understanding of the conditions under which the seals operate properly."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Initial data visualization"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 25,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAY0AAAEWCAYAAACaBstRAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMi4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvIxREBQAAIABJREFUeJzt3XmYHVWd//H3hwYkBBQBCWZhESMaQZbERAcHg7Lqg+DCADKgKERGEZkZGUH5OaA44iCOjqIhIC4goCJE1GhYpEVZEyBDAI0GiGRhhwCNEZLw/f1x6obqy+3b53a6um93Pq/n6advVZ1T9T213G9tt0oRgZmZWY71BjsAMzMbOpw0zMwsm5OGmZllc9IwM7NsThpmZpbNScPMzLKts0lD0lRJS9ai/nRJ/68/Y2owjZD02h6GHSHpqj6O9zRJF61ddGa96207k7SHpL9I6pJ0cC/j+kdJC0rdiyTt3Z/xNpn2CEm/kPSUpJ8OxDTb1ZBOGsVKs6JY4R6S9H1Jm1QwnQ9L+kO5X0QcFxFf7O9p5YqIH0XEvgM93aGQcIr14IzBjqM/rO3OzRDwBeBbEbFJRMxsVjAifh8ROw5QXPU+AIwCtoiIQwYphrYwpJNG4cCI2ATYFdgNOGWQ4zHLoqTSbVDS+lWOvx9sC9w9kBPs4zzZFvhzRKzqx3EOTRExZP+ARcDepe7/Bn5V6n4Z8FXgAeBhYDowohg2FVhSKnsycC/wDHAP8N6i/xuAvwOrgS5gedH/+8AZpfrHAguBJ4ArgdGlYQEcB/wFeBI4B1Ax7LXA74CngMeAH2fW+zDwh7qyJwD3FeM5C1ivh/l2GnAZ8OOivbcDu5SGjwZ+BjwK3A+cUPTfH3geWFnMi/8D9gLml+peA9xa6v4DcHCz8RbD1istg8eBnwCbF8O2K9r3oWJZPgZ8roe2TSvie76I8RcZ0z4N+ClwUTE/5gOvI+2APAIsBvYtle8EvgzcWiy3n9diLYa/BbgRWF7Mo6l1db8E3ACsKJb/0cAfi2nfB3ysKDuyKPNC0Zauoh3fp/u6N5Xu6/Ii4DPAncBzwPrN2t9gHr4buAN4umj7aaVhTZcFMKKI70nSdnRSOba66dxbtG1F0baX9TQvmrRz7x62x36bJ8DpdF/vP0ra/m4A/oe0zZ9RlP1IEf+TwGxg29J49gH+RFpnvkXa7o8prYMXNZjP6xfdrwC+CzwILAXOADrK3wWk77oni7YcUBrX5sD3gGXF8JlF/7tIO921chsUy3PXpt+7VX+xV/lXt9KMJW3s3ygN/zrpC3xzYFPgF8CXe1ipDilWovWAQ4FngVeXF0rdtNespMA7ipm9O2nF/yZwfalsAL8ENgO2KVbS/YthlwCfK6a7EfC2zHrdYirKXle0dRvgz7UVssF8O420AXygWFE+XaxoGxRx3AZ8HtgQeA1p492vh5V7I9JGvyVpQ3yoWDk3JX2BrAC2yBjvicDNxXJ8GXAucEndBnReMc5dSBv+G3po35plU3TntOnvwH5FG35YzI/PFfPkWOD+0vg6SRvuTqQv9p/V5gkwhpT03lVMd5+i+1Wlug8AbyymtQHpS3oHQMDbgb8BuzdaT1v4gpwHjCvmV9P2N5h/U4Gdi3pvIu1wHZyzLIAzgd+T1sNxpC+mhkmjfhsuurPnBa0njbWZJ6fRfb3/MLAK+GSxHEcAB5N2HN9Q9DsVuLEovyUpCde2uX8t6ucmjZmkbWIksBVph+VjpVhWktbTDuBfSNtgbQfzV6QdxFcW03570f8/6L6TehClHcAel9lgfNn311+xInSR9koCuBbYrBgm0hf/DqXyb6XY+OtXqgbjngccVFoozZLGd4H/Lg3bpFiI2xXdQfdk8BPg5OLzD4EZwNgGMTSr1y2mouz+pe6PA9c22QBuLnWvR9qD+UdgCvBAXflTgO81WrmLfr8H3kfaw76qiHN/0lHInUWZ3sb7R+CdpWGvLubh+ry4AY0tDb8VOKyH9q1ZNpnTPg24ujTswGK9qu3JbVpMv7ZudQJnlspPIO2JdpD2Zi+sm9Zs4EOlul/oZb2eCXyqp/W0Qfu6lSFtFx/JbX/GdvZ14H+Kz02XBemLt7weTquPv8E2vHeT4T3OC1pPGn2eJzROGvX1fw18tG67+hvp1NZRdN/mBCwhI2mQrqU8R3GWpBh+OHBdKZaFpWEbF3W3Jm1HLwCvbNCm0aTvzpcX3ZcB/9Hb+jAczsMdHBHXSHo7cDEpoy8HXkWaebdJqpUVacN+CUlHAf9GWliQvvi3zIxhNOkUDwAR0SXpcdJe56Ki90Ol8n8rxg8p238RuFXSk8DZEXFBqWxP9RpZXPr81yKuXstGxAvFxdbRpJVttKTlpbIdpMTQk99RbKTF5ydJe4nPFd2QNpxm490WuELSC6Xhq0kbTE0r86Kst2lD2puuWQE8FhGrS90U06uNo35eb0BaX7YFDpF0YGn4BqSjQBrURdIBwH+STomtR1pv52e1rGflaeS0vxzPFNIRw06kvfCXkU7flfW0LEbz0nmTraJ5UdPneZIxvto4vyHp7FI/kb4Hus2XiAhJ9fV7si1pHXqw9F22Xt301yyPiPhbUW4T0hHfExHxZP1II2KZpBuA90u6AjgA+FRvwQyHpAFARPxO0vdJ5/UOJp0uWgG8MSKWNqsraVvS4fY7gZsiYrWkeaQFDumLtJllpAVbG99I0imZptMt4n6IdFiJpLcB10i6PiIW9la3gXG8eFFxmyKuZmVr8a5HOi20jHTIfH9EjO8p5Ab9fgecTTrtciYpaZxHShrnFGUW9zLexaQ9wRvqB0jarkk7cmLsbdp9Ma70eRvSUdFjxbQujIhjc+KT9DLS6a2jgJ9HxEpJM2m+7j1L+jKt2brZNGi9/ReTzrkfEBF/l/R18negHuSl62GWjHnRTNXzpJFG69mXIuJH9QUljaf7Nie6r0PN4l9M2pa2jB4uxDexGNhc0mYRsbzB8B8Ax5BywU29fVfC8Lh7quzrwD6Sdo2IF0hfXP8jaSsASWMk7deg3kjSCvBoUe5o0l5WzcPAWEkb9jDdi4GjJe1arPj/BdwSEYt6C1jSIZLGFp1PFnGsblKlmZMkvVLSONIew4+blJ0o6X3FXR8nklbKm0mnGp6W9Jni3vQOSTtJenNR72Fgu7q7fm4EdgQmky6C301KolOA64syvY13OvClIoEj6VWSDurjfHiYdI66prdp98U/S5ogaWPSbaOXFUcmFwEHStqvmM5GxW2zY3sYT21P/lFgVbGnXb6V+mFgC0mvKPWbB7xL0uaStiYtv2Zabf+mpL3Tv0uaDHywl/GX/QQ4pVgPx5LO+efqbV40U/U8yTGd1PY3Akh6haTa7bm/At5Y2uZOoHtimAfsKWmbYlmvuQs0Ih4knfY9W9LLJa0naYfi7EpTRd1fA98ulskGkvYsFZlJuhb7KdKp8l4Nq6QREY+SGl770d1nSBembpb0NOnOnpfc5x0R95D2lG8ibaQ7k+6MqPktac/pIUmPNah/bTHNn5H2tHYADssM+83ALZK6SBftPxUR92fWrfdz0sW9eaSV9Lu9lD2UlKiOBN4XESuLL74DSbcw30/aez6fdPcGvHia4nFJtwNExLOk03N3R8TzxfCbgL9GxCNFmd7G+42i/VdJeoaUwKb0bTbwXWCCpOWSZmZMuy8uJJ1Hf4h0M8AJABGxmHRB8bOkL7/FpDuIGm5rEfFMUfcnpGXxQdJ8qA3/E+lmifuK9owupv1/pFOfV9F85yBn3tf7OPCFYjl8vogt1+mkU1L3F7FdmFuxt3nRi6rnSa8i4grgK8ClxffNXaRTPkTEY6Sbbc4k3RgxntJ3TERcXcR8J2kb/mXd6I8iJdV7SPPmMtL1ihxHko6E/0S6G3BNQo2IFaTvre2By3NGVru6bkOcpADG9/G0lrVAUifpouX5gx2LDV3tsh5J+jzwuoj455zyw+aahpmZtUbS5qTfnRyZW2dYnZ4yM7M8ko4lnT79dURc31v5NfV8esrMzHL5SMPMzLINq2saW265ZWy33XaDHUY3zz77LCNHjhzsMPqd2zX0DNe2uV1r57bbbnssIl6VW35YJY3tttuOuXPnDnYY3XR2djJ16tTBDqPfuV1Dz3Btm9u1diS19Kt9n54yM7NsThpmZpbNScPMzLI5aZiZWbaqXzW5v6QFkhZKOrnB8IMk3SlpnqS5xVNea8MWSZpfG1ZlnGZmlqeyu6ckdZAei70P6T0LcyRdWTwcsOZa4Mri2fJvIj2o7PWl4XsVD/oyM7M2UOWRxmTS26TuK558einp6Z9rRERXvPiT9Nrjyc3MrE1V9hgRSR8gvfbxmKL7SGBKRBxfV+69wJdJ7719d0TcVPS/nxffL3FuRMzoYTrTSK+UZNSoURMvvfTSStrTV11dXWyySe4L5oYOt2voGa5tc7vWzl577XVbREzKrtDb+2D7+kd6dvz5pe4jgW82Kb8ncE35/bXF/61Iz8nfs7dpTpw4MdrNddddN9ghVMLtGnqGa9vcrrUDzI0WvturPD21hO6vM6y9TrShSE9Z3EHSlkX3suL/I8AVpNNdZmY2iKpMGnOA8ZK2L16Tehh1b+GS9NriXblI2p30ZqrHJY2UtGnRfyTplY93VRirmZllqOzuqYhYJel4YDbQAVwQEXdLOq4YPh14P3CUpJXACuDQiAhJo4ArinyyPnBxRPymqljNzCxPpQ8sjIhZwKy6ftNLn79Ceqdufb37gF2qjM3MzFq3zv8ifOYdS9njzN+y/cm/Yo8zf8vMO5YOdkhmlqG27c5f+pS33QE0rB6N3qqZdyzllMvns2LlagCWLl/BKZfPB+Dg3cYMZmhm1kS3bXect92BtE4faZw1e8GahFGzYuVqzpq9YJAiMrMc3nYHzzqdNJYtX9FSfzNrD952B886nTRGbzaipf5m1h687Q6edTppnLTfjozYoKNbvxEbdHDSfjsOUkRmlsPb7uBZpy+E1y6YnTV7AcuWr2D0ZiM4ab8dfSHNrM2Vt114hjHedgfMOp00IK18XtHMhp7attvZ2cknj5g62OGsM9bp01NmZtYaJw0zM8vmpGFmZtmcNMzMLJuThpmZZXPSMDOzbE4aZmaWzUnDzMyyOWmYmVk2Jw0zM8tWadKQtL+kBZIWSjq5wfCDJN0paZ6kuZLellvXzMwGXmVJQ1IHcA5wADABOFzShLpi1wK7RMSuwEeA81uoa2ZmA6zKI43JwMKIuC8ingcuBQ4qF4iIroiIonMkELl1zcxs4FX5lNsxwOJS9xJgSn0hSe8FvgxsBby7lbpF/WnANIBRo0bR2dm5tnH3q66urraLqT+4XUPPcG2b2zWwqkwaatAvXtIj4grgCkl7Al8E9s6tW9SfAcwAmDRpUkydOrWv8Vais7OTdoupP7hdQ89wbZvbNbCqPD21BBhX6h4LLOupcERcD+wgactW65qZ2cCoMmnMAcZL2l7ShsBhwJXlApJeK0nF592BDYHHc+qamdnAq+z0VESsknQ8MBvoAC6IiLslHVcMnw68HzhK0kpgBXBocWG8Yd2qYjUzszyVvu41ImYBs+r6TS99/grwldy6ZmY2uPyLcDMzy+akYWZm2Zw0zMwsm5OGmZllc9IwM7NsThpmZpbNScPMzLI5aZiZWTYnDTMzy+akYWZm2Zw0zMwsm5OGmZllc9IwM7NsThpmZpbNScPMzLI5aZiZWTYnDTMzy5aVNCSNkLRj1cGYmVl76zVpSDoQmAf8pujeVdKVOSOXtL+kBZIWSjq5wfAjJN1Z/N0oaZfSsEWS5kuaJ2lufpPMzKwqOe8IPw2YDHQCRMQ8Sdv1VklSB3AOsA+wBJgj6cqIuKdU7H7g7RHxpKQDgBnAlNLwvSLisYwYzcxsAOScnloVEU/1YdyTgYURcV9EPA9cChxULhARN0bEk0XnzcDYPkzHzMwGSM6Rxl2SPgh0SBoPnADcmFFvDLC41L2E7kcR9T4K/LrUHcBVkgI4NyJmNKokaRowDWDUqFF0dnZmhDZwurq62i6m/uB2DT3DtW1u1wCLiKZ/wMbAl4A5xd8ZwEYZ9Q4Bzi91Hwl8s4eyewF/BLYo9Rtd/N8K+D9gz96mOXHixGg311133WCHUAm3a+gZrm1zu9YOMDd6+W4t/zU90iiuS5weEScBn2sxHy0BxpW6xwLLGkzjTcD5wAER8Xitf0QsK/4/IukK0umu61uMwczM+lHTaxoRsRqY2MdxzwHGS9pe0obAYUC3u64kbQNcDhwZEX8u9R8padPaZ2Bf4K4+xmFmZv0k55rGHcUttj8Fnq31jIjLm1WKiFWSjgdmAx3ABRFxt6TjiuHTgc8DWwDflgTpovskYBRwRdFvfeDiiPhNq40zM7P+lZM0NgceB95R6hekI4SmImIWMKuu3/TS52OAYxrUuw/Ypb6/mZkNrl6TRkQcPRCBmJlZ++s1aUj6HunIopuI+EglEZmZWdvKOT31y9LnjYD30uAuKDMzG/5yTk/9rNwt6RLgmsoiMjOzttWXR6OPB7bp70DMzKz95VzTeIZ0TUPF/4eAz1Qcl5mZtaGc01ObDkQgZmbW/npMGpJ2b1YxIm7v/3DMzKydNTvSOLvJsKD7j/3MzGwd0GPSiIi9BjIQMzNrf81OT70jIn4r6X2Nhvf27CkzMxt+mp2eejvwW+DABsOynj1lZmbDS7PTU/9Z/Pezp8zMDMh7jAiS3g28kfQYEQAi4gtVBWVmZu2p11+ES5oOHAp8kvQDv0OAbSuOy8zM2lDOY0T+ISKOAp6MiNOBt9L9Na5mZraOyEkaK4r/f5M0GlgJbF9dSGZm1q6yHo0uaTPgLOB20p1T51UalZmZtaWcZ099sfj4M0m/BDaKiKeqDcvMzNpRj6enJP1X6fM+ABHxXCsJQ9L+khZIWijp5AbDj5B0Z/F3o6RdcuuamdnAa3ZNY//S56+0OmJJHcA5wAHABOBwSRPqit0PvD0i3gR8EZjRQl0zMxtgfXkJU67JwMKIuC8ingcuBQ4qF4iIGyPiyaLzZmBsbl0zMxt4za5pbCXp30i/zah9XiMivtbLuMcAi0vdS4ApTcp/FPh1q3UlTQOmAYwaNYrOzs5ewhpYXV1dbRdTf3C7hp7h2ja3a2A1SxrnAZs2+JxLDfpFw4LSXqSk8bZW60bEDIrTWpMmTYqpU6e2GGa1Ojs7abeY+oPbNfQM17a5XQOr2bOnTl/LcS+h+48AxwLL6gtJehNwPnBARDzeSl0zMxtYVV7TmAOMl7S9pA2Bw4ArywUkbUN6Wu6REfHnVuqamdnAy3pgYV9ExCpJxwOzgQ7ggoi4W9JxxfDpwOeBLYBvSwJYFRGTeqpbVaxmZpan2UuYPhUR35C0R0Tc0JeRR8QsYFZdv+mlz8cAx+TWNTOzwdXs9FTtPRrfHIhAzMys/TU7PfVHSYuAV0m6s9RfQBQ/yDMzs3VIs7unDpe0Nem6wnsGLiQzM2tXTS+ER8RDwC7FHUyvK3oviIiVlUdmZmZtp9e7pyS9HfghsIh0amqcpA9FxPUVx2ZmZm0m55bbrwH7RsQCAEmvAy4BJlYZmJmZtZ+cH/dtUEsYAMWP8DaoLiQzM2tXOUcacyV9F7iw6D4CuK26kMzMrF3lJI1/AT4BnEC6pnE98O0qgzIzs/aU87rX50jXNXp7FLpVbOYdSzlr9gKWLV/B6M1GcNJ+O3LwbmPWuuxQM5zbVoVTZ87nklsWc+JOK/noKbM4fMo4zjh458EOy4aoyp49Zf1r5h1LOeXy+axYuRqApctXcMrl8wFe8oXZStmhZji3rQqnzpzPRTc/sKZ7dcSabicO64sqn3Jr/eis2QvWfFHWrFi5mrNmL1irskPNcG5bFS65ZXFL/c1602vSkLTTQARizS1bviK7fytlh5rh3LYqrI6G7y7rsb9Zb3KONKZLulXSxyVtVnlE1tDozUZk92+l7FAznNtWhQ41eglmz/3NetNr0oiIt5Fusx1Huv32Ykn7VB6ZdXPSfjsyYoOObv1GbNDBSfvtuFZlh5rh3LYqHD5lXEv9zXqTdSE8Iv4i6VRgLvC/wG5Kb036bERcXmWAltQu8ubcNdRK2aFmOLetCrWL3bVrGB2S756ytZLz7Kk3kd6t8W7gauDAiLhd0mjgJtLrWm0AHLzbmOwvx1bKDjXDuW1VOOPgnTnj4J3p7Ozk3iOmDnY4NsTlHGl8CziPdFSx5mpjRCwrjj7MzGwdkXMh/F3AxbWEIWk9SRsDRMSFzSpK2l/SAkkLJZ3cYPjrJd0k6TlJn64btkjSfEnzJM3Nb5KZmVUlJ2lcA5RvTdm46NeUpA7gHOAAYAJwuKQJdcWeID2e5Ks9jGaviNg1IiZlxGlmZhXLSRobRURXraP4vHFGvcnAwoi4LyKeBy4FDioXiIhHImIO4Jc6mZkNATnXNJ6VtHtE3A4gaSKQ80uqMUD5Z6dLgCktxBbAVZICODciZjQqJGkaMA1g1KhRdHZ2tjCJ6nV1dbVdTP3B7Rp6hmvb3K6BlZM0TgR+KmlZ0f1q4NCMeo1+PdTKz1D3KC62bwVcLelPjd4WWCSTGQCTJk2KqVOntjCJ6nV2dtJuMfUHt2voGa5tc7sGVs5TbudIej2wIykR/CnzHeFLSD8IrBkLLOuhbKPpLiv+PyLpCtLpLr9i1sxsEOU+sPDNwJuA3UgXtI/KqDMHGC9pe0kbAocBV+ZMTNJISZvWPgP7AndlxmpmZhXJ+XHfhcAOwDyg9njRAH7YrF5ErJJ0PDAb6AAuiIi7JR1XDJ8uaWvSr8xfDrwg6UTSnVZbAlekH52zPumW39/0oX1mZtaPcq5pTAImRLT+WMyImAXMqus3vfT5IdJpq3pPA7u0Oj0zM6tWzumpu4Ctqw7EzMzaX86RxpbAPZJuBZ6r9YyI91QWlZmZtaWcpHFa1UGYmdnQkHPL7e8kbQuMj4hriudOdfRWz8zMhp+c170eC1wGnFv0GgPMrDIoMzNrTzkXwj8B7EG6o4mI+AuwVZVBmZlZe8pJGs8VDxwEQNL6tPY4EDMzGyZyksbvJH0WGFG8G/ynwC+qDcvMzNpRTtI4GXgUmA98jPRjPb+xz8xsHZRz99QLpNe9nld9OGZm1s5ynj11Pw2uYUTEayqJyMzM2lbus6dqNgIOATavJhwzM2tnvV7TiIjHS39LI+LrwDsGIDYzM2szOaendi91rkc68ti0sojMzKxt5ZyeOrv0eRWwCPinSqIxM7O2lnP31F4DEYiZmbW/nNNT/9ZseER8rf/CMTOzdpZ799SbefH93gcC1wOLqwrKzMzaU84vwrcEdo+If4+IfwcmAmMj4vSIOL1ZRUn7S1ogaaGkkxsMf72kmyQ9J+nTrdQ1M7OBl5M0tgGeL3U/D2zXWyVJHcA5wAHABOBwSRPqij0BnAB8tQ91zcxsgOWcnroQuFXSFaRfhr8X+GFGvcnAwoi4D0DSpcBBwD21AhHxCPCIpHe3WtfMzAZezt1TX5L0a+Afi15HR8QdGeMeQ/frHkuAKZlxrU1dMzOrSM6RBsDGwNMR8T1Jr5K0fUTc30sdNeiX+x6O7LqSpgHTAEaNGkVnZ2fmJAZGV1dX28XUH9yuoWe4ts3tGlg5t9z+J+kOqh2B7wEbABeR3ubXzBJgXKl7LLAsM67suhExA5gBMGnSpJg6dWrmJAZGZ2cn7RZTf3C7hp7h2ja3a2DlXAh/L/Ae4FmAiFhG3mNE5gDjJW0vaUPgMF68bbfKumZmVpGc01PPR0RICgBJI3NGHBGrJB0PzAY6gAsi4m5JxxXDp0vaGpgLvBx4QdKJwISIeLpR3ZZbZ2Zm/SonafxE0rnAZpKOBT5C5guZImIW6U1/5X7TS58fIp16yqprZmaDK+fuqa8W7wZ/mnRd4/MRcXXlkZmZWdtpmjSKH9nNjoi9AScKM7N1XNML4RGxGvibpFcMUDxmZtbGcq5p/B2YL+lqijuoACLihMqiMjOztpSTNH5V/JmZ2Tqux6QhaZuIeCAifjCQAZmZWftqdk1jZu2DpJ8NQCxmZtbmmiWN8vOfXlN1IGZm1v6aJY3o4bOZma2jml0I30XS06QjjhHFZ4ruiIiXVx6dmZm1lR6TRkR0DGQgZmbW/nKecmtmZgY4aZiZWQucNMzMLJuThpmZZXPSMDOzbE4aZmaWzUnDzMyyOWmYmVm2SpOGpP0lLZC0UNLJDYZL0v8Ww++UtHtp2CJJ8yXNkzS3yjjNzCxPzvs0+qR4Vew5wD7AEmCOpCsj4p5SsQOA8cXfFOA7xf+avSLisapiNDOz1lR5pDEZWBgR90XE88ClwEF1ZQ4CfhjJzcBmkl5dYUxmZrYWKjvSAMYAi0vdS+h+FNFTmTHAg6Qn614lKYBzI2JGo4lImgZMAxg1ahSdnZ39Enx/6erqaruY+oPbNfQM17a5XQOryqShBv3qH7HerMweEbFM0lbA1ZL+FBHXv6RwSiYzACZNmhRTp05di5D7X2dnJ+0WU39wu4ae4do2t2tgVXl6agkwrtQ9FliWWyYiav8fAa4gne4yM7NBVGXSmAOMl7S9pA2Bw4Ar68pcCRxV3EX1FuCpiHhQ0khJmwJIGgnsC9xVYaxmZpahstNTEbFK0vHAbKADuCAi7pZ0XDF8OjALeBewEPgbcHRRfRRwhaRajBdHxG+qitXMzPJUeU2DiJhFSgzlftNLnwP4RIN69wG7VBmbmZm1zr8INzOzbE4aZmaWzUnDzMyyOWmYmVk2Jw0zM8vmpGFmZtmcNMzMLJuThpmZZXPSMDOzbE4aZmaWzUnDzMyyOWmYmVk2Jw0zM8vmpGFmZtmcNMzMLJuThpmZZXPSMDOzbE4aZmaWrdLXvUraH/gG6R3h50fEmXXDVQx/F+kd4R+OiNtz6poNlpl3LOWs2QtYtnwFozcbwUn77cjBu41pWPaI827ihnufWNO9xw6b86Nj37rW422l7Kkz53PJLYs5caeVfPSUWRw+ZRxnHLxzi60emFj7EsNh457hc2f+NmserI6gQxp286C/x9tMZUcakjqAc4ADgAnA4ZIm1BU7ABhf/E0DvtNCXbMBN/OOpZxy+XyWLl/aFGmJAAAJdklEQVRBAEuXr+CUy+cz846lLylbnzAAbrj3CY4476a1Gm8rZU+dOZ+Lbn6A1REArI7gopsf4NSZ8/vU/ipj7WsM9DLedWEe9Od4e1Pl6anJwMKIuC8ingcuBQ6qK3MQ8MNIbgY2k/TqzLpmA+6s2QtYsXJ1t34rVq7mrNkLXlK2PmE069/KeFspe8ktixvG0FP/HFXFWlUMngf9S1Fk334fsfQBYP+IOKboPhKYEhHHl8r8EjgzIv5QdF8LfAbYrre6pXFMIx2lMGrUqImXXnppJe3pq66uLjbZZJPBDqPfravtmr/0qR6H7TzmFW1ddtQIeHhFz2VzVRVrX2PorV1VxDAQ86A/18Vm9tprr9siYlJu+SqvaahBv/oM1VOZnLqpZ8QMYAbApEmTYurUqS2EWL3Ozk7aLab+sK6263Nn/nbNKZGyMZuN4JNHdK/34ZN/1eN4FtWVbWW8rZT96Cmz1pyW+fedV3H2/LTJd0jcW1c2V1Wx9jWGcrt6mwdl7T4P+nNd7E9Vnp5aAowrdY8FlmWWyalrNuBO2m9HRmzQ0a3fiA06OGm/HV9Sdo8dNm84jkb9WxlvK2UPnzLuJf2a9c9RVaxVxeB50L+qPNKYA4yXtD2wFDgM+GBdmSuB4yVdCkwBnoqIByU9mlHXbMDV7kzJuWPlR8e+NfvuqVbG20rZ2h1CtfP3/XHnUFWx9jUGeIYxmfOgv+6eard5MJB3TxERlf2RbqX9M3Av8Lmi33HAccVnke6SuheYD0xqVre3v4kTJ0a7ue666wY7hEq4XUPPcG2b27V2gLnRwvd6pb/TiIhZwKy6ftNLnwP4RG5dMzMbXP5FuJmZZXPSMDOzbE4aZmaWzUnDzMyyOWmYmVk2Jw0zM8vmpGFmZtmcNMzMLJuThpmZZXPSMDOzbE4aZmaWzUnDzMyyVfbmvsFQPFL9r4MdR50tgccGO4gKuF1Dz3Btm9u1draNiFflFh5WSaMdSZobLbxKcahwu4ae4do2t2tg+fSUmZllc9IwM7NsThrVmzHYAVTE7Rp6hmvb3K4B5GsaZmaWzUcaZmaWzUnDzMyyOWn0I0mLJM2XNE/S3KLfaZKWFv3mSXrXYMfZKkmbSbpM0p8k/VHSWyVtLulqSX8p/r9ysOPsix7aNqSXmaQdS7HPk/S0pBOH+jJr0q4hvbwAJP2rpLsl3SXpEkkbtevy8jWNfiRpETApIh4r9TsN6IqIrw5WXGtL0g+A30fE+ZI2BDYGPgs8ERFnSjoZeGVEfGZQA+2DHtp2IkN8mdVI6gCWAlOATzAMlhm8pF1HM4SXl6QxwB+ACRGxQtJPgFnABNpweflIw5qS9HJgT+C7ABHxfEQsBw4CflAU+wFw8OBE2HdN2jacvBO4NyL+yjBYZiXldg0H6wMjJK1P2nFZRpsuLyeN/hXAVZJukzSt1P94SXdKuqBdDjFb8BrgUeB7ku6QdL6kkcCoiHgQoPi/1WAG2Uc9tQ2G9jIrOwy4pPg8HJZZTbldMISXV0QsBb4KPAA8CDwVEVfRpsvLSaN/7RERuwMHAJ+QtCfwHWAHYFfSCnH2IMbXF+sDuwPfiYjdgGeBkwc3pH7TU9uG+jIDoDjd9h7gp4MdS39q0K4hvbyKJHcQsD0wGhgp6Z8HN6qeOWn0o4hYVvx/BLgCmBwRD0fE6oh4ATgPmDyYMfbBEmBJRNxSdF9G+qJ9WNKrAYr/jwxSfGujYduGwTKrOQC4PSIeLrqHwzKDunYNg+W1N3B/RDwaESuBy4F/oE2Xl5NGP5E0UtKmtc/AvsBdtYVeeC9w12DE11cR8RCwWNKORa93AvcAVwIfKvp9CPj5IIS3Vnpq21BfZiWH0/0UzpBfZoVu7RoGy+sB4C2SNpYk0nr4R9p0efnuqX4i6TWkowtIpz0ujogvSbqQdNgcwCLgY7XzlEOFpF2B84ENgftId6usB/wE2Ia00h8SEU8MWpB91EPb/pehv8w2BhYDr4mIp4p+WzDEl1kP7RoO29jpwKHAKuAO4BhgE9pweTlpmJlZNp+eMjOzbE4aZmaWzUnDzMyyOWmYmVk2Jw0zM8vmpGFDlqQtSk82fajuSacbDnZ8jUj6iKStKxz/SEmdktYruj8t6XZJHyi615e0uu5pseMk7Srp/KrisuFj/cEOwKyvIuJx0v35bfU0YUkdEbG6h8EfAW4HHmphfOtHxKrM4scAP42IFyS9AtiF9Avpy4o/gGciYte6eosl7SBpTPEsJLOGfKRhw5KkD0m6tdiT/rak9Yq97OWSzir2vmdLmiLpd5Luq72HQdIxkq4ohi+QdGrmeM+QdCswWdLpkuYovR9hupJDSUnux7WjIUlLJG1WjPstkq4pPp8h6VxJV5MeqLi+pK8V075T0jE9NP0IXvzlcKvb9y9JPzAz65GThg07knYiPU7iH4o96vVJT0UFeAVwVfFgyeeB00iPbTgE+EJpNJOLOrsDHyxO3/Q23tsjYnJE3AR8IyLeDOxcDNs/In4MzAMOjYhdI+L5XpqyG3BgRBwJTAMeiYjJwJtJD8Tcpq7dGwFjI2IJQEQ8CSwA5tD9cSKblk5NXVbqPxf4x15isnWcT0/ZcLQ36Yt1bnqUDyNIj54AWBERVxef55MeQ71K0nxgu9I4ZhdfukiaCbyNtL30NN7nefExMgDvlHQSsBGwJXAb8OsW2/HziPh78Xlf4A2SyklqPOnxEjVbAd0eMxERZwBn1I230ekpSA/EG91ijLaOcdKw4UjABRHx/7r1TC+4Ke/dvwA8V/pc3h7qn68TvYx3RRTP5Cmej/Qt0hNzl0o6g5Q8GlnFi0f89WWerWvTxyPi2h7GA7CiyXRybFSMw6xHPj1lw9E1wD9J2hLW3GW1TS916u2r9P7wjUnvOrihhfGOICWhx4onH7+/NOwZYNNS9yJgYvG5XK7ebODjRYKqvS97RLlARDwKbLQWd469jqH3hFgbYD7SsGEnIuYXTw29prj1dCVwHOkVmrn+AFxMernPhRExD9Y8jbTpeCPicaV3j98F/BW4pTT4e8D5klaQrpucBpwn6SHg1ibxnEt62um84tTYI6RkVu9a0rsYOltoa81etMnjt619+Sm3ZnWKO5N2iogTBzuWVkl6M+k01tEt1hsBXEd6+2RPtwub+fSU2XASEXOAP9R+3NeCbYD/cMKw3vhIw8zMsvlIw8zMsjlpmJlZNicNMzPL5qRhZmbZnDTMzCzb/wcByW/4iZuKdAAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "plt.scatter(data[\"Temperature\"], data[\"Frequency\"])\n",
    "plt.xlabel(\"Temperature (°F)\")\n",
    "plt.ylabel(\"Frequency of Failure\")\n",
    "plt.grid(True)\n",
    "plt.title(\"Relationship between temperature and failure frequency\")\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "This graph shows the relationship between flight temperature and failure frequency. Each point represents a flight and its associated temperature, with the failure frequency on the ordinate (vertical axis). This helps visualize whether there is a visible link between a drop in temperature and an increase in failures."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Although there is an apparent trend of increasing malfunctions at lower temperatures, this needs to be confirmed by a statistical model (logistic regression). The inclusion of all flights, including those without failures, is essential to obtain a complete picture of the relationship between temperature and malfunctions."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Estimation of the influence of temperature\n",
    "This part of the code uses a logistic regression model to estimate the impact of temperature on the probability of O-ring failure. The model takes \"temperature\" as the explanatory (or predictive) variable and \"failure frequency\" as the dependent (or response) variable."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Logistic regression on the entire data set"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 26,
   "metadata": {},
   "outputs": [],
   "source": [
    "data[\"Intercept\"] = 1\n",
    "logmodel = sm.GLM(data['Frequency'], data[['Intercept', 'Temperature']], family=sm.families.Binomial(sm.families.links.logit())).fit()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 27,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "                 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:                Mon, 21 Oct 2024   Deviance:                       3.0144\n",
      "Time:                        09:32:40   Pearson chi2:                     5.00\n",
      "No. Iterations:                     6   Covariance Type:             nonrobust\n",
      "===============================================================================\n",
      "                  coef    std err          z      P>|z|      [0.025      0.975]\n",
      "-------------------------------------------------------------------------------\n",
      "Intercept       5.0850      7.477      0.680      0.496      -9.570      19.740\n",
      "Temperature    -0.1156      0.115     -1.004      0.316      -0.341       0.110\n",
      "===============================================================================\n"
     ]
    }
   ],
   "source": [
    "print(logmodel.summary())"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "*Interpretation :*\n",
    "\n",
    "* Intercept Coefficient: 5.0850\n",
    "The intercept represents the predicted value of the failure frequency when the temperature is at 0°F. This is not very physically relevant, since temperatures in this setting are never that low, but the intercept is part of the model to fit the logistic curve.\n",
    "\n",
    "* Temperature Coefficient: -0.1156\n",
    "This negative coefficient indicates that as the temperature increases, the probability of failure decreases slightly. However, the effect is small and not significant.\n",
    "\n",
    "* Temperature P-value: 0.316\n",
    "With a p-value of 0.316 (greater than 0.05), this means that the effect of temperature on failure frequency is not statistically significant in this model. In other words, the analysis does not allow us to conclude that temperature has a significant effect on O-ring failures, probably due to the small sample size or the variability of the data.\n",
    "\n",
    "* Temperature confidence interval: [-0.341, 0.110]\n",
    "This interval contains 0, which confirms that the effect of temperature could be zero.\n",
    "\n",
    "*Conclusion :*\n",
    "The logistic regression model shows a trend suggesting that higher temperatures may reduce the probability of failure, but this relationship is not statistically significant in this case. It would be necessary to explore other variables or increase the sample size to obtain a more precise analysis."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Estimation of the probability of O-ring malfunction"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Model predictions for a range of temperatures\n",
    "Here, predictions are generated for a temperature range from 30°F to 90°F. The model predicts the failure frequency for each of these temperatures, allowing us to see the predicted failure probability for extreme temperatures such as 31°F (the temperature on the day of the disaster)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 28,
   "metadata": {},
   "outputs": [],
   "source": [
    "data_pred = pd.DataFrame({'Temperature': np.linspace(start=30, stop=90, num=121), 'Intercept': 1})\n",
    "data_pred['Frequency'] = logmodel.predict(data_pred[['Intercept', 'Temperature']])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "An increase in failure probability is expected for low temperatures, if the model is well fitted."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Visualization of predictions with confidence interval"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 29,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYUAAAEWCAYAAACJ0YulAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMi4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvIxREBQAAIABJREFUeJzs3Xd8FWXWwPHfyU1CSKFDEKlC6CUkgKCgIIpYgBVFVERlV1Gwrb2w+6rviq+ufRVQLKgLiij2CkqxINJEioB0iPSahBDSzvvHTMIFU27Kzb1Jzvfzmc+9M3fKeXKTOZnnmXkeUVWMMcYYgJBAB2CMMSZ4WFIwxhiTx5KCMcaYPJYUjDHG5LGkYIwxJo8lBWOMMXksKZhyJSIjRGRWCbddLSJ9yzikoCQiqSJyWqDjKA4ReUNEHvVx3S0icq6/YzLFZ0nBFMgff7iqOk1VB/hw7D+dYFS1g6rOK8t4gpWqRqvqpkDHYaoeSwqmShGR0GDenzGBZknBlIiI3CAiG0TkgIh8IiKNvD4bICLrROSwiEwUkfkicr372XUi8oP7XkTkWRHZ4667QkQ6ishoYARwr1uN8qm7ft6Vi4h4RORBEdkoIikislREmuQTZ3MRURH5m4hsA+a4y3uKyAIROSQiv3pXS4lICxH5zt3vNyIyQUSmlmJ/14nIJnd/m0VkhLu8lfuzOSwi+0TkXa9tVERaue9rishbIrJXRLaKyD9EJMT75ykiT4nIQXf/FxTyvW0RkXvcn/UREXlNRGJF5Euv8tb2Wn+wW213SETmiUg7r8+6isgyd7t3gYiTjnWxiCx3t10gIp0LissEEVW1yaZ8J2ALcG4+y88B9gEJQDXgBeA797N6QDIwFAgFbgcygevdz68DfnDfnw8sBWoBArQDTnE/ewN4tKB4gHuAlUAbd9suQN18Ym0OKPAWEAVUB04F9gMX4vxjdJ47X9/d5ifgKSAc6O2WZ2pJ9ueukwy0cbc/Bejgvn8HGOduEwH09opbgVbu+7eAj4EY9/i/A3/z+nlmAjcAHmAMsAOQQr7ThUCsG/ceYBnQ1f0u5wAPueu2Bo645QkD7gU2uD+XcGArcIf72WVuHI+62ya4+z7djeta99jVCvvdsinwk10pmJIYAbyuqstU9RjwANBLRJrjnBhXq+oHqpoF/AfYVcB+MnFOdG1xTmJrVHWnjzFcD/xDVdep41dV3V/I+g+r6hFVPQpcDXyhql+oao6qzgaWABeKSFOgO/A/qpqhqj8An5R0f+66OUBHEamuqjtVdbVX+ZsBjVQ13T3WCUTEAwwHHlDVFFXdAjwNjPRabauqvqKq2cCbOIkntpCfxQuqultV/wC+B35W1V/c7/JDnASBe9zPVXW2qmbiJMrqwBlAT5xk8JyqZqrq+8Bir2PcALysqj+raraqvgkcc7czQcySgimJRjj/JQKgqqk4/xmf6n623eszBZLy24mqzgFeBCYAu0VksojU8DGGJsDGYsS83et9M2CYW61xSEQO4VwRnOLGf0BV0wrYtlj7U9UjOCfXm4CdIvK5iLR1t7sX5ypnkVtF89d8jlOP4/+V59qK87POlZd0veKOzmdfuXZ7vT+az3zutid/zzk45c79nv9wv1/vuHI1A+466WfSxN3OBDFLCqYkduD80QMgIlFAXeAPYCfQ2Osz8Z4/mar+R1UTgQ441RX35H5URAzbgZbFiNl7f9uB/6pqLa8pSlUfd+OvIyKRXuv/qa2iGPtDVb9W1fNwks5a4BV3+S5VvUFVGwE3AhNz2xG87OP4FUWupjg/a387+XsWnJ9F7vd8qrvMO65c24HxJ/1MIlX1nXKI25SCJQVTlDARifCaQoG3gVEiEi8i1YDHcKogtgCfA51E5C/uujcDDfPbsYh0F5HTRSQMp+46Hch2P94NFHaf/qvAv0QkThydRaSuj2WaCgwSkfPFabCOEJG+ItJYVbfiVP08LCLhItILGFTS/bmNuIPdxHkMSM0to4gME5HchHkQJ9Fke+/YrRKaAYwXkRgRaQbc6R7T32YAF4lIf/c7usstwwKcdpcs4DYRCRWRoUAPr21fAW5yv18RkSgRuUhEYsohblMKlhRMUb7AqVLInR5W1W+BfwIzcf5jbAlcAaCq+4BhwL9xqpTa45xkj+Wz7xo4J4+DOFUP+3HqrQFeA9q7VQ8f5bPtMzgnrVk4Dbmv4dR3F0lVtwNDgAeBvTj/1d7D8b+HEUAvN55HgXcLiN+X/YXgnEx3AAeAs4Gx7qbdgZ9FJBWn3eJ2Vd2czyFuxUmam4AfcJLy676UtTRUdR1Oe8kLOFcsg4BBbltLBs7NBNfhfH/DgQ+8tl2C067wovv5BnddE+TkxCpBY8qWe+tkEjBCVecGOp6ScG+3XKuqDwU6FmP8za4UTJlzq1FquVVLD+I0pi4McFg+c6u1WopIiIgMxLkKyO9qxZhKx57GNP7QC6eKIxz4DfiLe+tmRdEQpyqkLs5VzhhV/SWwIRlTPqz6yBhjTB6rPjLGGJOnwlUf1atXT5s3b16ibY8cOUJUVFTZBhQgVpbgVFnKUlnKAVaWXEuXLt2nqvWLWq/CJYXmzZuzZMmSEm07b948+vbtW7YBBYiVJThVlrJUlnKAlSWXiGwtei2rPjLGGOPFkoIxxpg8lhSMMcbkqXBtCsaYyiczM5OkpCTS09PLfN81a9ZkzZo1Zb7fQPClLBERETRu3JiwsLASHcOSgjEm4JKSkoiJiaF58+ac2PFq6aWkpBATUzn64SuqLKrK/v37SUpKokWLFiU6hlUfGWMCLj09nbp165Z5QqhqRIS6deuW6orLkoIxJihYQigbpf05VpmkcOBIBm+vOcaRY1mBDsUYY4JWlUkKP2zYx+ytWQx+8QfW7koOdDjGmCDj8XiIj4+nY8eODBs2jLS0tKI3KsC8efO4+OKLAfjkk094/PHHC1z30KFDTJw4sdjHePjhh3nqqaeKXrGYqkxSGNylEfd2jyA5PYshL/7IjMX5DbtrjKmqqlevzvLly1m1ahXh4eG89NJLJ3yuquTk5BR7v4MHD+b+++8v8POSJgV/qTJJAaBdXQ9f3NaHbs1rc+/MFdw/cwXpmdlFb2iMqVL69OnDhg0b2LJlC+3atWPs2LEkJCSwfft2Zs2aRa9evUhISGDYsGGkpqYC8NVXX9G2bVt69+7NBx/kDULHG2+8wS233ALA7t27ueSSS+jSpQtdunRhwYIF3H///WzcuJH4+HjuuccZovzJJ5+ke/fudO7cmYceOj6205NPPkmbNm0499xzWbdunV/KXuVuSa0fU423/no6z8xex4S5G1m9I5mXRiZyai2fRnI0xvjZI5+u5rcdZVfFm52dTacmtXloUAef1s/KyuLLL79k4MCBAKxbt44pU6YwceJE9u3bx6OPPso333xDVFQUTzzxBM888wz33nsvN9xwA3PmzKFVq1YMHz48333fdtttnH322Xz44YdkZ2eTmprK448/zqpVq1i+fDkAs2bNYv369SxatAhVZfDgwXz33XdERUUxc+ZMfvnlF7KyskhISCAxMbFsfkheqlxSAPCECPec35b4JrW5893lDH7hByaMSKDnab6O+26MqWyOHj1KfHw84Fwp/O1vf2PHjh00a9aMnj17ArBw4UJ+++03zjzzTAAyMjLo1asXa9eupUWLFsTFxQFw9dVXM3ny5D8dY86cObz11luA04ZRs2ZNDh48eMI6s2bNYtasWXTt2hWA1NRU1q9fT0pKChdffDGRkZGAUy3lD35NCu5Qhs8DHuBVVX38pM9rAlOBpm4sT6nqFH/G5O289rF8ePOZjP7vEq5+9WceGtSekb2al9fhjTH58PU/el/5+vBabpvCyby7qlZVzjvvPN55550T1lm+fHmZ3VKrqjzwwAPceOONJyx/7rnnyuW2Xb+1KYiIB5gAXAC0B64UkfYnrXYz8JuqdgH6Ak+LSLi/YspPqwbRfHTzmZzVuj7//Hg1//xoFZnZxW9MMsZUfj179uTHH39kw4YNAKSlpfH777/Ttm1bNm/ezMaNGwH+lDRy9e/fn0mTJgFOtVZycjIxMTGkpKTkrXP++efz+uuv57VV/PHHH+zZs4ezzjqLzz77jKNHj5KSksKnn37qlzL6s6G5B7BBVTepagYwHWcAdG8KxIiT/qKBA0C5P0hQIyKMV67pxo1nncZ/F25l1JTFHD6aWd5hGGOCXP369XnjjTe48sor6dy5Mz179mTt2rVEREQwefJkLrroInr37k2zZs3y3f75559n7ty5dOrUicTERFavXk3dunU588wz6dixI/fccw8DBgzgqquuolevXnTq1InLLruMlJQUEhISGDp0KPHx8Vx66aX06dPHL2X02xjNInIZMFBVr3fnRwKnq+otXuvEAJ8AbYEYYLiqfp7PvkYDowFiY2MTp0+fXqKYUlNTiY6OLnSd75MyeWN1BrGRwh2JEdSPDM4btHwpS0VhZQk+5V2OmjVr0qpVK7/sOzs7G4/H45d9lzdfy7JhwwYOHz58wrJ+/fotVdVuRW3rzzaF/Cq/Ts5A5wPLgXOAlsBsEfleVU+49UBVJwOTAbp166YlHXnIl1GL+gLnbtzPjf9dwhPLcnjt2gS6NKlVouP5k40mFZwqS1nKuxxr1qzxW6d1ValDvFwRERF5DdXF5c9/g5OAJl7zjYEdJ60zCvhAHRuAzThXDQHVq2VdPhh7BhFhIVwxeSFz1+4JdEjGGFMu/JkUFgNxItLCbTy+AqeqyNs2oD+AiMQCbYBNfozJZ60axPDh2DNp2SCK699awowl9gS0Maby81tSUNUs4Bbga2ANMENVV4vITSJyk7vav4AzRGQl8C1wn6ru81dMxVU/phrTR/fijJZ1uff9FUyctwF/tcEYY0ww8OtzCqr6BfDFScte8nq/AxjgzxhKK7paKK9d25173v+Vf3+1jkNpmTxwQVvr5tcYUylVySeaiys8NIRnL4+nVvUwJn+3iYNHMnj80s54QiwxGGMql+C83zIIhYQIDw/uwG3943hvaRK3Tf/FHnIzppJJSkpiyJAhxMXF0bJlS26//XYyMjJO6NQumPjjtmFLCsUgItx5XmsevLAtn6/YyZipS62XVWMqCVVl6NCh/OUvf2H9+vX8/vvvpKamMm7cOL8cLysrOAf8sqRQAqPPasm//tKRb9bs4Ya3llhiMKa8TZsGzZtDSIjzOm1aqXc5Z84cIiIiGDVqFOB0WPfss8/y+uuvk5aWxvbt2xk4cCBt2rThkUceAeDIkSNcdNFFdOnShY4dO/Luu+8CsHTpUs4++2wSExM5//zz2blzJwB9+/blwQcf5Oyzz2b8+PE0b948b4yGtLQ0mjRpQmZmJhs3bmTgwIEkJibSp08f1q5dC8CWLVvo1asX3bt355///Gepy5wfa1MooZE9m1EtNIT7Zq7gb28u5tVrulM9vHI8NWlMUJs2DUaPhtyR0bZudeYBRowo8W5Xr179p66oa9SoQdOmTcnKymLRokWsWrWKyMhIunfvzkUXXcTWrVtp1KgRn3/udMRw+PBhMjMzufXWW/n444+pX78+7777LuPGjeP1118HnEF15s+fD8CyZcuYP38+/fr149NPP+X8888nLCyM0aNH89JLLxEXF8fPP//M2LFjmTNnDvfddx9jxozhmmuuYcKECSUua2HsSqEULu/WhKeHdeGnjfsZ9cYi0jKC83LQmEpl3LjjCSFXWpqzvBRUNd+7CnOXn3feedStW5fq1aszdOhQfvjhBzp16sQ333zDfffdx/fff0/NmjVZt24dq1at4rzzziM+Pp5HH32UpKSkvP15j7UwfPjwvKuL6dOnM3z4cFJTU1mwYAHDhg0jPj6eG2+8Me9KY+HChVx55ZUAjBw5slTlLYhdKZTS0ITGeEKEO95dzvVvLuG1a+2KwRi/2rateMt91KFDB2bOnHnCsuTkZLZv347H4/lTwhARWrduzdKlS/niiy944IEHGDBgAJdccgkdOnTgp59+yvc43l1xDx48mAceeIADBw6wdOlSzjnnHI4cOUKtWrXy7cY797j+ZFcKZWBI/Kk8fXkXftq039oYjPG3pk2Lt9xH/fv3Jy0tLW8QnOzsbO666y6uu+46IiMjmT17NgcOHODo0aN89NFHnHnmmezYsYPIyEiuvvpq7r77bpYtW0abNm3Yu3dvXlLIzMxk9erV+R4zOjqaHj16cPvtt3PxxRfj8XioUaMGLVq04L333gOcK5Vff/0VcLruzu0QdFoZtKPkx5JCGbmka2OeuqwLP27cx+j/LuVYliUGY/xi/HhwRx/LExnpLC8FEeHDDz/kvffeIy4ujtatWxMREcFjjz0GQO/evRk5cmRe19XdunVj5cqV9OjRg/j4eMaPH88//vEPwsPDef/997nvvvvo0qUL8fHxLFiwoMDjDh8+nKlTp55QrTRt2jRee+01unTpQocOHfj4448BeOKJJ5gwYQLdu3f/Uy+oZUZVK9SUmJioJTV37twSb+ur6Yu2arP7PtMb3lysGVnZfjtOeZSlvFhZgk95l+O3334r3gZTp6o2a6Yq4rxOnVrgqsnJyaWKLZj4Wpb8fp7AEvXhHGttCmVsePemHM3I5uFPf+POGb/y3PB4e/LZmLI2YkSp7jQyBSsyKYhIJHAX0FRVbxCROKCNqn7m9+gqqOvObEF6Vg6Pf7mW6GoeHrukk/WVZIypEHy5UpgCLAV6ufNJwHuAJYVC3HR2S1LSM5kwdyM1IsK43zrRM6ZQWsAtoaZ4tJQ9OfuSFFqq6nARudI94FGxb84ndw9oQ/LRLF7+bhM1qodxcz//DDdoTEUXERHB/v37qVu3riWGUlBV9u/fT0RERIn34UtSyBCR6rhDaYpIS+BYiY9YhYgIjwzuQHJ6Jk9+vY66UeFc0aN0t80ZUxk1btyYpKQk9u7dW+b7Tk9PL9VJMpj4UpaIiAgaN25c4mP4khQeAr4CmojINOBM4LoSH7GKCQkRnrysC4fSMnnww5XUjgrn/A4NAx2WMUElLCyMFi1a+GXf8+bNK/F4xcGmPMpS6HMKbjXRWmAoTiJ4B+imqvP8GlUlEx4awqSrE+jcuBa3vvMLizYfCHRIxhiTr0KTgntv60equl9VP1fVzzSIhsusSCLDQ5lyXXca167O9W8uZv3ulECHZIwxf+LLE80LRaS73yOpAmpHhfPmqB5UC/Nw7euL2HU4PdAhGWPMCXxJCv2An0Rko4isEJGVIrLC34FVVk3qRPLGqO4kp2dx3ZRFpKRnBjokY4zJ40tSuABoCZwDDAIudl9NCXVoVJNJVyewYU8qY6cts2E9jTFBw5ekoAVMphT6xNXnsaGd+H79Ph78YGWpHzgxxpiy4MstqZ/jJAEBIoAWwDqggx/jqhIu79aEpANp/GfOBprXi7KH24wxAVdkUlDVTt7zIpIA3Oi3iKqYO85rzdYDaTz59Tqa1olkUJdGgQ7JGFOFFXs8BVVdBtjdSGVERHji0s50a1abu977laVbDwY6JGNMFVZkUhCRO72mu0XkbaDsn0WvwiLCPEy+phun1Ixg9FtLSDqYVvRGxhjjB75cKcR4TdVw2hiG+DOoqqhOVDivXdudjOwcrn9zCanHsgIdkjGmCvKlTeGR8gjEQKsG0UwckcB1UxZz+zu/MPmabjZAjzGmXBWYFETkUwq59VRVB/sloiquT1x9Hh7Unn9+vJp/f72WBy5oF+iQjDFVSGFXCk+VWxTmBCN7NWfd7hRenr+JNrExDE0oeTe4xhhTHAUmBVWdX56BmBM9NKgDG/cc4f4PVtKiXhRdm9YOdEjGmCqgwIZmEZnhvq50+zw6YSq/EKumME8IE0ckEFujGqP/u5TdydZ5njHG/wq7++h29zW3r6OTJ+NntaPCefWa7hw5lsWN/11KemZ2oEMyxlRyBSYFVd3pvm7Nbyq/EKu2Ng1jeObyeJZvP8S4D1dZH0nGGL/y5eG1niKyWERSRSRDRLJFJLk8gjOOgR0bclv/OGYuS+KNBVsCHY4xphLzpUO8F4ErgPeAbsA1gPXcVs7+3j+O33Yk8+jna2h3So1Ah2OMqaR86vtIVTcAHlXNVtUpOAPvmHIUEiI8M7wLzepGcvO0Zew/amMwGGPKni9JIU1EwoHlIvJvEbkDiPJzXCYfNSLCmDyyG8eycnjxl2PW8GyMKXO+JIWR7nq3AEeAJsClvuxcRAaKyDoR2SAi9xewTl8RWS4iq0XEno0oQqsG0TxzeRc2J+fw8CerAx2OMaaSKaybi6aqus3rTqN0wOd+kETEA0wAzgOSgMUi8omq/ua1Ti1gIjBQVbeJSIOSFKKqGdChIYNOC2P64u10aVKLK3s0DXRIxphKorArhY9y34jIzBLsuwewQVU3qWoGMJ0/9656FfCBqm4DUNU9JThOlXRJXBhnta7PQx+v5pdtNgaDMaZsSEH3vYvIL6ra9eT3Pu9Y5DKcK4Dr3fmRwOmqeovXOs8BYThDe8YAz6vqW/nsazQwGiA2NjZx+vTpxQklT2pqKtHR0SXaNtikpqZCeBQPLTiKAg+fUZ0a4RWzR9XK9r1UhrJUlnKAlSVXv379lqpqt6LWK+yWVC3gva/yO0OdvJ9QIBHoD1QHfhKRhar6+wkbqU4GJgN069ZN+/btW4JwYN68eZR022CTW5Zm7Q9z6UsLeHdbdd766+kVsqvtyvi9VHSVpRxgZSmuwqqPuohIsoikAJ3d98kikuLjw2tJOI3SuRoDO/JZ5ytVPaKq+4DvgC7FKUBV16lxTf41pAM/btjP07PWBTocY0wFV1g3Fx5VraGqMaoa6r7Pnffl6anFQJyItHBvab0C+OSkdT4G+ohIqIhEAqcDa0pamKpqePemXNG9CRPnbWT2b7sDHY4xpgLz6eG1klDVLJzbWL/GOdHPUNXVInKTiNzkrrMG+ApYASwCXlXVVf6KqTJ7eHAHOp5agztnLGfbfhvj2RhTMn5LCgCq+oWqtlbVlqo63l32kqq+5LXOk6raXlU7qupz/oynMosI8zBpRCICjJlmPaoaY0qmsPEUqpVnIKb0mtSJ5Lkr4lm9I9kebDPGlEhhVwo/AYjIf8spFlMGzmkby839WjJ98XbeX5oU6HCMMRVMYbekhovItcAZIjL05A9V9QP/hWVK445zW7Ns6yH+8dFKOp5ag7YNrVdVY4xvCrtSuAnoCdTiz6OuXez/0ExJhXpCeP7KeGpEhDF26jJS0jMDHZIxpoIo8EpBVX8AfhCRJar6WjnGZMpAg5gIXriyK1e9+jMPfLCSF67sikjFe7DNGFO+fLn76L8icpuIvO9Ot4pImN8jM6V2+ml1uXtAGz5bsZP/LrQRVI0xRfMlKUzE6YpiojslAJP8GZQpOzeedRr92zbgX5/9xq/bDwU6HGNMkPMlKXRX1WtVdY47jQK6+zswUzZCQoSnL+9Cg5gIbn57GYfTrH3BGFMwX5JCtoi0zJ0RkdMAezKqAqkVGc6LV3Vld3I6d7//KwX1jGuMMb4khXuAuSIyzx0ZbQ5wl3/DMmWta9PaPHBBO2b/tpvXftgc6HCMMUGqsOcUAFDVb0UkDmiD0x32WlU95vfITJkbdWZzFm0+wONfrqVr09okNqsd6JCMMUHGp76PVPWYqq5Q1V8tIVRcIsK/h3WmUa3q3Pr2Mg4eyQh0SMaYIOPXDvFM8KkREcaEqxLYl5rBnTOWk5Nj7QvGmOMsKVRBnRrX5B8Xt2Puur28/N2mQIdjjAkiRSYFEZkpIheJiCWQSmRkz2Zc1OkUnpq1jsVbDgQ6HGNMkPDlRD8JuApYLyKPi0hbP8dkyoGI8PilnWhcuzq3vv0L+1OtqcgY40NSUNVvVHUEzpPMW4DZIrJAREZZdxcVW4zbvnAgLYM7Z/xq7QvGGN/aFESkLnAdcD3wC/A8TpKY7bfITLnoeGpN/ufi9sz/fS+T5m8MdDjGmADzpU3hA+B7IBIYpKqDVfVdVb0ViPZ3gMb/RpzelIs7n8LTs9axaLO1LxhTlflypfCiO4by/6nqTu8PVLWbn+Iy5UhE+L+hnWhWN4pb31nGPmtfMKbK8iUptBORWrkzIlJbRMb6MSYTALntC4fSMrnjXXt+wZiqypekcIOq5vW5rKoHgRv8F5IJlPaNavDw4A58v34fE+ZuCHQ4xpgA8CUphIjXkF0i4gHC/ReSCaQrujdhSHwjnv3mdxZs3BfocIwx5cyXpPA1MENE+ovIOcA7wFf+DcsEiojw2CWdaF4vitunL2dPSnqgQzLGlCNfksJ9ON1ljwFuBr4F7vVnUCawoqqFMnFEAinpmfx9+nKyrX3BmCrDl4fXclR1kqpepqqXqurLqmqD7FRybRvW4H+HdGTBxv08/+36QIdjjCknvjyncKaIzBaR30Vkk4hsFhHrRa0KuLxbEy5NaMwLc9bz/fq9gQ7HGFMOfKk+eg14BuiNMzZzN2yM5irjX3/pQFyDaG6fvpxdh619wZjKzpekcFhVv1TVPaq6P3fye2QmKESGO+0L6ZnZ3PrOMjKzcwIdkjHGj3xJCnNF5EkR6SUiCbmT3yMzQaNVgxj+b2gnFm85yFNfrwt0OMYYPypyjGbgdPfVu0sLBc4p+3BMsBoSfyqLNh/g5e820a15Hc5rHxvokIwxflBkUlDVfuURiAl+/7y4Pb8mHeKuGcv57NY+NK0bGeiQjDFlzJe7j2JF5DUR+dKdby8if/N/aCbYRIR5mDQiEYCxby8lPdPuTDamsvGlTeENnKeaG7nzvwN/91dAJrg1qRPJ05fHs+qPZB759LdAh2OMKWO+JIV6qjoDyAFQ1SzA/kWsws5rH8tNZ7fknUXbmLk0KdDhGGPKkC9J4Yg78poCiEhP4LBfozJB7+4BrTm9RR3GfbSStbuSAx2OMaaM+JIU7gQ+AVqKyI/AW8Ctfo3KBL1QTwgvXNWVGhFhjJm6jOT0zECHZIwpA770fbQMOBs4A7gR6KCqK3zZuYgMFJF1IrJBRO4vZL3uIpItIpf5GrgJvAYxEUwYkcC2A2nc896vqFrHecZUdL7cfXQNcBWQCCQAV7rLitrOA0wALgBaRQd5AAAgAElEQVTau9u1L2C9J3Aas00F0715HR64oC1fr97Ny99Zl1jGVHS+PLzm3c9RBNAfWIZTjVSYHsAGVd0EICLTgSHAybes3ArMxPpTqrD+1rsFv2w/xL+/WkvnU2tyRqt6gQ7JGFNCvlQf3eo13QB0xbeR104FtnvNJ7nL8ojIqcAlwEu+h2yCjYjw70s7c1r9aG595xd2Hj4a6JCMMSUkxa0HFpEwYIWqtitivWHA+ap6vTs/Euihqrd6rfMe8LSqLhSRN4DPVPX9fPY1GhgNEBsbmzh9+vRixZwrNTWV6OjoEm0bbIKxLDtSc/jfn47SKDqEB06PICxEit6I4CxLSVWWslSWcoCVJVe/fv2Wqmq3otYrsvpIRD7FvR0V58qiPTDDhxiSgCZe842BHSet0w2Y7g4BXQ+4UESyVPUj75VUdTIwGaBbt27at29fHw7/Z/PmzaOk2wabYC1Lnea7uGnqUuYcqsf/De3k0zbBWpaSqCxlqSzlACtLcfnSpvCU1/ssYKuq+vLE0mIgTkRaAH8AV+A0WOdR1Ra5772uFE5ICKZiGdixIWP7tmTivI10blyTK3s0DXRIxphi8KVDvPkl2bGqZonILTh3FXmA11V1tYjc5H5u7QiV1F0D2rDyj8M89PFq2jSMIaFp7UCHZIzxkS+3pKaISHI+U4qIFPooq6p+oaqtVbWlqo53l72UX0JQ1evya08wFY8nRHjhyq7E1qzGmKlL2ZNiI7YZU1H48kTzs8D9OHcONQbuAx5V1RhVreHP4EzFVSsynMkju5F8NIuxU5eRkWUjthlTEfiSFM5X1YmqmqKqyao6CbjU34GZiq/dKTX492WdWbL1IA9/ujrQ4RhjfOBLUsgWkREi4hGREBEZgfWSanw0qEsjxvRtyds/b2Pqwq2BDscYUwRfksJVwOXAbncaxkl3ERlTmLsHtKFvm/o8/MlqFm0+EOhwjDGF8OWJ5i2qOkRV66lqfVX9i6puKYfYTCXhCRGev6IrTetEMmbqUpIOpgU6JGNMAXy5+6i1iHwrIqvc+c4i8g//h2Yqk5rVw3jl2m5kZOdw/ZtLOHIsK9AhGWPy4Uv10SvAA0AmgNtt9hX+DMpUTi3rR/PiVQn8vjuFO2csJyfHuto2Jtj4khQiVXXRScvs3zxTIme3rs+DF7bj69W7eXr2upLvaNo0aN4cQkKc12nTyipEUx7s+wtavnRzsU9EWnJ8OM7LgJ1+jcpUan/r3YINe1KZMHcjrRpEU+znnadNg9GjIc1tm9i61ZkHGDGiLEM1/mDfX1Dz5UrhZuBloK2I/AH8HbjJr1GZSk1E+N8hHel5Wh3ue38l6w8W8w7nceOOn1BypaU5y03ws+8vqBWaFEQkBOimqucC9YG2qtpbVe2Gc1Mq4aEhvHR1Io1qRfCfZels21+MO5K2bSvechNc7PsLaoUmBVXNAW5x3x9R1ZRyicpUCbUiw3n9uu7kAKPeWMThtEzfNmxaQM+rBS03wcW+v6DmS/XRbBG5W0SaiEid3MnvkZkq4bT60dzaNYJtB9IYM22pb30kjR8PkZEnLouMdJab4GffX1DzJSn8Fadd4TtgqTst8WdQpmppW8fD40M7s2Djfh74YCVFjgY4YgRMngzNmoGI8zp5sjVSVhT2/QU1X8ZTaFHUOsaU1qWJjdl+MI3nvllPkzrV+fu5rQvfYMQIO4lUZPb9Ba0CrxRE5DGv9+eVTzimKru9fxyXJTbmuW/W896S7YEOx5gqqbDqo4Fe75/wdyDGiAiPXdKJ3q3qcf8HK5m3bk+gQzKmyvGlTcGYchMeGsKkqxNoExvD2GnLWJF0KNAhGVOlFJYUGojInSJyl9f7vKm8AjRVT0xEGG+M6k6dqHBGTVnMln1HAh2SMVVGYUnhFSAGiPZ67z0Z4zcNakTw5l97kKPKNa8vYk+yjfNsTHko8O4jVX2kPAMx5mQt60czZVQPrnplIddOWcz00T2pWT0s0GEZU6lZm4IJavFNavHyyEQ27EnhhjeXcDTDRoI1xp8sKZig1yeuPs9cHs/irQd8f+rZGFMilhRMhTCoSyMeu6QT89bt5c4Zy8m2AXqM8YsC2xSKusNIVZ8p+3CMKdiVPZqSkp7JY1+sJTLc6RojJEQCHZYxlUph3Vzk3mHUBugOfOLOD8LpB8mYcjf6rJakHsvmP9+uJyLMwyODOyBiicGYslLk3UciMgtIyO02W0QeBt4rl+iMyccd58ZxLDObl7/bRLXQEB68sJ0lBmPKiC/DcTYFMrzmM4DmfonGGB+ICPdf0Jb0zGxe+X4znpAQ7hvYxhKDMWXAl6TwX2CRiHyIM07zJcBbfo3KmCKICA8P7kBWjvLS/I14QuDuAZYYjCktX7rOHi8iXwJ93EWjVPUX/4ZlTNFEhH8N6UiOKhPmbgQsMRhTWr5cKQBEAsmqOkVE6otIC1Xd7M/AjPFFSIgw/i+dAJgwdyNZOcr9A9taYjCmhIpMCiLyENAN5y6kKUAYMBU407+hGeOb3MTgCRFenr+JrGzlHxdZ47MxJeHLlcIlQFdgGYCq7hAR6xDPBJWQEKcqKTQkhNd+2Ex6Zjb/GtLRnmMwpph8SQoZqqoiogAiEuXnmIwpERHhoUHtqR7uYdK8jRzNyObfl3Um1GMP7hvjK1+SwgwReRmoJSI3AH8FXvVvWMaUjIhw38C2RFcL5cmv15FyLIsXruxKRJgn0KEZUyEU+S+Uqj4FvA/MxGlX+B9V/Y+/AzOmNG7u14qHB7Vn9m+7uW7KIlLSMwMdkjEVQpFJQUSeUNXZqnqPqt6tqrNFxMZsNkHvujNb8NzweBZvOciVryxkb8qxQIdkTNDzpbL1vHyWXVDWgRjjD3/peiqvXJPIhj2pXDppgQ3taUwRCkwKIjJGRFYCbUVkhde0GVjpy85FZKCIrBORDSJyfz6fj/Da7wIR6VLyohiTv3PaxvL2DT1JSc/k0kkLWL79UKBDMiZoFXal8DZOj6gfu6+5U6KqjihqxyLiASbgXFW0B64UkfYnrbYZOFtVOwP/AiYXuwTG+CChaW1mjjmDyGoerpj8E1+v3hXokIwJSgUmBVU9rKpbgOeBA6q6VVW3ApkicroP++4BbFDVTaqaAUwHhpx0jAWqetCdXQg0LkkhjPHFafWj+WDMmbRpWIObpi7l1e83oWqD9RjjTYr6oxCRX3C6zs59TiEEWKKqCUVsdxkwUFWvd+dHAqer6i0FrH830DZ3/ZM+Gw2MBoiNjU2cPn16kQXLT2pqKtHR0SXaNthYWUruWLYyecUxlu7Opm/jUK5uH05oGT3kVlm+l8pSDrCy5OrXr99SVe1W5IqqWugELM9n2QofthsGvOo1PxJ4oYB1+wFrgLpF7TcxMVFLau7cuSXeNthYWUonOztHH/9yjTa77zO94uWf9OCRY2Wy38ryvVSWcqhaWXLh/DNf5Dnfl7uPNonIbSIS5k63A5t82C4JaOI13xjYcfJKItIZ52G4Iaq634f9GlNqISHOQ27PXN6FpVsPMvjFH1m7KznQYRkTcL4khZuAM4A/cE70p+NW5RRhMRAnIi1EJBy4guNDegIgIk2BD4CRqvp7cQI3piwMTWjM9Bt7kp6ZzdCJC/hy5c5Ah2RMQPnyRPMeVb1CVRuoaqyqXqWqe3zYLgu4Bfgap2pohqquFpGbROQmd7X/AeoCE0VkuYgsKUVZjCmRhKa1+fTW3rRpGMOYact4/Mu1ZGXnBDosYwKisOcU7nVfXxCR/5w8+bJzVf1CVVuraktVHe8ue0lVX3LfX6+qtVU13p2KbgQxxg9iP53J9Keu5arlX/LS/I2MfOwT9r35NjRvDiEhzuu0af45+LRp5XMcX4wdC6GhsHSp8zp2bOBiMQFRWId4a9xX++/dVG7TpsHo0VRLS+OxzRNI+GMN486/mQsXZ/OfnBh6qsLWrTDarTUdUeRjOsU+Nmlpzry/juOLsWNh0qTj89nZx+cnTizfWEzAFPacwqfu65v5TeUXojF+Nm7c8ZMycNmqOXz01l1EH0vjqivG858zriBbQpx1xo3z67EB/xzHF5MLeHa0oOWmUirwSkFEPgUKfIhBVQf7JSJjytu2bX9a1G7vFj556w7GDbiZZ/pczYKmnXn286c5JZ91y/rYhS73p+zs4i03lVJhDc1PAU/jdEVxFHjFnVKBVf4PzZhy0rRpvoujM47y3GdP8eTnz7LilDgGjnqRr84YVC7HLnC5P3kKGHOioOWmUiqs+mi+qs4HuqrqcFX91J2uAnqXX4jG+Nn48RAZeeKysDAID0eAYau+5fM3bqdZ8m5u6j2aO2csJ7msxmfI79iRkc7y8ja6gDvNC1puKiVfnlOoLyKn5c6ISAugvv9CMqacjRjh1Js3awYizuuUKfD663nLWtQIY2b/etx2Tis+Xr6DC577nh/W7/PPsSdPLv9GZnAak8eMOX5l4PE489bIXKX4MhznHcA8Ecl9irk5cKPfIjImEEaMyP9E7LUsDLgT6Ne2AXfN+JWrX/uZq05vygMXtCUmIqzsjx0IEyc607x5kJUV6GhMABSZFFT1KxGJA9q6i9aqqg1hZaqsrk1r88XtfXhm9u+8+v0m5q3dw7/+0pH+7WIDHZoxpebLcJyRwD3ALar6K9BURC72e2TGBLGIMA8PXtiO98ecQXREKH97cwk3v72MQ8fsSWhTsfnSpjAFyAB6ufNJwKN+i8iYCiShaW0+u7UPd53Xmtmrd/PA90d5c8EWsnNsnAZTMfmSFFqq6r+BTABVPQqUTefzxlQC4aEh3No/jq/vOIvTaobw0CerGfziDyzZciDQoRlTbL4khQwRqY77IJuItASsTcGYk7SoF8Xd3SJ48aquHDiSwWUv/cTt039h5+GjgQ7NGJ/5cvfRQ8BXQBMRmQacCVznz6CMqahEhIs7N+Kctg2YNG8jL3+3ia9X72J0n9O48eyWRFXz5U/OmMAp9EpBRARYCwzFSQTvAN1UdZ7fIzOmAosMD+WuAW349s6zObddLP+Zs4G+T81j6sKtZFq33CaIFZoU3CHcPlLV/ar6uap+pqpl8MSOMVVDkzqRvHhVAjPHnEHzupH846NVDHj2Oz5e/gc51hhtgpAvbQoLRaS73yMxphJLbFabGTf24rVru1EtNITbpy/ngue/56tVOy05mKDiSwVnP+AmEdkCHMG580hVtbM/AzOmshER+reLpV+bBny2cifPffM7N01dRpvYGG45pxUXdjoFT4jd2GcCy5ekcIHfozCmCgkJEQZ3acRFnU7hsxU7eGHOBm595xeenf07o886jUsSTqVaqPVMagKjsPEUIoCbgFbASuA1d9xlY0wZ8IQIQ+JPZVDnRny1eheT5m3k/g9W8szs37n2jOaMOL0ptSLDAx2mqWIKu1J4E+eBte9xrhbaA7eXR1DGVCUhIcKFnU7hgo4N+XHDfl7+biNPfr2OF+dsYGjCqVx7RnNax8YEOkxTRRSWFNqraicAEXkNWFQ+IRlTNYkIvePq0TuuHmt3JfP6D5t5b2kS037exhkt63J1z2ac1z6WMI8v94cYUzKFJYW8UURUNct5ZMEYUx7aNqzBvy/rwv0XtGP64m1MW7iNsdOWUT+mGsMSG3N5tyY0rxcV6DBNJVRYUugiIsnuewGqu/O5dx/V8Ht0xlRxdaLCGdu3FTee1ZL5v+9h6sJtvDR/IxPnbeT0FnW4LLExF3Q6hWh7UtqUkQJ/k1TVbn8wJkh4QoRz2sZyTttYdh1OZ+ayJGYs2c4976/gnx+v4vwODRkS34g+cfWtesmUiv17YUwF07BmBDf3a8XYvi1Ztu0QM5cl8cXKnXy8fAe1IsMY2KEhF3Y6hV4t61qCMMVmScGYCkpESGxWm8RmtXl4UAe+X7+Xj5fv4NNfdzB98XZqRYbRv20s53eIpU9cfaqH28W/KZolBWMqgfDQEPq3i6V/u1jSM7OZ//tevl61i9m/7WLmsiSqhYbQu1U9+reLpW+b+jSqVT3QIZsgZUnBmEomIszD+R0acn6HhmRm57Bw036+XbOHb9bs5tu1ewBoHRvN2a3r0zuuPj2a17GrCJPHkoIxlViYJ4Q+cfXpE1efhwa1Z8OeVOat28vcdXt4c8FWXvl+M+GeEBKa1eKMlvXoeVpdujSpGeiwTQBZUjCmihAR4mJjiIuN4YazTuNoRjY/b97Pjxv28dOm/Tz7ze+oOlVRLWJgybF1dGtem65Na1OzeligwzflxJKCMVVU9XAPfds0oG+bBgAcSstg0eYD/Lz5AN+u2Mqk+RvJnut0692qQTQJTWvRuXEt4pvUok3DGLuzqZKypGCMAaBWZDgDOjRkQIeG9IneQ/devfl1+yGWbTvIsm2HmP3bbmYsSQIg3BNC21Ni6NCoJh0a1aB9oxq0bRhDZLidUio6+waNMfmKqhbKGa3qcUaregCoKkkHj7J8+yFW/XGYlX8c5vMVO3hn0TYARKB53SjaxMbQumEMrWOjiWsQQ/N6kdYVeAViScEY4xMRoUmdSJrUiWRQl0aAkyj+OHSU1TuSWbszhbW7klm7K4VZv+0id0A5T4jQtE4kLetHcVr9aFrUi6J53Sia14skNiaCEBtYKKhYUjDGlJiI0Lh2JI1rR3J+h4Z5y9Mzs9m4N5UNe45Pm/Ye4bv1+8jIyslbr1poCE3rRNLUTTaNa1ence3qnFrLeV8rMgzrjLN8WVIwxpS5iDCP295w4u2t2TnKjkNH2bo/jc37j7Bt/xG27k9j24E0ft58gNRjWSftJ4RGtarTqGZ1GtaMoGGNCGJrRhAbU40GNSJoEFON+jHVrNG7DFlSMMaUG0/I8Sqo3nH1TvhMVTmUlskfh46SdPAoSQfT2Hk4nR2HjrIrOZ0fN+xjT8oxsnPrpbzUiQqnXnQ49aKrUS+6Wt58nahq7NydReTmA9SODKNWZDi1IsMsiRTCkoIxJiiICLWjwqkdFU7HU/N/gC47R9mfeoxdyensST7GnpRj7ElJZ1/qMfamHGNfagYrkg6xLzXjhKuOF3756YT9xFQLpWZkGDWr/3mqUT2MmIhQZ6oWRnREKNHVnPmoas77aqEhlbZay69JQUQGAs8DHuBVVX38pM/F/fxCIA24TlWX+TMmYyqVadNg3DjYtg2aNoXx42HKFPj22+Pr9O8Po0b9eT3487Iff4TJk+GJJ+Dcc2H0aJg40bfjjhhR8HJfts89dnY2eDz5HtsTIjT4dCYN8ivLP09cln75FRxMy2DW/AW0bNeFA2kZHJr/Iwe/nsPBYzkk143lUOcEDoc1YP2eVA4fzSQlPZP0zJw/x3sST4gQFe4hulookdVCiQr3UD3cQ2R4qPMa5iEy3ENEuIfqYe4U7iEi1EO1sBAiwjxUC3Vevd9XCw0h3J2qhYYQ7in/5OO3pCAiHmACcB6QBCwWkU9U9Tev1S4A4tzpdGCS+2qMKcq0ac6JMy3Nmd+6Fa6++s/rffvtiUli61b4619BFTIzjy+75hrI8TohZmfDpEnOe++Tc37HHT3aOam/+eafl8OJiSG/7Utz7FGjnPthMzJOOG4EcMqIETSr4XGqqqZNgwe8tgX4LNJJRGOOx3csK5uU9Cx3yiQ1PYvk9CyOHMviSIaz/MgxZ0rLyCYtI5vUY1kczchmT0o6aRnZHHWXp2dmcyyr6CRTmHBPCGEeITw0hL6NoG/fUu2uSP68UugBbFDVTQAiMh0YAngnhSHAW6qqwEIRqSUip6jqTj/GZUzlMG7ciSe44sg9gXrLKeDkNXnyiSfm/I6blnb8v/yTl48bd2JSyG/70hw7M5M/8fW4+axXLdRDtWgP9aKr5R9TMWXnKOmZToI4mplNemaOmyyyOZaZQ7r7mpHtLM/IyiHdnT+WlUNGVg6Z2c5rrWO7yiSmwohzPvbDjkUuAwaq6vXu/EjgdFW9xWudz4DHVfUHd/5b4D5VXXLSvkYDowFiY2MTp0+fXqKYUlNTiY6OLtG2wcbKEpzKtSxLl/pt16mNGxOdlHR8QWJi6Y5bmu1LuW3ed1LYtt7HCGKl+f3q16/fUlXtVuSKquqXCRiG046QOz8SeOGkdT4HenvNfwskFrbfxMRELam5c+eWeNtgY2UJTuValmbNVJ1KoDKf5j711PF5j8e343o8+S9v1qzkcft67EKOm/edFLTtyfEFsdL8fgFL1Idztz/vy0oCmnjNNwZ2lGAdY0x+xo+HyMiSbRseDmEn9XwaUsDpILddoLDjRkY66+W3PLchuLDtS3PssDCnPCU5bn7rVXH+TAqLgTgRaSEi4cAVwCcnrfMJcI04egKH1doTjPHNiBFOnXuzZk5Da7NmMHWqc7eRt/79neXe673+unOXkveyt96CMWOcO3/AeR0z5s93H+V33Ny6//yWn3z3UX7bl+bYU6Y45SnJcfNbr6rz5XKipBPOraa/AxuBce6ym4Cb3PeCc4fSRmAl0K2ofVr1kcPKEpwqS1kqSzlUrSy58LH6yK/PKajqF8AXJy17yeu9Ajf7MwZjjDG+s2e9jTHG5LGkYIwxJo8lBWOMMXksKRhjjMljScEYY0weSwrGGGPyWFIwxhiTx28d4vmLiOwFtpZw83rAvjIMJ5CsLMGpspSlspQDrCy5mqlq/aJWqnBJoTREZIn60ktgBWBlCU6VpSyVpRxgZSkuqz4yxhiTx5KCMcaYPFUtKUwOdABlyMoSnCpLWSpLOcDKUixVqk3BGGNM4aralYIxxphCWFIwxhiTp9ImBRGJEJFFIvKriKwWkUfc5XVEZLaIrHdfawc6Vl+IiEdEfhGRz9z5ilqOLSKyUkSWi8gSd1lFLUstEXlfRNaKyBoR6VURyyIibdzvI3dKFpG/V9Cy3OH+va8SkXfc80CFKweAiNzulmO1iPzdXeb3slTapAAcA85R1S5APDDQHfLzfuBbVY0DvnXnK4LbgTVe8xW1HAD9VDXe637rilqW54GvVLUt0AXn+6lwZVHVde73EQ8kAmnAh1SwsojIqcBtOCM4dgQ8OMMAV6hyAIhIR+AGoAfO79bFIhJHeZTFl+HZKvoERALLgNOBdcAp7vJTgHWBjs+H+Bu7vwDnAJ+5yypcOdxYtwD1TlpW4coC1AA2496sUZHLclL8A4AfK2JZgFOB7UAdIBT4zC1PhSqHG+cw4FWv+X8C95ZHWSrzlUJulctyYA8wW1V/BmJVdSeA+9ogkDH66DmcX4gcr2UVsRwACswSkaUiMtpdVhHLchqwF5jiVuu9KiJRVMyyeLsCeMd9X6HKoqp/AE8B24CdwGFVnUUFK4drFXCWiNQVkUic8e6bUA5lqdRJQVWz1bkkbgz0cC/JKhQRuRjYo6pLAx1LGTlTVROAC4CbReSsQAdUQqFAAjBJVbsCR6gA1RKFEZFwYDDwXqBjKQm3fn0I0AJoBESJyNWBjapkVHUN8AQwG/gK+BXIKo9jV+qkkEtVDwHzgIHAbhE5BcB93RPA0HxxJjBYRLYA04FzRGQqFa8cAKjqDvd1D069dQ8qZlmSgCT36hPgfZwkURHLkusCYJmq7nbnK1pZzgU2q+peVc0EPgDOoOKVAwBVfU1VE1T1LOAAsJ5yKEulTQoiUl9Earnvq+P8wqwFPgGudVe7Fvg4MBH6RlUfUNXGqtoc59J+jqpeTQUrB4CIRIlITO57nPreVVTAsqjqLmC7iLRxF/UHfqMClsXLlRyvOoKKV5ZtQE8RiRQRwflO1lDxygGAiDRwX5sCQ3G+G7+XpdI+0SwinYE3ce5ACAFmqOr/ikhdYAbQFOeXaJiqHghcpL4Tkb7A3ap6cUUsh4ichnN1AE71y9uqOr4ilgVAROKBV4FwYBMwCvd3jYpXlkicRtrTVPWwu6zCfS/urefDcapafgGuB6KpYOUAEJHvgbpAJnCnqn5bHt9JpU0Kxhhjiq/SVh8ZY4wpPksKxhhj8lhSMMYYk8eSgjHGmDyWFIwxxuSxpGCCmvuYf27vnbtE5A+v+fBAx5cfEfmriDT04/6jRGSeiIS483eLyDIRucydDxWR7JN6Pm0iIvEi8qq/4jKVQ2igAzCmMKq6H6eXW0TkYSBVVZ8KaFBOLB5VzS7g47/idMC4qxj7C1VVX7sxuB54T1VzRKQmTi+aPXCeqn7fXSfF7eLF23YRaSkip7r9BBnzJ3alYCosEblWnDEzlovIRBEJcf9LPiQiT7r/PX8tIqeLyHwR2SQiF7rbXi8iH7qfrxORf/i430dFZBFOX1qPiMhicfq8f0kcw3GS2Lu5VzMikuT1dH1PEfnGff+oiLwsIrNxOtYLFZFn3GOvEJHrCyj6CI4/yVrcv+HPcB7uMiZflhRMhSRO54aXAGe4/xGH4nQDAlATmOV2vJcBPIzT5cEw4H+9dtPD3SYBuMqtXilqv8tUtYeq/gQ8r6rdgU7uZwNV9V1gOTBcnTEKMoooSldgkKqOBEbjdH7YA+iO02Fg05PKHQE0VtUkAFU9iNOd8mJO7KIixqvq6H2v5UuAPkXEZKowqz4yFdW5OCfOJU43N1TH6aYB4Kiqznbfr8TpQjlLRFYCzb328bV7UkVEPgJ64/xNFLTfDI530wHQX0TuASKAesBS4MtiluNjVU133w8A2omIdxKKw+nOIFcDnM7R8qjqo8CjJ+03v+ojcDpQa1TMGE0VYknBVFQCvK76/+3dv0oDQRDH8e+ARSLYWacRrKwEH8DGFxDs7SSVlZ1g4SPYiIKI4DsICVrEQm0CvoAWivinEkmhOBa7t55H/nhYhMPfp0pyd5M7CDfs3mbGN398aDZBuHlnPgld+LLX+d98scaLj4jb81gXJtYK2gHm3f3OzLYJyaGfD75H5cV93grX1HT39oA4AL0h3/MbtRhDpC9NH0lVtYAVM5uGtEqpMeKYoiULfZYnCXX4z0vErROSzLOFyq/LuW2vwFTu/bvbWOwAAADXSURBVA2hzSWF/YpOgGZMQFnv5Hp+B3d/Amp/WHk1S6hMK9KXRgpSSe5+HStituLSzHdgDbgvEaYDHAMzwJG7dyFV2hwa191fzOyQcIO9BS5ymw+AfTPrEZ5bbAF7ZvYAXA45n11C9ctunLp6JCSrojahT8BZiWvNLFKR0tEyHqqSKv9SXNkz5+7r4z6XssxsgTDNtFryuDpwSuh+N2g5rfxzmj4SqRh3vwI62Z/XSmgAG0oIMoxGCiIikmikICIiiZKCiIgkSgoiIpIoKYiISKKkICIiyRcB+QfBIrGzDQAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "plt.plot(data_pred['Temperature'], data_pred['Frequency'], label='Predicted')\n",
    "plt.scatter(data[\"Temperature\"], data[\"Frequency\"], color='red', label='Observed')\n",
    "plt.xlabel(\"Temperature (°F)\")\n",
    "plt.ylabel(\"Predicted Frequency of Failure\")\n",
    "plt.legend()\n",
    "plt.grid(True)\n",
    "plt.title(\"Logistic regression model\")\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Uncertainty of estimates (Confidence interval)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 30,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Confidence interval for parameters :\n",
      "                     0          1\n",
      "Intercept   -9.569730  19.739685\n",
      "Temperature -0.341358   0.110156\n"
     ]
    }
   ],
   "source": [
    "# Taking uncertainty into account\n",
    "# Added confidence interval representation\n",
    "conf = logmodel.conf_int()\n",
    "print(\"Confidence interval for parameters :\\n\", conf)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "*Interpretation :*\n",
    "\n",
    "* Intercept : [-9.5697, 19.7397]\n",
    "The confidence interval for the intercept is very wide, which means there is a large uncertainty about this estimate. This indicates that the data do not provide a precise estimate of the intercept.\n",
    "\n",
    "* Temperature : [-0.3414, 0.1102]\n",
    "The confidence interval for the temperature coefficient is also wide and includes zero. This means that the effect of temperature could be zero or even slightly positive or negative. In other words, one cannot conclude with certainty that temperature has a significant impact on the probability of failure.\n",
    "\n",
    "*Conclusion :* These confidence intervals reveal a large uncertainty in the model estimates. In particular, the temperature interval being wide and containing zero, this reinforces the fact that the effect of temperature on failures is not statistically significant in this analysis. This suggests that a model revision or additional data would be necessary to draw reliable conclusions."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## II. Second attempt of Improvement"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Trying non-linear transformations\n",
    "Let's try adding a data transformation, for example by adding a quadratic variable \"*Temperature_squared* \" to the model to capture nonlinear effects:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 31,
   "metadata": {},
   "outputs": [],
   "source": [
    "data[\"Temperature_squared\"] = data[\"Temperature\"] ** 2\n",
    "logmodel = sm.GLM(data['Frequency'], data[['Intercept', 'Temperature', 'Temperature_squared']], family=sm.families.Binomial(sm.families.links.logit())).fit()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 32,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "                 Generalized Linear Model Regression Results                  \n",
      "==============================================================================\n",
      "Dep. Variable:              Frequency   No. Observations:                   23\n",
      "Model:                            GLM   Df Residuals:                       20\n",
      "Model Family:                Binomial   Df Model:                            2\n",
      "Link Function:                  logit   Scale:                          1.0000\n",
      "Method:                          IRLS   Log-Likelihood:                -3.8798\n",
      "Date:                Mon, 21 Oct 2024   Deviance:                       2.9319\n",
      "Time:                        09:32:48   Pearson chi2:                     4.48\n",
      "No. Iterations:                     6   Covariance Type:             nonrobust\n",
      "=======================================================================================\n",
      "                          coef    std err          z      P>|z|      [0.025      0.975]\n",
      "---------------------------------------------------------------------------------------\n",
      "Intercept              22.1261     58.284      0.380      0.704     -92.109     136.362\n",
      "Temperature            -0.6509      1.814     -0.359      0.720      -4.207       2.905\n",
      "Temperature_squared     0.0041      0.014      0.297      0.766      -0.023       0.031\n",
      "=======================================================================================\n"
     ]
    }
   ],
   "source": [
    "print(logmodel.summary())"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "*Interpretation :*\n",
    "\n",
    "* Intercept: 22.1261\n",
    "The intercept is very high but its wide confidence interval and p-value (0.704) indicate that this estimate is extremely uncertain and not significant.\n",
    "\n",
    "* Temperature coefficient: -0.6509\n",
    "This coefficient remains negative, which would indicate that, overall, an increase in temperature could reduce the probability of failure. However, with a p-value of 0.720, this effect is not statistically significant. The very wide confidence interval ([-4.207, 2.905]) confirms this uncertainty.\n",
    "\n",
    "* Temperature_squared coefficient: 0.0041\n",
    "Adding this quadratic component allows us to test whether the relationship between temperature and the probability of failure is curved. However, this coefficient is very small and, with a p-value of 0.766, it is not significant either. The confidence interval includes zero ([-0.023, 0.031]), which means that the quadratic effect might not exist.\n",
    "\n",
    "*Conclusion :*\n",
    "\n",
    "* Adding a quadratic component did not improve the significance of the model. The coefficients of temperature and Temperature_squared are not statistically significant (p-values ​​are greater than 0.05), meaning that in this analysis there is no clear evidence that temperature has a significant effect on the probability of O-ring failure, either linearly or nonlinearly.\n",
    "\n",
    "* The wide confidence intervals and high p-values ​​indicate that there is a large uncertainty in the estimates. This could be due to the small sample size and the lack of data for very low temperatures, where the risk of failure could be higher."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Exploring other explanatory variables\n",
    "We'll try to include additional variables that might influence the probability of o-ring failure, such as pressure. We will also leave out temperature and its quadratic transformation, and add pressure as a factor in the model."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 33,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "                 Generalized Linear Model Regression Results                  \n",
      "==============================================================================\n",
      "Dep. Variable:              Frequency   No. Observations:                   23\n",
      "Model:                            GLM   Df Residuals:                       19\n",
      "Model Family:                Binomial   Df Model:                            3\n",
      "Link Function:                  logit   Scale:                          1.0000\n",
      "Method:                          IRLS   Log-Likelihood:                -3.7883\n",
      "Date:                Mon, 21 Oct 2024   Deviance:                       2.7490\n",
      "Time:                        09:32:50   Pearson chi2:                     4.14\n",
      "No. Iterations:                     6   Covariance Type:             nonrobust\n",
      "=======================================================================================\n",
      "                          coef    std err          z      P>|z|      [0.025      0.975]\n",
      "---------------------------------------------------------------------------------------\n",
      "Intercept               8.4977     64.779      0.131      0.896    -118.467     135.462\n",
      "Temperature            -0.2808      1.962     -0.143      0.886      -4.126       3.564\n",
      "Temperature_squared     0.0014      0.015      0.093      0.926      -0.028       0.031\n",
      "Pressure                0.0079      0.020      0.396      0.692      -0.031       0.047\n",
      "=======================================================================================\n"
     ]
    }
   ],
   "source": [
    "# Adding a column for the quadratic term of temperature\n",
    "data[\"Temperature_squared\"] = data[\"Temperature\"] ** 2\n",
    "\n",
    "# Logistic regression model with temperature, quadratic temperature term, and pressure\n",
    "logmodel_extended = sm.GLM(data['Frequency'], \n",
    "                           data[['Intercept', 'Temperature', 'Temperature_squared', 'Pressure']], \n",
    "                           family=sm.families.Binomial(sm.families.links.logit())).fit()\n",
    "\n",
    "# Summary of results\n",
    "print(logmodel_extended.summary())"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "*Interpretation :*\n",
    "\n",
    "* Intercept: 8.4977\n",
    "The intercept coefficient remains very high, but with a p-value of 0.896, it is not statistically significant. This indicates a large uncertainty in this estimate.\n",
    "\n",
    "* Temperature coefficient: -0.2808\n",
    "The temperature coefficient remains negative, but it is still insignificant with a p-value of 0.886. This means that there is still no clear evidence that temperature (as a linear variable) significantly affects the probability of o-ring failure.\n",
    "\n",
    "* Temperature_squared coefficient: 0.0014\n",
    "The addition of the quadratic component still does not seem significant, with a p-value of 0.926. This means that adding a nonlinear term for temperature did not improve the modeling of the temperature effect.\n",
    "\n",
    "* Pressure coefficient: 0.0079\n",
    "The coefficient for pressure is positive, which would indicate that in this model, a small increase in pressure could increase the probability of failure. However, the p-value is 0.692, which means that this effect is not statistically significant. The confidence interval ([-0.031, 0.047]) includes zero, which confirms that the effect of pressure is uncertain.\n",
    "\n",
    "*Conclusion :*\n",
    "\n",
    "* Neither temperature, nor its quadratic component, nor pressure are significant in this model. The high p-values ​​indicate that none of these variables has a clear effect on the probability of O-ring failure based on the current data.\n",
    "\n",
    "* Adding pressure did not reduce uncertainty in the model. The wide confidence intervals and high p-values ​​show that it is difficult to isolate a factor that significantly influences O-ring failures."
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.6.4"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}