{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "# Analysis of the risk of failure of the O-rings on the Challenger shuttle" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "On January 27, 1986, the day before the takeoff of the shuttle *Challenger*, had\n", "a three-hour teleconference was held between \n", "Morton Thiokol (the manufacturer of one of the engines) and NASA. The\n", "discussion focused on the consequences of the\n", "temperature at take-off of 31°F (just below\n", "0°C) for the success of the flight and in particular on the performance of the\n", "O-rings used in the engines. Indeed, no test\n", "had been performed at this temperature.\n", "\n", "The following study takes up some of the analyses carried out that\n", "night with the objective of assessing the potential influence of\n", "the temperature and pressure to which the O-rings are subjected\n", "on their probability of malfunction. Our starting point is \n", "the results of the experiments carried out by NASA engineers\n", "during the six years preceding the launch of the shuttle\n", "Challenger." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Loading the data\n", "We start by loading this data:" ] }, { "cell_type": "code", "execution_count": 1, "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", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \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": 1, "metadata": {}, "output_type": "execute_result" } ], "source": [ "import numpy as np\n", "import pandas as pd\n", "import warnings\n", "warnings.filterwarnings(\"ignore\", category=DeprecationWarning) \n", "warnings.filterwarnings(\"ignore\", category=FutureWarning) \n", "data = pd.read_csv(\"../data/shuttle.csv\")\n", "data" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The data set shows us the date of each test, the number of O-rings (there are 6 on the main launcher), the temperature (in Fahrenheit) and pressure (in psi), and finally the number of identified malfunctions." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Graphical inspection\n", "Flights without incidents do not provide any information\n", "on the influence of temperature or pressure on malfunction.\n", "We thus focus on the experiments in which at least one O-ring\n", "was defective." ] }, { "cell_type": "code", "execution_count": 2, "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", "
DateCountTemperaturePressureMalfunction
111/12/81670501
82/03/846572001
94/06/846632001
108/30/846702001
131/24/856532002
2010/30/856752002
221/12/866582001
\n", "
" ], "text/plain": [ " Date Count Temperature Pressure Malfunction\n", "1 11/12/81 6 70 50 1\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", "13 1/24/85 6 53 200 2\n", "20 10/30/85 6 75 200 2\n", "22 1/12/86 6 58 200 1" ] }, "execution_count": 2, "metadata": {}, "output_type": "execute_result" } ], "source": [ "data = data[data.Malfunction>0]\n", "data" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "We have a high temperature variability but\n", "the pressure is almost always 200, which should\n", "simplify the analysis.\n", "\n", "How does the frequency of failure vary with temperature?" ] }, { "cell_type": "code", "execution_count": 3, "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYIAAAEKCAYAAAAfGVI8AAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjQuMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/MnkTPAAAACXBIWXMAAAsTAAALEwEAmpwYAAAVc0lEQVR4nO3de5CddX3H8fd3kxCCSTGFdnWSIFgilUqKsICUWtd6KdoxqYVa6FiqtcS20HrpVKh1KEPrTMEqrQ5eIqUKVSOXqmublkrtqrUCCRrDrdEtItloQWO4LIZc2G//OM/qye7Z5Gyyzzns/t6vmZ09z/V898vD+eS5nOeJzESSVK6ebhcgSeoug0CSCmcQSFLhDAJJKpxBIEmFMwgkqXC1BUFEXBMRD0XEXZNMj4h4b0QMRcSmiDiprlokSZOrc4/gI8CZ+5j+CmB59bMa+ECNtUiSJlFbEGTmF4Ef7GOWVcC12XAr8PSIeGZd9UiSWpvbxfdeAmxpGh6uxn13/IwRsZrGXgMLFiw4edmyZR0p8GCNjo7S0+NpmGb2ZCJ7MpE9ae1g+vKNb3zj+5n5U62mdTMI2paZa4A1AH19fblhw4YuV9SewcFB+vv7u13GU4o9mcieTGRPWjuYvkTEtyeb1s3I3Qo0/9N+aTVOktRB3QyCAeC86uqhFwCPZOaEw0KSpHrVdmgoIj4B9ANHRsQw8BfAPIDM/CCwDnglMAT8EHh9XbVIkiZXWxBk5rn7mZ7ABXW9vySpPZ6Wl6TCGQSSVDiDQJIKZxBIUuEMAkkqnEEgSYUzCCSpcAaBJBXOIJCkwhkEklQ4g0CSCmcQSFLhDAJJKpxBIEmFMwgkqXAGgSQVziCQpMIZBJJUOINAkgpnEEhS4QwCSSqcQSBJhTMIJKlwBoEkFc4gkKTCGQSSVDiDQJIKZxBIUuEMAkkqnEEgSYUzCCSpcAaBJBXOIJCkwhkEklQ4g0CSCldrEETEmRGxOSKGIuLiFtOPioj/jIivRcSmiHhlnfVIkiaqLQgiYg5wFfAK4Hjg3Ig4ftxs7wCuz8znA+cA76+rHklSa3XuEZwKDGXmfZm5C1gLrBo3TwI/Ub0+HPhOjfVIklqIzKxnxRFnA2dm5u9Vw78NnJaZFzbN80zg34HFwNOAl2bmHS3WtRpYDdDb23vy2rVra6l5uo2MjLBw4cJul/GUYk8msicT2ZPWDqYvL37xi+/IzL5W0+YeVFUH71zgI5n57og4HbguIp6XmaPNM2XmGmANQF9fX/b393e+0gMwODjITKm1U+zJRPZkInvSWl19qfPQ0FZgWdPw0mpcszcA1wNk5leAQ4Eja6xJkjROnUGwHlgeEcdExCE0TgYPjJvnAeAlABHxXBpB8L0aa5IkjVNbEGTmHuBC4GbgXhpXB90dEZdFxMpqtj8Bzo+IrwOfAF6XdZ20kCS1VOs5gsxcB6wbN+6Sptf3AGfUWYMkad/8ZrEkFc4gkKTCGQSSVDiDQJIKZxBIUuEMAkkqnEEgSYUzCCSpcAaBJBXOIJCkwhkEklQ4g0CSCmcQSFLhDAJJKpxBIEmFMwgkqXAGgSQVziCQpMIZBJJUOINAkgpnEEhS4QwCSSqcQSBJhTMIJKlwBoEkFc4gkKTCGQSSVDiDQJIKZxBIUuEMAkkqnEEgSYUzCCSpcAaBJBXOIJCkwtUaBBFxZkRsjoihiLh4knleExH3RMTdEfHxOuuRJE00t64VR8Qc4CrgZcAwsD4iBjLznqZ5lgN/BpyRmdsj4qfrqkeS1FpbewQRccIBrPtUYCgz78vMXcBaYNW4ec4HrsrM7QCZ+dABvI8k6SC0u0fw/oiYD3wE+FhmPtLGMkuALU3Dw8Bp4+Z5DkBEfBmYA1yamf82fkURsRpYDdDb28vg4GCbZXfXyMjIjKm1U+zJRPZkInvSWl19aSsIMvOF1WGc3wXuiIjbgX/IzM9Nw/svB/qBpcAXI+KEzHx43PuvAdYA9PX1ZX9//0G+bWcMDg4yU2rtFHsykT2ZyJ60Vldf2j5ZnJnfBN4BXAS8CHhvRPxPRPz6JItsBZY1DS+txjUbBgYyc3dmfgv4Bo1gkCR1SLvnCFZExJXAvcAvA6/KzOdWr6+cZLH1wPKIOCYiDgHOAQbGzfNpGnsDRMSRNA4V3TfFv0GSdBDaPUfwPuBq4O2ZuWNsZGZ+JyLe0WqBzNwTERcCN9M4/n9NZt4dEZcBGzJzoJr28oi4B3gS+NPM3HYQf48kaYraDYJfBXZk5pMAEdEDHJqZP8zM6yZbKDPXAevGjbuk6XUCb61+JEld0O45gluABU3Dh1XjJEkzXLtBcGhmjowNVK8Pq6ckSVIntRsEj0fESWMDEXEysGMf80uSZoh2zxG8GbghIr4DBPAM4DfrKkqS1DntfqFsfUT8LHBcNWpzZu6uryxJUqdM5aZzpwBHV8ucFBFk5rW1VCVJ6pi2giAirgN+BthI43p/gAQMAkma4drdI+gDjq+u+5ckzSLtXjV0F40TxJKkWabdPYIjgXuqu47uHBuZmStrqUqS1DHtBsGldRYhSeqedi8f/UJEPAtYnpm3RMRhNG4kJ0ma4dq9DfX5wI3Ah6pRS2jcQlqSNMO1e7L4AuAM4FH40UNqfNC8JM0C7QbBzuoB9ABExFwa3yOQJM1w7QbBFyLi7cCCiHgZcAPw2frKkiR1SrtBcDHwPeBO4I00HjbT8slkkqSZpd2rhkaBD1c/kqRZpN17DX2LFucEMvPZ016RJKmjpnKvoTGHAr8B/OT0lyNJ6rS2zhFk5ramn62Z+bc0HmgvSZrh2j00dFLTYA+NPYSpPMtAkvQU1e6H+bubXu8B7gdeM+3VSJI6rt2rhl5cdyGSpO5o99DQW/c1PTPfMz3lSJI6bSpXDZ0CDFTDrwJuB75ZR1GSpM5pNwiWAidl5mMAEXEp8C+Z+dq6CpMkdUa7t5joBXY1De+qxkmSZrh29wiuBW6PiE9Vw78GfLSWiiRJHdXuVUPvjIh/BV5YjXp9Zn6tvrIkSZ3S7qEhgMOARzPz74DhiDimppokSR3U7qMq/wK4CPizatQ84B/rKkqS1Dnt7hG8GlgJPA6Qmd8BFtVVlCSpc9oNgl2ZmVS3oo6Ip9VXkiSpk9oNgusj4kPA0yPifOAWfEiNJM0K+w2CiAjgk8CNwE3AccAlmfm+NpY9MyI2R8RQRFy8j/nOioiMiL7J5pEk1WO/l49mZkbEusw8AfhcuyuOiDnAVcDLgGFgfUQMZOY94+ZbBLwJuG1KlUuSpkW7h4a+GhGnTHHdpwJDmXlfZu4C1gKrWsz3l8DlwBNTXL8kaRq0+83i04DXRsT9NK4cCho7Cyv2scwSYEvT8HC1nh+pHnizLDP/JSL+dLIVRcRqYDVAb28vg4ODbZbdXSMjIzOm1k6xJxPZk4nsSWt19WWfQRARR2XmA8CvTPcbR0QP8B7gdfubNzPXAGsA+vr6sr+/f7rLqcXg4CAzpdZOsScT2ZOJ7ElrdfVlf3sEn6Zx19FvR8RNmXnWFNa9FVjWNLy0GjdmEfA8YLBxPppnAAMRsTIzN0zhfSRJB2F/5wii6fWzp7ju9cDyiDgmIg4BzuHHzzMgMx/JzCMz8+jMPBq4FTAEJKnD9hcEOcnr/crMPcCFwM3AvcD1mXl3RFwWESunVqYkqS77OzT08xHxKI09gwXVa/jxyeKf2NfCmbkOWDdu3CWTzNvfVsWSpGm1zyDIzDmdKkSS1B1TuQ21JGkWMggkqXAGgSQVziCQpMIVEwTbRnby9S0Ps21kZ7dLkaQp2zaykx27n6zlM6yIIPjMxq2ccfnnee3Vt3HG5Z9nYOPW/S8kSU8RY59h3/re47V8hs36INg2spOLbtrEE7tHeWznHp7YPcrbbtrknoGkGaH5M+zJzFo+w2Z9EAxv38G8nr3/zHk9PQxv39GliiSpfZ34DJv1QbB08QJ2j47uNW736ChLFy/oUkWS1L5OfIbN+iA4YuF8rjhrBYfO62HR/LkcOq+HK85awREL53e7NEnar+bPsDkRtXyGtftgmhlt5YlLOOPYIxnevoOlixcYApJmlLHPsNu/8l98eeUvTvtnWBFBAI1UNQAkzVRHLJzPgnlzavkcm/WHhiRJ+2YQSFLhDAJJKpxBIEmFMwgkqXAGgSQVziCQpMIZBJJUOINAkgpnEEhS4QwCSSqcQSBJhTMIJKlwBoEkFc4gkKTCGQSSVDiDQJIKZxBIUuEMAkkqnEEgSYUzCCSpcLUGQUScGRGbI2IoIi5uMf2tEXFPRGyKiP+IiGfVWY8kaaLagiAi5gBXAa8AjgfOjYjjx832NaAvM1cANwJX1FWPJKm1OvcITgWGMvO+zNwFrAVWNc+Qmf+ZmT+sBm8FltZYjySphbk1rnsJsKVpeBg4bR/zvwH411YTImI1sBqgt7eXwcHBaSqxXiMjIzOm1k6xJxPZk4nsSWt19aXOIGhbRLwW6ANe1Gp6Zq4B1gD09fVlf39/54o7CIODg8yUWjvFnkxkTyayJ63V1Zc6g2ArsKxpeGk1bi8R8VLgz4EXZebOGuuRJLVQ5zmC9cDyiDgmIg4BzgEGmmeIiOcDHwJWZuZDNdYiSZpEbUGQmXuAC4GbgXuB6zPz7oi4LCJWVrO9C1gI3BARGyNiYJLVSZJqUus5gsxcB6wbN+6SptcvrfP9S7JtZCfD23ewdPECjlg4v7ZlSmBfumvowcfY/sPdDD34GMf2Lup2OUV4Spws1sH5zMatXHTTJub19LB7dJQrzlrByhOXTPsyJbAv3XXJp+/k2lsf4E9O2MNbrvwi551+FJetOqHbZc163mJihts2spOLbtrEE7tHeWznHp7YPcrbbtrEtpHJz7sfyDIlsC/dNfTgY1x76wN7jbv2Kw8w9OBjXaqoHAbBDDe8fQfzevb+zzivp4fh7TumdZkS2Jfu2rjl4SmN1/QxCGa4pYsXsHt0dK9xu0dHWbp4wbQuUwL70l0nLnv6lMZr+hgEM9wRC+dzxVkrOHReD4vmz+XQeT1ccdaKfZ7kPJBlSmBfuuvY3kWcd/pRe4077/SjPGHcAZ4sngVWnriEM449ckpXuhzIMiWwL9112aoTOO8FR3PnHbdyy1teYAh0iEEwSxyxcP6UP7QOZJkS2JfuOrZ3EcOHzTMEOshDQ5JUOINAkgpnEEhS4QwCSSqcQSBJhTMIJKlwBoEkFc4gkKTCGQSSVDiDQJIKZxBIUuEMAkkqnEEgSYUzCCSpcAaBJBXOIJCkwhkEklQ4g0CSCmcQSFLhDAJJKpxBIEmFMwgkqXAGgSQVziCQpMIZBJJUOINAkgpnEEhS4QwCSSpcrUEQEWdGxOaIGIqIi1tMnx8Rn6ym3xYRR9dZjyRpotqCICLmAFcBrwCOB86NiOPHzfYGYHtmHgtcCVxeVz2SpNbq3CM4FRjKzPsycxewFlg1bp5VwEer1zcCL4mIqLEmSdI4c2tc9xJgS9PwMHDaZPNk5p6IeAQ4Avh+80wRsRpYXQ2ORMTmWiqefkcy7m+RPWnBnkxkT1o7mL48a7IJdQbBtMnMNcCabtcxVRGxITP7ul3HU4k9mcieTGRPWqurL3UeGtoKLGsaXlqNazlPRMwFDge21ViTJGmcOoNgPbA8Io6JiEOAc4CBcfMMAL9TvT4b+HxmZo01SZLGqe3QUHXM/0LgZmAOcE1m3h0RlwEbMnMA+HvguogYAn5AIyxmkxl3OKsD7MlE9mQie9JaLX0J/wEuSWXzm8WSVDiDQJIKZxBMk4i4PyLujIiNEbGhGndpRGytxm2MiFd2u85Oi4inR8SNEfE/EXFvRJweET8ZEZ+LiG9Wvxd3u85OmqQnxW4rEXFc09+9MSIejYg3l7yd7KMntWwnniOYJhFxP9CXmd9vGncpMJKZf9OturotIj4KfCkzr66uHjsMeDvwg8z86+oeVIsz86KuFtpBk/TkzRS+rcCPbk2zlcaXTy+g4O1kzLievJ4athP3CFSbiDgc+CUaV4eRmbsy82H2vrXIR4Ff60Z93bCPnqjhJcD/Zua3KXg7Gae5J7UwCKZPAv8eEXdUt8QYc2FEbIqIa0rata0cA3wP+IeI+FpEXB0RTwN6M/O71Tz/B/R2rcLOm6wnUPa2MuYc4BPV65K3k2bNPYEathODYPr8YmaeRONuqxdExC8BHwB+BjgR+C7w7u6V1xVzgZOAD2Tm84HHgb1uR159gbCk45OT9aT0bYXqMNlK4Ibx0wrcToCWPallOzEIpklmbq1+PwR8Cjg1Mx/MzCczcxT4MI07spZkGBjOzNuq4RtpfAg+GBHPBKh+P9Sl+rqhZU/cVoDGP6K+mpkPVsMlbydj9upJXduJQTANIuJpEbFo7DXwcuCusY248mrgrm7U1y2Z+X/Alog4rhr1EuAe9r61yO8An+lCeV0xWU9K31Yq57L3IZBit5Mme/Wkru3Eq4amQUQ8m8ZeADR2/T+eme+MiOto7MIlcD/wxqZjnkWIiBOBq4FDgPtoXPXQA1wPHAV8G3hNZv6gWzV22iQ9eS8FbyvVP6AeAJ6dmY9U446g7O2kVU9q+UwxCCSpcB4akqTCGQSSVDiDQJIKZxBIUuEMAkkq3Ix4eL3UruqSw/+oBp8BPEnjlg7Q+JLfrq4U1kJE9AO7MvO/u1yKCmcQaFbJzG00rrN+Stz9NSLmZuaeSSb3AyNA20Gwn/VJB8RDQ5r1IuLkiPhCdUPAm5tuWzAYEVdGxIbquQCnRMQ/Vfe//6tqnqOr5wZ8rJrnxog4rI31/m00nkvxpoh4VUTcVt1k7paI6I2Io4HfB95S3Vf+hRHxkYg4u6nukep3f0R8KSIGaHwLeU5EvCsi1lc3H3tjRxuqWccg0GwXwPuAszPzZOAa4J1N03dlZh/wQRq3MLgAeB7wuuowE8BxwPsz87nAo8AfRsS8/az3kMzsy8x3A/8FvKC6ydxa4G2ZeX/1nldm5omZ+aX9/B0nAW/KzOcAbwAeycxTgFOA8yPimKm3Rmrw0JBmu/k0Ptg/FxEAc2jctXHMQPX7TuDusa/rR8R9wDLgYWBLZn65mu8fgT8G/m0/6/1k0+ulwCerPYZDgG8dwN9xe2aOLfdyYEXT3sPhwPIDXK9kEGjWCxof8KdPMn1n9Xu06fXY8Nj/H+Pvw5JtrPfxptfvA96TmQPVCeJLJ1lmD9VeekT00AiNVusL4I8y8+ZJ1iNNiYeGNNvtBH4qIk4HiIh5EfFzU1zHUWPLA79F41DP5ims93AajxqEH99NE+AxYFHT8P3AydXrlcC8SdZ3M/AH1eEpIuI5TQ+3kabMINBsNwqcDVweEV8HNgK/MMV1bKbxsKF7gcU0HiqzawrrvRS4ISLuAL7fNP6zwKvHThbTuL/8i6r1nc7eewHNrqZxO++vRsRdwIdw714HwbuPSvtQXd3zz5n5vG7XItXFPQJJKpx7BJJUOPcIJKlwBoEkFc4gkKTCGQSSVDiDQJIK9//DiEewM8RHmgAAAABJRU5ErkJggg==\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": [ "At first glance, the dependence does not look very important, but let's try to\n", "estimate the impact of temperature $t$ on the probability of O-ring malfunction." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Estimation of the temperature influence\n", "\n", "Suppose that each of the six O-rings is damaged with the same\n", "probability and independently of the others and that this probability\n", "depends only on the temperature. If $p(t)$ is this probability, the\n", "number $D$ of malfunctioning O-rings during a flight at\n", "temperature $t$ follows a binomial law with parameters $n=6$ and\n", "$p=p(t)$. To link $p(t)$ to $t$, we will therefore perform a\n", "logistic regression." ] }, { "cell_type": "code", "execution_count": 4, "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", "
Generalized Linear Model Regression Results
Dep. Variable: Frequency No. Observations: 7
Model: GLM Df Residuals: 5
Model Family: Binomial Df Model: 1
Link Function: logit Scale: 1.0000
Method: IRLS Log-Likelihood: -2.5250
Date: Fri, 22 Oct 2021 Deviance: 0.22231
Time: 13:13:23 Pearson chi2: 0.236
No. Iterations: 4 Pseudo R-squ. (CS): 1.926e-05
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 -1.3895 7.828 -0.178 0.859 -16.732 13.953
Temperature 0.0014 0.122 0.012 0.991 -0.238 0.240
" ], "text/plain": [ "\n", "\"\"\"\n", " Generalized Linear Model Regression Results \n", "==============================================================================\n", "Dep. Variable: Frequency No. Observations: 7\n", "Model: GLM Df Residuals: 5\n", "Model Family: Binomial Df Model: 1\n", "Link Function: logit Scale: 1.0000\n", "Method: IRLS Log-Likelihood: -2.5250\n", "Date: Fri, 22 Oct 2021 Deviance: 0.22231\n", "Time: 13:13:23 Pearson chi2: 0.236\n", "No. Iterations: 4 Pseudo R-squ. (CS): 1.926e-05\n", "Covariance Type: nonrobust \n", "===============================================================================\n", " coef std err z P>|z| [0.025 0.975]\n", "-------------------------------------------------------------------------------\n", "Intercept -1.3895 7.828 -0.178 0.859 -16.732 13.953\n", "Temperature 0.0014 0.122 0.012 0.991 -0.238 0.240\n", "===============================================================================\n", "\"\"\"" ] }, "execution_count": 4, "metadata": {}, "output_type": "execute_result" } ], "source": [ "import statsmodels.api as sm\n", "\n", "data[\"Success\"]=data.Count-data.Malfunction\n", "data[\"Intercept\"]=1\n", "\n", "logmodel=sm.GLM(data['Frequency'], data[['Intercept','Temperature']], family=sm.families.Binomial(link=sm.families.links.logit())).fit()\n", "\n", "logmodel.summary()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The most likely estimator of the temperature parameter is 0.0014\n", "and the standard error of this estimator is 0.122, in other words we\n", "cannot distinguish any particular impact and we must take our\n", "estimates with caution." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Estimation of the probability of O-ring malfunction\n", "\n", "The expected temperature on the take-off day is 31°F. Let's try to\n", "estimate the probability of O-ring malfunction at\n", "this temperature from the model we just built:" ] }, { "cell_type": "code", "execution_count": 5, "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXQAAAEKCAYAAAACS67iAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjQuMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/MnkTPAAAACXBIWXMAAAsTAAALEwEAmpwYAAAbNklEQVR4nO3df5RV5X3v8fd3ZkCGH4JBpeAQIauItSo/ZvgVrBmSKJikqC1XJcbU3BBye0NimkqXrOZGa3Wtpthoa62VKjWtVwfisgSzWIHEMjX1Rh0QUIHyo2YiMxpREn6MDsrMfO8fe59hz2FmzpmZczjnPH5ea83i7L2f/eznOZv9mX2ec84z5u6IiEjpKyt0A0REJDcU6CIigVCgi4gEQoEuIhIIBbqISCAU6CIigcgY6Ga22swOmtmrPWw3M/s7M9tvZi+b2fTcN1NERDLJ5g79UWBBL9uvAibFP0uBBwfeLBER6auMge7uzwK/7qXI1cC/eOR5YJSZjc1VA0VEJDsVOajjPOBAYrkpXvdmekEzW0p0F09lZWX1+PHj+3XAjo4OysrCGP5XX4pPKP0A9aVYDaQve/fufcfdz+luWy4CPWvuvgpYBVBTU+NbtmzpVz319fXU1tbmsGWFo74Un1D6AepLsRpIX8zslz1ty8Wvu2YgeatdFa8TEZHTKBeBvh74Yvxpl9nAEXc/ZbhFRETyK+OQi5k9AdQCZ5tZE3A7MAjA3f8R2AB8BtgPvAd8KV+NFRGRnmUMdHdfnGG7A1/LWYtEpGScOHGCpqYmjh8/nvdjjRw5kt27d+f9OKdDNn0ZMmQIVVVVDBo0KOt6T+uboiISlqamJkaMGMGECRMws7we69ixY4wYMSKvxzhdMvXF3Tl06BBNTU1MnDgx63rD+AyQiBTE8ePHGT16dN7D/MPGzBg9enSfX/ko0EVkQBTm+dGf51WBLiISCI2hi0hJKy8v55JLLulcXrduHRMmTChcgwpIgS4iJa2yspLt27d3u83dcfdgpgzI5MPRSxH50GhsbGTy5Ml88Ytf5OKLL+bAgQOsXLmSGTNmcOmll3L77bd3lr377ru54IILuOyyy1i8eDH33HMPALW1taSmJnnnnXc67/jb29tZvnx5Z10PPfQQcPKr/IsWLeLCCy/kxhtvJPpENzQ0NPDxj3+cKVOmMHPmTI4dO8aCBQu6/BK67LLL2LFjx4D7rjt0EcmJv3h6J7veOJrTOi8adya3//7v9lqmtbWVqVOnAjBx4kTuvfde9u3bx/e//31mz57Npk2b2LdvHy+++CLuzsKFC3n22WcZNmwYdXV1bN++nba2NqZPn051dXWvx3rkkUcYOXIkDQ0NvP/++8ydO5crr7wSgG3btrFz507GjRvH3Llzee6555g5cybXX389a9asYcaMGRw9epTKykpuuukmHn30Ue677z727t3L8ePHmTJlyoCfLwW6iJS09CGXxsZGzj//fGbPng3Apk2b2LRpE9OmTQOgpaWFffv2cezYMa699lqGDh0KwMKFCzMea9OmTbz88ss8+eSTABw5coR9+/YxePBgZs6cSVVVFQBTp06lsbGRkSNHMnbsWGbMmAHAmWeeCcC1117L3LlzWblyJatXr+bmm2/OyXOhQBeRnMh0J306DRs2rPOxu7NixQq++tWvdilz33339bh/RUUFHR0dAF0+C+7u3H///cyfP79L+fr6es4444zO5fLyctra2nqsf+jQoVxxxRX88Ic/ZO3atWzdujWrfmWiMXQRCdr8+fNZvXo1LS0tADQ3N3Pw4EEuv/xy1q1bR2trK8eOHePpp5/u3GfChAmdIZu6G0/V9eCDD3LixAkA9u7dy7vvvtvjsSdPnsybb75JQ0MDEH1DNBX0S5Ys4Rvf+AYzZszgrLPOyklfdYcuIkG78sor2b17N3PmzAFg+PDhPPbYY0yfPp3rr7+eKVOmcO6553YOiwDceuutXHfddaxatYrPfvazneuXLFlCY2Mj06dPx90555xzWLduXY/HHjx4MGvWrOHrX/86ra2tVFZW8tOf/hSA6upqzjzzTL70pRzOZ5j6WM/p/qmurvb+2rx5c7/3LTbqS/EJpR/u+e/Lrl278lp/0tGjR/Na/+233+4rV67M6zFSjh496s3NzT5p0iRvb2/vsVx3zy+wxXvIVQ25iIicZo8//jizZs3i7rvvzuln5DXkIiIC3HHHHaftWJ///OdPeZM2F3SHLiID4vEXaCS3+vO8KtBFpN+GDBnCoUOHFOo55vF86EOGDOnTfhpyEZF+q6qqoqmpibfffjvvxzp+/HifA65YZdOX1F8s6gsFuoj026BBg/r0F3UGor6+vvPbnqUuX33RkIuISCAU6CIigVCgi4gEQoEuIhIIBbqISCAU6CIigVCgi4gEQoEuIhIIBbqISCAU6CIigVCgi4gEQoEuIhIIBbqISCAU6CIigVCgi4gEQoEuIhKIrALdzBaY2R4z229mt3Wz/aNmttnMtpnZy2b2mdw3VUREepMx0M2sHHgAuAq4CFhsZhelFfs2sNbdpwE3AP+Q64aKiEjvsrlDnwnsd/fX3P0DoA64Oq2MA2fGj0cCb+SuiSIikg3L9Ne6zWwRsMDdl8TLNwGz3H1ZosxYYBNwFjAM+LS7b+2mrqXAUoAxY8ZU19XV9avRLS0tDB8+vF/7Fhv1pfiE0g9QX4rVQPoyb968re5e0+1Gd+/1B1gEPJxYvgn4+7Qy3wL+NH48B9gFlPVWb3V1tffX5s2b+71vsVFfik8o/XBXX4rVQPoCbPEecjWbIZdmYHxiuSpel/RlYG38C+LnwBDg7CzqFhGRHMkm0BuASWY20cwGE73puT6tzOvApwDM7HeIAv3tXDZURER6lzHQ3b0NWAZsBHYTfZplp5ndaWYL42J/CnzFzHYATwA3xy8NRETkNKnIppC7bwA2pK37TuLxLmBubpsmIiJ9oW+KiogEQoEuIhIIBbqISCAU6CIigVCgi4gEQoEuIhIIBbqISCAU6CIigVCgi4gEQoEuIhIIBbqISCAU6CIigVCgi4gEQoEuIhIIBbqISCAU6CIigVCgi4gEQoEuIhIIBbqISCAU6CIigVCgi4gEQoEuIhIIBbqISCAU6CIigVCgi4gEQoEuIhIIBbqISCAU6CIigVCgi4gEQoEuIhIIBbqISCAU6CIigVCgi4gEQoEuIhKIrALdzBaY2R4z229mt/VQ5joz22VmO83s8dw2U0REMqnIVMDMyoEHgCuAJqDBzNa7+65EmUnACmCuu//GzM7NV4NFRKR72dyhzwT2u/tr7v4BUAdcnVbmK8AD7v4bAHc/mNtmiohIJubuvRcwWwQscPcl8fJNwCx3X5Yosw7YC8wFyoE73P3H3dS1FFgKMGbMmOq6urp+NbqlpYXhw4f3a99io74Un1D6AepLsRpIX+bNm7fV3Wu625ZxyCVLFcAkoBaoAp41s0vc/XCykLuvAlYB1NTUeG1tbb8OVl9fT3/3LTbqS/EJpR+gvhSrfPUlmyGXZmB8YrkqXpfUBKx39xPu/guiu/VJuWmiiIhkI5tAbwAmmdlEMxsM3ACsTyuzjujuHDM7G7gAeC13zRQRkUwyBrq7twHLgI3AbmCtu+80szvNbGFcbCNwyMx2AZuB5e5+KF+NFhGRU2U1hu7uG4ANaeu+k3jswLfiHxERKQB9U1REJBAKdBGRQCjQRUQCoUAXEQmEAl1EJBAKdBGRQCjQRUQCoUAXEQmEAl1EJBAKdBGRQCjQRUQCoUAXEQmEAl1EJBAKdBGRQCjQRUQCoUAXEQmEAl1EJBAKdBGRQCjQRUQCoUAXEQmEAl1EJBAKdBGRQCjQRUQCoUAXEQmEAl1EJBAKdBGRQCjQRUQCoUAXEQmEAl1EJBAKdBGRQCjQRUQCoUAXEQmEAl1EJBAKdBGRQCjQRUQCkVWgm9kCM9tjZvvN7LZeyv2hmbmZ1eSuiSIiko2MgW5m5cADwFXARcBiM7uom3IjgFuAF3LdSBERySybO/SZwH53f83dPwDqgKu7KfeXwHeB4zlsn4iIZMncvfcCZouABe6+JF6+CZjl7ssSZaYDf+7uf2hm9cCt7r6lm7qWAksBxowZU11XV9evRre0tDB8+PB+7Vts1JfiE0o/QH0pVgPpy7x587a6e7fD2hUDahVgZmXA94CbM5V191XAKoCamhqvra3t1zHr6+vp777FRn0pPqH0A9SXYpWvvmQz5NIMjE8sV8XrUkYAFwP1ZtYIzAbW641REZHTK5tAbwAmmdlEMxsM3ACsT2109yPufra7T3D3CcDzwMLuhlxERCR/Mga6u7cBy4CNwG5grbvvNLM7zWxhvhsoIiLZyWoM3d03ABvS1n2nh7K1A2+WiIj0lb4pKiISCAW6iEggFOgiIoFQoIuIBEKBLiISiAF/U1RkINZta2blxj28cbiVcaMqWT5/MtdMO6/QzZIs6fwVFwW6FMy6bc2seOoVWk+0A9B8uJUVT70CoFAoATp/xUdDLlIwKzfu6QyDlNYT7azcuKdALZK+0PkrPgp0KZg3Drf2ab0UF52/4qNAl4IZN6qyT+uluOj8FR8FuhTM8vmTqRxU3mVd5aByls+fXKAWSV/o/BUfvSkqBZN640yfkihNOn/FR4EuBXXNtPMUACVM56+4aMhFRCQQCnQRkUAo0EVEAqFAFxEJhAJdRCQQCnQRkUAo0EVEAqFAFxEJhAJdRCQQCnQRkUAo0EVEAqFAFxEJhAJdRCQQCnQRkUAo0EVEAqFAFxEJhAJdRCQQCnQRkUAo0EVEAqG/KSoindydDocOd9o7HHdod6fDnY6OaF27x+s7Uuvj8qkyyX1TZeJ62ztOlulw6Ii3R+VIPE6tp7P8f71+ggPP/7KzHV3KxMfuSLW38xipdpMon3bszvJ06WeqPck2JfdtTxyza7vp9nk4Wd75g48ZtXk4fwp0KRqnXASJC7Kjm3BIBkxPZaK6Tg2MKJS6D4NX32jj1y81nXJhdiQu/G7rTB3rlLAi0aauIdVdECX3SR03GQre2Ye0dqSHmDvvvXecwT9/pts6k8dN1lPUdr2aVbEyg/Iyo8yin+gxlJUZ5WaUxcvlZli8vbzMsHhd575ldKmjPF5XUVHWWSZVt5l17mvx8bscK1H+3Pa38vL0ZBXoZrYA+FugHHjY3f8qbfu3gCVAG/A28D/d/Zc5bmvOpe5Guv62j9advGhOvRvJFDRdL/5T71RSF5a78/Kv2ji6443OC7m7C/PUNma+I+kSZr0ETPodRG/PSU93VqnwbWl5jzMaNne9a0v2K7VcCkHy8o5+75q6kK1LqBBf2KmLmsTF3jVMyuL16WXKLAqXivIyzqjoGkrJ0EgGy8G33mLc2LO71NkloOLjltnJ8CmP22eJsOoSjKmgSvQpvd7UPieDsWtbk89L6nF6CJ58LqL1L/z858ydO/fkPqljp8rYyXrMLIf/GXKvvv6dvNSbMdDNrBx4ALgCaAIazGy9u+9KFNsG1Lj7e2b2x8BfA9fno8Frtxzg7/7zPYa89B9dQrOjg87wak887vElURy0RWH7tgFXkX4xd96NJEKky11L2gXQ9W4j7WJNuyuJjpF2vDJ4p6yVsb816tTjZbgb6vmC7yYkEmGZDKVkUCb3sfTAOCWUTt6FpY6xpaGBObNndRsyyecjPYhTz0cxqa+vp7Z2SqGbkROjhpRxzogzCt2MopbNHfpMYL+7vwZgZnXA1UBnoLv75kT554Ev5LKRSSMrBzFmWBljzh1+6kuctLubk8GSWE7cXXS5A+lysRrGqWGUse5uXm71fIxon61btjB71owuYdr1zq73/VPhVgxBEoXHtEI3Y8CahpVx/uhhhW6GSJ+ZZ7hNNbNFwAJ3XxIv3wTMcvdlPZT/e+BX7n5XN9uWAksBxowZU11XV9evRre0tDB8+PB+7Vts1JfiE0o/QH0pVgPpy7x587a6e01323L6pqiZfQGoAT7R3XZ3XwWsAqipqfHa2tp+HSe6E+zfvsVGfSk+ofQD1Jdila++ZBPozcD4xHJVvK4LM/s08OfAJ9z9/dw0T0REspXNF4sagElmNtHMBgM3AOuTBcxsGvAQsNDdD+a+mSIikknGQHf3NmAZsBHYDax1951mdqeZLYyLrQSGAz8ws+1mtr6H6kREJE+yGkN39w3AhrR130k8/nSO2yXSL+u2NbNy4x7eONzKuFGVLJ8/GeCUdddMO++0HDsfx8nGt9e9whMvHOCbF5/gyys2sHjWeO665pKCtEVOH31TVIKxblszK556hdYT7QA0H25l+Q92gMGJdu9ct+KpVwByGrbdHTsfx8nGt9e9wmPPv9653O7euaxQD5sm55JgrNy4pzNQU050eGeYp7SeaGflxj15P3Y+jpONJ1440Kf1Eg4FugTjjcOteSk7kPpyfZxstPfw3ZKe1ks4FOgSjHGjKvNSdiD15fo42Sjv4VvDPa2XcCjQJRjL50+mclB5l3WDyoxB5V2DrHJQeeebpfk8dj6Ok43Fs8b3ab2EQ2+KSjBSbz4W4lMuPR27EJ9ySb3xmRozLzfTp1w+JBToEpRrpp3XbYiejmDt6diFcNc1l3DXNZdQX1/Pf99YW+jmyGmiIRcRkUAo0EVEAqFAFxEJhAJdRCQQCnQRkUAo0EVEAqFAFxEJhAJdRCQQCnQRkUAo0EVEAqFAFxEJhAJdRCQQCnQRkUAo0EVEAqFAFxEJhAJdRCQQCnQRkUAo0EVEAqFAFxEJhAJdRCQQCnQRkUAo0EVEAqFAFxEJhAJdRCQQCnQRkUAo0EVEAqFAFxEJhAJdRCQQWQW6mS0wsz1mtt/Mbutm+xlmtibe/oKZTch5S0VEpFcZA93MyoEHgKuAi4DFZnZRWrEvA79x998G7gW+m+uGiohI77K5Q58J7Hf319z9A6AOuDqtzNXA9+PHTwKfMjPLXTNFRCSTiizKnAccSCw3AbN6KuPubWZ2BBgNvJMsZGZLgaXxYouZ7elPo4Gz0+suYepL8QmlH6C+FKuB9OX8njZkE+g54+6rgFUDrcfMtrh7TQ6aVHDqS/EJpR+gvhSrfPUlmyGXZmB8YrkqXtdtGTOrAEYCh3LRQBERyU42gd4ATDKziWY2GLgBWJ9WZj3wR/HjRcC/u7vnrpkiIpJJxiGXeEx8GbARKAdWu/tOM7sT2OLu64FHgH81s/3Ar4lCP58GPGxTRNSX4hNKP0B9KVZ56YvpRlpEJAz6pqiISCAU6CIigSj6QDezIWb2opntMLOdZvYX8fqJ8TQD++NpBwYXuq3ZMLNyM9tmZj+Kl0u1H41m9oqZbTezLfG6j5jZT8xsX/zvWYVuZzbMbJSZPWlm/2Vmu81sTin2xcwmx+cj9XPUzL5Zon35k/h6f9XMnohzoFSvlVvifuw0s2/G6/JyToo+0IH3gU+6+xRgKrDAzGYTTS9wbzzdwG+Iph8oBbcAuxPLpdoPgHnuPjXxedrbgGfcfRLwTLxcCv4W+LG7XwhMITo/JdcXd98Tn4+pQDXwHvBvlFhfzOw84BtAjbtfTPRhjBsowWvFzC4GvkL0jfspwOfM7LfJ1zlx95L5AYYCLxF9U/UdoCJePwfYWOj2ZdH+qvjkfRL4EWCl2I+4rY3A2Wnr9gBj48djgT2FbmcW/RgJ/IL4AwKl3Je09l8JPFeKfeHkN88/QvRJvB8B80vxWgH+B/BIYvn/AH+Wr3NSCnfoqWGK7cBB4CfAfwOH3b0tLtJE9J+g2N1HdDI74uXRlGY/ABzYZGZb4ykdAMa4+5vx418BYwrTtD6ZCLwN/HM8FPawmQ2jNPuSdAPwRPy4pPri7s3APcDrwJvAEWArpXmtvAr8npmNNrOhwGeIvoSZl3NSEoHu7u0evYysInrpcmFhW9R3ZvY54KC7by10W3LkMnefTjQL59fM7PLkRo9uPUrhM7EVwHTgQXefBrxL2svfEuoLAPHY8kLgB+nbSqEv8Xjy1US/bMcBw4AFBW1UP7n7bqKhok3Aj4HtQHtamZydk5II9BR3PwxsJnq5NSqeZgC6n46g2MwFFppZI9GMlZ8kGrsttX4AnXdRuPtBonHamcBbZjYWIP73YOFamLUmoMndX4iXnyQK+FLsS8pVwEvu/la8XGp9+TTwC3d/291PAE8RXT+leq084u7V7n450dj/XvJ0Too+0M3sHDMbFT+uBK4getNqM9E0AxBNO/DDgjQwS+6+wt2r3H0C0cvhf3f3GymxfgCY2TAzG5F6TDRe+ypdp4Aoib64+6+AA2Y2OV71KWAXJdiXhMWcHG6B0uvL68BsMxsaT8OdOicld60AmNm58b8fBf4AeJw8nZOi/6aomV1KNNd6OdEvoLXufqeZfYzoTvcjwDbgC+7+fuFamj0zqwVudffPlWI/4jb/W7xYATzu7neb2WhgLfBR4JfAde7+6wI1M2tmNhV4GBgMvAZ8ifj/GqXXl2FEgfgxdz8Sryu58xJ/PPl6oI3oulhCNGZeUtcKgJn9jOj9shPAt9z9mXydk6IPdBERyU7RD7mIiEh2FOgiIoFQoIuIBEKBLiISCAW6iEggTusfiRbJVvyxrmfixd8i+nbd2/HyTHf/oCAN60b8MdQP3P3/Fbgp8iGnQJei5O6HiGbXxMzuAFrc/Z5CtcfMKhLziKSrBVqArAM9Q30i/aIhFykZZlZtZv8RTwi2MfHV6Xozu9fMtsTzmc8ws6fiuabvistMiOc7/79xmSfjyZIy1XufRfO932Jmvx/Px73NzH5qZmPMbALwv4A/iecg/z0ze9TMFiXa3RL/W2tmPzOz9cCueNK5lWbWYGYvm9lXT+sTKsFRoEupMOB+YJG7VwOrgbsT2z/waF72fyT6GvXXgIuBm+PhG4DJwD+4++8AR4H/bWaDMtQ72N1r3P1vgP8EZseTeNUBf+bujfEx7/VoLvKfZejHdOAWd7+AaD7vI+4+A5gBfMXMJvb9qRGJaMhFSsUZRAH9k2h6D8qJplZNWR//+wqwMzU1qZm9RjRd6WHggLs/F5d7jOiPKPw4Q71rEo+rgDXxHfxgonnU++pFd0/tdyVwaeJufiQwqZ/1iijQpWQYUVDP6WF7ak6PjsTj1HLq/3n6PBeeRb3vJh7fD3zP3dfHb4Te0cM+bcSvfs2sjCj8u6vPgK+7+8Ye6hHpEw25SKl4HzjHzOYAmNkgM/vdPtbx0dT+wOeJhlD29KHekZycsvWPEuuPASMSy41EfwIOonnJB/VQ30bgj+NhH8zsgnhyLZF+UaBLqeggmjr1u2a2g+gPBXy8j3XsIfpjHLuBs4j+qMUHfaj3DuAHZraV6M+hpTwNXJt6UxT4J+ATcX1z6HpXnvQw0bSwL5nZq8BD6FWzDIBmW5QPhfjTKD/y6I8OiwRJd+giIoHQHbqISCB0hy4iEggFuohIIBToIiKBUKCLiARCgS4iEoj/D4/WYtJsTILdAAAAAElFTkSuQmCC\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[['Intercept','Temperature']])\n", "data_pred.plot(x=\"Temperature\",y=\"Frequency\",kind=\"line\",ylim=[0,1])\n", "plt.scatter(x=data[\"Temperature\"],y=data[\"Frequency\"])\n", "plt.grid(True)" ] }, { "cell_type": "markdown", "metadata": { "hideCode": false, "hidePrompt": false, "scrolled": true }, "source": [ "As expected from the initial data, the\n", "temperature has no significant impact on the probability of failure of the\n", "O-rings. It will be about 0.2, as in the tests\n", "where we had a failure of at least one joint. Let's get back\n", "to the initial dataset to estimate the probability of failure:" ] }, { "cell_type": "code", "execution_count": 6, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "0.06521739130434782\n" ] } ], "source": [ "data = pd.read_csv(\"../data/shuttle.csv\")\n", "print(np.sum(data.Malfunction)/np.sum(data.Count))" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "This probability is thus about $p=0.065$. Knowing that there is\n", "a primary and a secondary O-ring on each of the three parts of the\n", "launcher, the probability of failure of both joints of a launcher\n", "is $p^2 \\approx 0.00425$. The probability of failure of any one of the\n", "launchers is $1-(1-p^2)^3 \\approx 1.2%$. That would really be\n", "bad luck.... Everything is under control, so the takeoff can happen\n", "tomorrow as planned.\n", "\n", "But the next day, the Challenger shuttle exploded and took away\n", "with her the seven crew members. The public was shocked and in\n", "the subsequent investigation, the reliability of the\n", "O-rings was questioned. Beyond the internal communication problems\n", "of NASA, which have a lot to do with this fiasco, the previous analysis\n", "includes (at least) a small problem.... Can you find it?\n", "You are free to modify this analysis and to look at this dataset\n", "from all angles in order to to explain what's wrong." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Comments" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "All the statistics look reasonable at a first glance. However, the sample is quite small (22 observations) and the considered temperatures are quite far from the launch temperature (30 °F), which can lead to non accurate predictions. In addition, excluding the observations without failures makes it difficult for the model to learn how the temperature impacts on having a failure or not, and not only on the number of failures. Pressure can also be an important predictor that is not taken into account for the analysis.\n", "\n", "The same analysis can therefore be repeated, considering all observations and pressure as an additional predictor." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## The new analysis" ] }, { "cell_type": "code", "execution_count": 7, "metadata": { "scrolled": true }, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAA3gAAAFBCAYAAAAlhA0CAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjQuMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/MnkTPAAAACXBIWXMAAAsTAAALEwEAmpwYAAAlmElEQVR4nO3de5hddX3v8fc3F0JCUkCigZIgKEjLEYoQEQ5ag7eCniZtoRZPFW1FejzSKrVVvDyUQ2ufgi20Wqri5VjUGhFajJUWRI1Yj1wChABBcOSWBAgYAyQQcmG+54+9hu5MZ5I9k732rLX2+/U882Svtdes/f2tffnms9dlIjORJEmSJNXfpIkuQJIkSZLUHQY8SZIkSWoIA54kSZIkNYQBT5IkSZIawoAnSZIkSQ1hwJMkSZKkhigt4EXEFyLi0Yi4Y5T7IyI+EREDEbEiIo4qqxZJkqrEHilJKkuZe/C+CJy4g/tPAg4pfs4APlViLZIkVckXsUdKkkpQWsDLzOuAn+9gkUXApdlyPbBXROxXVj2SJFWFPVKSVJaJPAdvf2BV2/TqYp4kSf3OHilJGpcpE11AJyLiDFqHqDB9+vSj582bN8EV7dzg4CCTJvXfNWwcd//oxzGD4+61e+6552eZ+fyeP3CN1KlH1vX9Y929Zd29V9fa+7nuHfXHiQx4a4D2LjS3mPdfZOYlwCUA8+fPz2XLlpVf3S5aunQpCxYsmOgyes5x949+HDM47l6LiAd6/qDV0MgeWdf3j3X3lnX3Xl1r7+e6d9QfJzLyLgFOK64UdizwRGY+PIH1SJJUFfZISdK4lLYHLyK+CiwAZkfEauDPgKkAmflp4CrgjcAA8DTwe2XVIklSldgjJUllKS3gZeZbdnJ/Au8p6/ElSaoqe6QkqSz1OytRkiRJkjQiA54kSZIkNYQBT5IkSZIawoAnSZIkSQ1hwJMkSZKkhjDgSZIkSVJDGPAkSZIkqSEMeJIkSZLUEAY8SZIkSWoIA54kSZIkNYQBT5IkSZIawoAnSZIkSQ1hwJMkSZKkhjDgSZIkSVJDGPAkSZIkqSEMeJIkSZLUEAY8SZIkSWoIA54kSZIkNYQBT5IkSZIawoAnSZIkSQ1hwJMkSZKkhjDgSZIkSVJDGPAkSZIkqSEMeJIkSZLUEAY8SZIkSWoIA54kSZIkNYQBT5IkSZIawoAnSZIkSQ1hwJMkSZKkhjDgSZIkSVJDGPAkSZIkqSEMeJIkSZLUEAY8SZIkSWoIA54kSZIkNYQBT5IkSZIawoAnSZIkSQ1hwJMkSZKkhjDgSZIkSVJDGPAkSZIkqSEMeJIkSZLUEAY8SZIkSWoIA54kSZIkNYQBT5IkSZIawoAnSZIkSQ1hwJMkSZKkhjDgSZIkSVJDlBrwIuLEiLg7IgYi4uwR7j8gIr4XEbdGxIqIeGOZ9UiSVAX2R0lSWUoLeBExGbgYOAk4DHhLRBw2bLGPApdl5suAU4F/KKseSZKqwP4oSSpTmXvwjgEGMvPezNwCLAYWDVsmgV8obu8JPFRiPZIkVYH9UZJUmsjMclYccQpwYmaeXky/DXhFZp7Ztsx+wDXA3sAewOsy8+YR1nUGcAbAnDlzjl68eHEpNXfTxo0bmTlz5kSX0XOOu3/045jBcffaCSeccHNmzu/5A5eom/2xWLY2PbKu7x/r7i3r7r261t7Pde+oP07ZpTXvurcAX8zMv4mI44AvRcRLM3OwfaHMvAS4BGD+/Pm5YMGC3lc6RkuXLqUOdXab4+4f/ThmcNzqmY76I9SrR9b1dWTdvWXdvVfX2q17ZGUeorkGmNc2PbeY1+6dwGUAmfkjYHdgdok1SZI00eyPkqTSlBnwbgIOiYiDImI3WieJLxm2zIPAawEi4pdpNbDHSqxJkqSJZn+UJJWmtICXmduAM4GrgbtoXQ3szog4LyIWFou9H3hXRNwGfBV4R5Z1UqAkSRVgf5QklanUc/Ay8yrgqmHzzmm7vRI4vswaJEmqGvujJKkspf6hc0mSJElS7xjwJEmSJKkhDHiSJEmS1BAGPEmSJElqCAOeJEmSJDWEAU+SJEmSGsKAJ0mSJEkNYcCTJEmSpIYw4EmSJElSQxjwJEmSJKkhDHiSJEmS1BAGPEmSJElqCAOeJEmSJDWEAU+SJEmSGsKAJ0mSJEkNYcCTJEmSpIYw4EmSJElSQxjwJEmSJKkhDHiSJEmS1BAGPEmSJElqCAOeJEmSJDWEAU+SJEmSGsKAJ0mSJEkNYcCTJEmSpIYw4EmSJElSQxjwJEmSJKkhDHiSJEmS1BAGPEmSJElqCAOeJEmSJDWEAU+SJEmSGsKAJ0mSJEkNYcCTJEmSpIYw4EmSJElSQxjwJEmSJKkhDHiSJEmS1BAGPEmSJElqCAOeJEmSJDWEAU+SJEmSGsKAJ0mSJEkNYcCTJEmSpIYw4EmSJElSQxjwJEmSJKkhDHiSJEmS1BAGPEmSJElqCAOeJEmSJDVEqQEvIk6MiLsjYiAizh5lmTdHxMqIuDMi/qnMeiRJqgL7oySpLFPKWnFETAYuBl4PrAZuioglmbmybZlDgA8Bx2fm+oh4QVn1SJJUBfZHSVKZytyDdwwwkJn3ZuYWYDGwaNgy7wIuzsz1AJn5aIn1SJJUBfZHSVJpygx4+wOr2qZXF/PavQR4SUT8MCKuj4gTS6xHkqQqsD9KkkoTmVnOiiNOAU7MzNOL6bcBr8jMM9uW+VdgK/BmYC5wHXB4Zj4+bF1nAGcAzJkz5+jFixeXUnM3bdy4kZkzZ050GT3nuPtHP44ZHHevnXDCCTdn5vyeP3CJutkfi2Vr0yPr+v6x7t6y7t6ra+39XPeO+mNp5+ABa4B5bdNzi3ntVgM3ZOZW4L6IuAc4BLipfaHMvAS4BGD+/Pm5YMGCsmrumqVLl1KHOrvNcfePfhwzOG51Rdf6I9SrR9b1dWTdvWXdvVfX2q17ZGUeonkTcEhEHBQRuwGnAkuGLXMlsAAgImbTOiTl3hJrkiRpotkfJUml6SjgRcThY11xZm4DzgSuBu4CLsvMOyPivIhYWCx2NbAuIlYC3wP+NDPXjfWxJEmaKGPtkfZHSVKZOj1E8x8iYhrwReArmflEJ7+UmVcBVw2bd07b7QT+uPiRJKmOxtwj7Y+SpLJ0tAcvM18F/C6tcwZujoh/iojXl1qZJEk1YI+UJFVJx+fgZeZPgI8CHwReDXwiIn4cEb9VVnGSJNWBPVKSVBWdnoN3RERcROtcgdcAv56Zv1zcvqjE+iRJqjR7pCSpSjo9B++TwOeAD2fmpqGZmflQRHy0lMokSaoHe6QkqTI6DXhvAjZl5rMAETEJ2D0zn87ML5VWnSRJ1WePlCRVRqfn4F0LTG+bnlHMkySp39kjJUmV0WnA2z0zNw5NFLdnlFOSJEm1Yo+UJFVGpwHvqYg4amgiIo4GNu1geUmS+oU9UpJUGZ2eg/c+4OsR8RAQwL7A75RVlCRJNfI+7JGSpIroKOBl5k0R8UvAocWsuzNza3llSZJUD/ZISVKVdLoHD+DlwIHF7xwVEWTmpaVUJUlSvdgjJUmV0FHAi4gvAS8GlgPPFrMTsHlJkvqaPVKSVCWd7sGbDxyWmVlmMZIk1ZA9UpJUGZ1eRfMOWieNS5Kk7dkjJUmV0ekevNnAyoi4Edg8NDMzF5ZSlSRJ9WGPlCRVRqcB79wyi5AkqcbOnegCJEka0umfSfh+RLwQOCQzr42IGcDkckuTJKn67JGSpCrp6By8iHgXcDnwmWLW/sCVJdUkSVJt2CMlSVXS6UVW3gMcDzwJkJk/AV5QVlGSJNWIPVKSVBmdBrzNmbllaCIiptD6Gz+SJPU7e6QkqTI6DXjfj4gPA9Mj4vXA14FvlleWJEm1YY+UJFVGpwHvbOAx4HbgD4CrgI+WVZQkSTVij5QkVUanV9EcBD5b/EiSpII9UpJUJR0FvIi4jxHOJ8jMF3W9IkmSasQeKUmqkk7/0Pn8ttu7A78NPK/75UiSVDv2SElSZXR0Dl5mrmv7WZOZfwu8qdzSJEmqPnukJKlKOj1E86i2yUm0vq3sdO+fJEmNZY+UJFVJpw3ob9pubwPuB97c9WokSaofe6QkqTI6vYrmCWUXIklSHdkjJUlV0ukhmn+8o/sz88LulCNJUr3YIyVJVTKWq2i+HFhSTP86cCPwkzKKkiSpRuyRkqTK6DTgzQWOyswNABFxLvCtzHxrWYVJklQT9khJUmV09GcSgDnAlrbpLcU8SZL6nT1SklQZne7BuxS4MSL+pZj+DeAfS6lIkqR6sUdKkiqj06tofiwi/g14VTHr9zLz1vLKkiSpHuyRkqQq6fQQTYAZwJOZ+XfA6og4qKSaJEmqG3ukJKkSOgp4EfFnwAeBDxWzpgJfLqsoSZLqwh4pSaqSTvfg/SawEHgKIDMfAmaVVZQkSTVij5QkVUanAW9LZiaQABGxR3klSZJUK/ZISVJldBrwLouIzwB7RcS7gGuBz5ZXliRJtWGPlCRVxk6vohkRAXwN+CXgSeBQ4JzM/HbJtUmSVGn2SElS1ew04GVmRsRVmXk4YMOSJKlgj5QkVU2nh2jeEhEvL7USSZLqyR4pSaqMjv7QOfAK4K0RcT+tq4QFrS8ujyirMEmSasIeKUmqjB0GvIg4IDMfBH6tR/VIklQL9khJUhXt7BDNKwEy8wHgwsx8oP1nZyuPiBMj4u6IGIiIs3ew3MkRkRExf0zVS5I0ca6E8fVI+6MkqSw7C3jRdvtFY1lxREwGLgZOAg4D3hIRh42w3CzgvcANY1m/JEkTbFw90v4oSSrTzgJejnK7E8cAA5l5b2ZuARYDi0ZY7s+B84Fnxrh+SZIm0nh7pP1RklSanQW8X4mIJyNiA3BEcfvJiNgQEU/u5Hf3B1a1Ta8u5j0nIo4C5mXmt8ZcuSRJE2u8PdL+KEkqzQ4vspKZk8t64IiYBFwIvKODZc8AzgCYM2cOS5cuLausrtm4cWMt6uw2x90/+nHM4Lj1n8rqkWPpj8XytemRdX0dWXdvWXfv1bV26x5Zp38mYTzWAPPapucW84bMAl4KLI0IgH2BJRGxMDOXta8oMy8BLgGYP39+LliwoMSyu2Pp0qXUoc5uc9z9ox/HDI5bXdG1/gj16pF1fR1Zd29Zd+/VtXbrHlmnf+h8PG4CDomIgyJiN+BUYMnQnZn5RGbOzswDM/NA4HpgxOYlSVKD2B8lSaUpLeBl5jbgTOBq4C7gssy8MyLOi4iFZT2uJElVZn+UJJWpzEM0ycyrgKuGzTtnlGUXlFmLJElVYX+UJJWlzEM0JUmSJEk9ZMCTJEmSpIYw4EmSJElSQxjwJEmSJKkhDHiSJEmS1BAGPEmSJElqCAOeJEmSJDWEAU+SJEmSGsKAJ0mSJEkNYcCTJEmSpIbou4C3buNmblv1OOs2bp7oUiTVzLqNm9m09Vk/PyRJ0rgMrN3A+qe3MrB2Q2mP0VcB7xvL13D8+d/lrZ+7gePP/y5Llq+Z6JIk1cTQ58d9jz3l54ckSRqzc668nddddB2r1z/N6y66jnO+cXspj9M3AW/dxs188IoVPLN1kA2bt/HM1kE+cMUKv4mXtFPtnx/PZvr5IUmSxmRg7QYuvf7B7eZd+qMHS9mT1zcBb/X6TUydtP1wp06axOr1myaoIkl14eeHJEnaFctXPT6m+buibwLe3L2ns3VwcLt5WwcHmbv39AmqSFJd+PkhSZJ2xZHz9hrT/F3RNwFvn5nTuODkI9h96iRmTZvC7lMnccHJR7DPzGkTXZqkimv//Jgc4eeHJEkak4PnzOK04w7Ybt5pxx3AwXNmdf2xpnR9jRW28Mj9Of7g2axev4m5e0/3P2eSOjb0+XHjj/6DHy58pZ8fkiRpTM5bdDinHXsgt998PdeedWwp4Q76LOBB65t4/2MmaTz2mTmN6VMn+xkiSZLG5eA5s1g9Y2pp4Q766BBNSZIkSWo6A54kSZIkNYQBT5IkSZIawoAnSZIkSQ1hwJMkSZKkhjDgSZIkSVJDGPAkSZIkqSEMeJIkSZLUEAY8SZIkSWoIA54kSZIkNYQBT5IkSZIawoAnSZIkSQ1hwJMkSZKkhjDgSZIkSVJDGPAkSZIkqSEMeJIkSZLUEAY8SZIkSWoIA54kSZIkNYQBT5IkSZIawoAnSZIkSQ1hwJMkSZKkhjDgSZIkSVJDGPAkSZIkqSEMeJIkSZLUEAY8SZIkSWoIA54kSZIkNYQBT5IkSZIawoAnSZIkSQ1RasCLiBMj4u6IGIiIs0e4/48jYmVErIiI70TEC8usR5KkKrA/SpLKUlrAi4jJwMXAScBhwFsi4rBhi90KzM/MI4DLgQvKqkeSpCqwP0qSylTmHrxjgIHMvDcztwCLgUXtC2Tm9zLz6WLyemBuifVIklQF9kdJUmkiM8tZccQpwImZeXox/TbgFZl55ijL/z3wSGb+xQj3nQGcATBnzpyjFy9eXErN3bRx40Zmzpw50WX0nOPuH/04ZnDcvXbCCSfcnJnze/7AJepmfyzur02PrOv7x7p7y7p7r66193PdO+qPU3ZpzV0SEW8F5gOvHun+zLwEuARg/vz5uWDBgt4VN05Lly6lDnV2m+PuH/04ZnDc6q2d9UeoV4+s6+vIunvLunuvrrVb98jKDHhrgHlt03OLeduJiNcBHwFenZmbS6xHkqQqsD9KkkpT5jl4NwGHRMRBEbEbcCqwpH2BiHgZ8BlgYWY+WmItkiRVhf1RklSa0gJeZm4DzgSuBu4CLsvMOyPivIhYWCz2cWAm8PWIWB4RS0ZZnSRJjWB/lCSVqdRz8DLzKuCqYfPOabv9ujIfX5KkKrI/SpLKUuofOpckSZIk9Y4Br6bWbdzMbaseZ93Gzs67H+vyTdGv4y7bwNoNrH96KwNrN0x0KZK6yM9MSSrXd1Y+wpr1m/jOykdKe4xK/JkEjc03lq/hg1esYOqkSWwdHOSCk49g4ZH7d235pujXcZftnCtv59LrH+T9h2/jrIuu47TjDuC8RYdPdFmSdpGfmZJUrjdctJR71j7F+w/fxkcuvZlD5+zB1Wct6PrjuAevZtZt3MwHr1jBM1sH2bB5G89sHeQDV6wY9dvWsS7fFP067rINrN3Apdc/uN28S3/0oHvypJrzM1OSyvWdlY9wz9qntpt399qnStmTZ8CrmdXrNzF10vZP29RJk1i9flNXlm+Kfh132ZavenxM8yXVg5+ZklSua1auHdP8XWHAq5m5e09n6+DgdvO2Dg4yd+/pXVm+Kfp13GU7ct5eY5ovqR78zJSkcr3hsDljmr8rDHg1s8/MaVxw8hHsPnUSs6ZNYfepk7jg5CPYZ+a0rizfFP067rIdPGcWpx13wHbzTjvuAA6eM2uCKpLUDX5mSlK5XnvYvhw6Z4/t5h06Zw9ee9i+XX8sL7JSQwuP3J/jD57N6vWbmLv39J024LEu3xT9Ou6ynbfocE479kBuv/l6rj3rWMOd1BB+ZkpSua4+awHfWfkIj9x9K58/7WWlhDsw4NXWPjOnjan5jnX5pujXcZft4DmzWD1jquFOahg/MyWpXK89bF+WPjqdBSWFO/AQTUmSJElqDAOeJEmSJDWEAU+SJEmSGsKAJ0mSJEkNYcCTJEmSpIYw4EmSJElSQxjwJEmSJKkhDHiSJEmS1BAGPEmSJElqCAOeJEmSJDWEAU+SJEmSGsKAJ0mSJEkNYcCTJEmSpIYw4EmSJElSQxjwJEmSJKkhDHiSJEmS1BAGPEmSJElqCAOeJEmSJDWEAU+SJEmSGsKAJ0mSJEkNYcCTJEmSpIYw4EmSJElSQxjwJEmSJKkhDHiSJEmS1BAGPEmSJElqCAOeJEmSJDWEAU+SJEmSGsKAJ0mSJEkNYcCTJEmSpIYw4EmSJElSQxjwJEmSJKkhDHiSJEmS1BAGPEmSJElqCAOeJEmSJDWEAU+SJEmSGsKAJ0mSJEkNUWrAi4gTI+LuiBiIiLNHuH9aRHytuP+GiDiwzHokSaoC+6MkqSylBbyImAxcDJwEHAa8JSIOG7bYO4H1mXkwcBFwfln1SFWxbuNmblv1OOs2bu5o+WX3rePCa+5m2X3rSnuMsS4/sHYD65/eysDaDR3XNFZlj2G8NW3a+uyYttPly1ZVajuN93fKrqmfVLk/XnnLKk7/x5u48pZVvXg4lagXn9OSxu5j37yDHz+ygY99847SHmNKaWuGY4CBzLwXICIWA4uAlW3LLALOLW5fDvx9RERmZol1SRPmG8vX8MErVjB10iS2Dg5ywclHsPDI/Udd/q2fu57/GGgFu098d4BXHbwPXzr92K4+xliXP+fK27n0+gd5/+HbOOui6zjtuAM4b9HhOxn52JQ9hl2p6Y9+eStnnf/djrfTkCpsp/H+Ttk19aFK9sdj//LbPPLkFgCuvetRzv/3H/OjD7++rIdTiXrxOS1p7F509rcYBN5/+CCf/eEDfP6HD3DvX72p649T5iGa+wPtXwGuLuaNuExmbgOeAPYpsSZpwqzbuJkPXrGCZ7YOsmHzNp7ZOsgHrlgx6l6OZfetey7cDfnBwLod7skb62OMdfmBtRu2Cy0Al/7owa5+Q1z2GHa1pmcza7mdxvs7ZdfUpyrXH6+8ZdVz4W7Iw09ucU9eDfXi80fS2H3sm3cwOGzeYDG/26KsLwMj4hTgxMw8vZh+G/CKzDyzbZk7imVWF9M/LZb52bB1nQGcUUweCtxdStHdNRv42U6Xah7HPYqYOn3GlL33e0lMmjR5aF4ODj67bf3D9+TWTU8PX37yrNm/OHmPvfYbPv/Zpx5/+NkNP3uoG48x1uUnzdhznym/8PwDAZ59+gkmz9gTgG1PPnb/4NNPdH4M6Q6UPYZdrWlo3J1up3YTuZ3G+zttuv4a79ALM/P54/zdSupmfyzu2+UeOWXvX3zxpGkz9ho+f3Dz049vW//QT8e6vh2oa4+oTd29+Jzugdps72HqWjfUt/ba1D31+QceHpOn7Abbvzfz2W1btj52/+3jWOWo/bHMQzTXAPPapucW80ZaZnVETAH2BP7Lh09mXgJcUlKdpYiIZZk5f6Lr6DXH3T8iYtm2Jx7tqzFDf4+7317jJepaf4R69ci6vo7qXHcdP6/qvL3rWDfUt/Y6113me7PMQzRvAg6JiIMiYjfgVGDJsGWWAG8vbp8CfNfz7yRJDWd/lCSVprQ9eJm5LSLOBK4GJgNfyMw7I+I8YFlmLgE+D3wpIgaAn9NqcpIkNZb9UZJUpjIP0SQzrwKuGjbvnLbbzwC/XWYNE6gWh8uUwHH3j34cMzhudUEf98e6vo6su7esu/fqWrt1j6C0i6xIkiRJknqrzHPwJEmSJEk9ZMDrkoi4PyJuj4jlEbGsmHduRKwp5i2PiDdOdJ3dFBF7RcTlEfHjiLgrIo6LiOdFxLcj4ifFv3tPdJ3dNsq4m/5cH9o2tuUR8WREvK/pz/cOxt305/usiLgzIu6IiK9GxO7FBUFuiIiBiPhacXEQaVR17RF1ev1HxBci4tHiz2oMzRtxG0fLJ4oxrIiIoypW98eL18qKiPiXiNir7b4PFXXfHRG/NiFFM3Ldbfe9PyIyImYX05Xe3sX8Pyy2+Z0RcUHb/Mpu74g4MiKuL3rvsog4pphfpe09LyK+FxEri2373mJ+796bmelPF36A+4HZw+adC/zJRNdW4pj/ETi9uL0bsBdwAXB2Me9s4PyJrrNH4270cz1s/JOBR4AX9sPzPcq4G/t80/oD2/cB04vpy4B3FP+eWsz7NPDuia7Vn2r/1LFH1O31D/wqcBRwR9u8Ebcx8Ebg34AAjgVuqFjdbwCmFLfPb6v7MOA2YBpwEPBTYHJV6i7mz6N10aQHhv4vWIPtfQJwLTCtmH5BHbY3cA1wUts2XlrB7b0fcFRxexZwT7Fde/bedA+exiUi9qT1xvs8QGZuyczHgUW0mjrFv78xEfWVZQfj7ievBX6amQ/Q8Od7mPZxN90UYHq0/v7aDOBh4DXA5cX9TX+utYtq3iNq8/rPzOtoXWW13WjbeBFwabZcD+wVEfv1pNBhRqo7M6/JzG3F5PW0/j4ktOpenJmbM/M+YAA4pmfFbl/jSNsb4CLgA0D7hS0qvb2BdwN/lZmbi2UeLeZXfXsn8AvF7T2Bh4rbVdreD2fmLcXtDcBdtL486tl704DXPQlcExE3R8QZbfPPLHa3fqGKh6LsgoOAx4D/GxG3RsTnImIPYE5mPlws8wgwZ8IqLMdo44bmPtfDnQp8tbjd9Oe7Xfu4oaHPd2auAf4aeJDWf2yfAG4GHm/7z9dqWs1KGk0te0RDXv+jbeP9gVVty1V5HL9Pa48GVLzuiFgErMnM24bdVem6gZcAryoOPf5+RLy8mF/1ut8HfDwiVtF6r36omF/JuiPiQOBlwA308L1pwOueV2bmUcBJwHsi4leBTwEvBo6k1Sj+ZuLK67optHabfyozXwY8RWt383Oytd+5aZdpHW3cTX6un1Ocd7IQ+Prw+xr6fAMjjruxz3cRVhfR+g/6LwJ7ACdOaFGqo1r2iKa9/qu4jXcmIj4CbAO+MtG17ExEzAA+DJyzs2UraArwPFqHBP4pcFlExMSW1JF3A2dl5jzgLIqjBKooImYCVwDvy8wn2+8r+71pwOuS4lu/oV3c/wIck5lrM/PZzBwEPssE7eIuyWpgdWbeUExfTquZrx3arVz8++gov19XI4674c91u5OAWzJzbTHd9Od7yHbjbvjz/Trgvsx8LDO3Av8MHE/rkJGhv506F1gzUQWqFuraI5rw+h9tG6+hda7YkMqNIyLeAfwP4HeL/wBDtet+Ma0vA26LiPtp1XZLROxLteuG1nv0n4vDAm8EBoHZVL/ut9N6X0LrS9eh/lupuiNiKq1w95XMHKq3Z+9NA14XRMQeETFr6DatE4XvGHb87G8C/+WqS3WVmY8AqyLi0GLWa4GVwBJabz6Kf78xAeWVZrRxN/m5HuYtbH+YYqOf7zbbjbvhz/eDwLERMaP4Nnfovf094JRimSY/1+qCGveIJrz+R9vGS4DTiiv2HQs80Xa42ISLiBNpnce2MDOfbrtrCXBqREyLiIOAQ4AbJ6LG4TLz9sx8QWYemJkH0gpNRxWv/0pvb+BKWhdaISJeQutCSD+jwtu78BDw6uL2a4CfFLcrs72Lz47PA3dl5oVtd/XuvTneq7P4s93Vcl5E64pDtwF3Ah8p5n8JuB1YUTx5+010rV0e95HAsmJ8VwJ7A/sA36H1hrsWeN5E19mjcTf6uS7GvQewDtizbV4/PN8jjbvRzzfwf4Af0wquX6J1NbUX0WryA7S+NZ020XX6U+2fuvaIOr3+aX3x9DCwlVa4eOdo25jWFfoupnVVxNuB+RWre4DWeUjLi59Pty3/kaLuuymuoFiVuofdfz//eRXNqm/v3YAvF6/zW4DX1GF7A6+kdV7sbbTOazu6gtv7lbQOv1zR9np+Yy/fm1GsWJIkSZJUcx6iKUmSJEkNYcCTJEmSpIYw4EmSJElSQxjwJEmSJKkhDHiSJEmS1BBTdr6IJICIGLq8LcC+wLPAY8X0MZm5ZUIKG0FELAC2ZOb/m+BSJEl6TkQ8S+tS8FOAu4C35/Z/907SLnIPntShzFyXmUdm5pHAp4GLhqYnItxFxI6+oFkA/Pcurk+SpG7YVPTNlwJbgP/Vfmcve5F9T01lwJN2QUQcHRHfj4ibI+LqiNivmL80Ii6KiGURcVdEvDwi/jkifhIRf1Esc2BE/DgivlIsc3lEzOhgvX8bEcuA90bEr0fEDRFxa0RcGxFzIuJAWg3zrIhYHhGviogvRsQpbXVvLP5dEBE/iIglwMqImBwRH4+ImyJiRUT8QU83qCSpn/wAOLjTXhQR+0XEdUVvu6Pob5OLHndHRNweEWcVyy6NiPnF7dkRcX9x+x0RsSQivgt8JyL2iIgvRMSNRS9dNDGbQuoev7mQxi+ATwKLMvOxiPgd4GPA7xf3b8nM+RHxXuAbwNHAz4GfRsRFxTKHAu/MzB9GxBeA/x0Rf7eT9e6WmUNNa2/g2MzMiDgd+EBmvj8iPg1szMy/LpZ75w7GcRTw0sy8LyLOAJ7IzJdHxDTghxFxTWbet+ubS5KklmLv2UnAvxezdtqLgN8Crs7Mj0XEZGAGcCSwf7FHkIjYq4OHPwo4IjN/HhF/CXw3M3+/+N0bI+LazHyqe6OVesuAJ43fNOClwLcjAmAy8HDb/UuKf28H7szMhwEi4l5gHvA4sCozf1gs92Xgj2g1ux2t92ttt+cCXyv28O0GjCeI3dgW4N4AHNG2t29P4JBxrleSpOGmR8Ty4vYPgM/TOqWgk150E/CFiJgKXJmZy4ue+qKI+CTwLeCaDmr4dmb+vO2xFkbEnxTTuwMH0Do/UKolA540fkEruB03yv2bi38H224PTQ+993LY72QH623/VvGTwIWZuaS4sMq5o/zONopDsiNiEq0wONL6AvjDzLx6lPVIkrQrNhXnsj+n+DKzo14UEb8KvAn4YkRcmJmXRsSvAL9G6/SEN9M64uW5vkcrtLUb/lgnZ+bd4x6RVDGegyeN32bg+RFxHEBETI2I/zbGdRww9PvA/wT+A7h7DOvdE1hT3H572/wNwKy26ftpHSIKsBCYOsr6rgbeXXw7SkS8JCL26Hw4kiTtshF7UUS8EFibmZ8FPgccFRGzgUmZeQXwUVqHX8L2fe8URnc18IdRpMyIeFnXRyP1mAFPGr9BWk3j/Ii4DVjOGK9cSSvMvSci7gL2Bj5VXJGz0/WeC3w9Im4GftY2/5vAbw5dZAX4LPDqYn3Hsf23l+0+B6wEbomIO4DP4J5+SVJvjdaLFgC3RcStwO8AfwfsDywtDvv8MvChYh1/TSsk3grM3sFj/TmtLz1XRMSdxbRUa5E5/AgxSb1QXO3yX4dODJckSZJ2lXvwJEmSJKkh3IMnSZIkSQ3hHjxJkiRJaggDniRJkiQ1hAFPkiRJkhrCgCdJkiRJDWHAkyRJkqSGMOBJkiRJUkP8f6Yn9B/Vuw/tAAAAAElFTkSuQmCC\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "# Reload data\n", "data = pd.read_csv(\"../data/shuttle.csv\")\n", "#Graphical inspection\n", "data[\"Frequency\"]=data.Malfunction/data.Count\n", "_, axs = plt.subplots(1, 2, figsize=(15, 5))\n", "data.plot(x=\"Temperature\",y=\"Frequency\",kind=\"scatter\",ylim=[0,1],ax=axs[0])\n", "data.plot(x=\"Pressure\",y=\"Frequency\",kind=\"scatter\",ylim=[0,1],ax=axs[1])\n", "axs[0].grid(True)\n", "axs[1].grid(True)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Only by graphical inspection, it is clear that temperature has indeed a big impact on failures: they are unlikely over 65 °F. On the other hand, no contribution of pressure can be observed." ] }, { "cell_type": "code", "execution_count": 8, "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", "
Generalized Linear Model Regression Results
Dep. Variable: Frequency No. Observations: 23
Model: GLM Df Residuals: 20
Model Family: Binomial Df Model: 2
Link Function: logit Scale: 1.0000
Method: IRLS Log-Likelihood: -3.7926
Date: Fri, 22 Oct 2021 Deviance: 2.7576
Time: 13:13:24 Pearson chi2: 4.19
No. Iterations: 6 Pseudo R-squ. (CS): 0.05416
Covariance Type: nonrobust
\n", "\n", "\n", " \n", "\n", "\n", " \n", "\n", "\n", " \n", "\n", "\n", " \n", "\n", "
coef std err z P>|z| [0.025 0.975]
Intercept 2.5202 8.541 0.295 0.768 -14.220 19.260
Temperature -0.0983 0.110 -0.894 0.371 -0.314 0.117
Pressure 0.0085 0.019 0.451 0.652 -0.028 0.045
" ], "text/plain": [ "\n", "\"\"\"\n", " 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.7926\n", "Date: Fri, 22 Oct 2021 Deviance: 2.7576\n", "Time: 13:13:24 Pearson chi2: 4.19\n", "No. Iterations: 6 Pseudo R-squ. (CS): 0.05416\n", "Covariance Type: nonrobust \n", "===============================================================================\n", " coef std err z P>|z| [0.025 0.975]\n", "-------------------------------------------------------------------------------\n", "Intercept 2.5202 8.541 0.295 0.768 -14.220 19.260\n", "Temperature -0.0983 0.110 -0.894 0.371 -0.314 0.117\n", "Pressure 0.0085 0.019 0.451 0.652 -0.028 0.045\n", "===============================================================================\n", "\"\"\"" ] }, "execution_count": 8, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# Estimation of temperature and pressure influence\n", "data[\"Intercept\"]=1\n", "\n", "logmodel=sm.GLM(data['Frequency'], data[['Intercept','Temperature','Pressure']], family=sm.families.Binomial(sm.families.links.logit())).fit()\n", "\n", "logmodel.summary()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Unlike the previous fit of logistic regression, the most likely estimator of the temperature parameter is not negligible and has negative value, implying a slightly negative correlation between temperature and probability of failure (the lower the temperature, the higher the probability). The p-value corresponding to the null hypothesis that its parameter can be set to 0 without loss of information is also lower (0.371). The estimated parameter for the pressure is instead very close to 0. A new model is therefore fit without taking pressure into account." ] }, { "cell_type": "code", "execution_count": 9, "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", "
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: Fri, 22 Oct 2021 Deviance: 3.0144
Time: 13:13:24 Pearson chi2: 5.00
No. Iterations: 6 Pseudo R-squ. (CS): 0.04355
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: Fri, 22 Oct 2021 Deviance: 3.0144\n", "Time: 13:13:24 Pearson chi2: 5.00\n", "No. Iterations: 6 Pseudo R-squ. (CS): 0.04355\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": 9, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# Estimation of temperature influence\n", "data[\"Intercept\"]=1\n", "\n", "logmodel=sm.GLM(data['Frequency'], data[['Intercept','Temperature']], family=sm.families.Binomial(sm.families.links.logit())).fit()\n", "\n", "logmodel.summary()" ] }, { "cell_type": "code", "execution_count": 10, "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXQAAAEKCAYAAAACS67iAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjQuMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/MnkTPAAAACXBIWXMAAAsTAAALEwEAmpwYAAAs7klEQVR4nO3deXhU5fn/8fc9k8lGFiBAWMKmhCCyJ2ERqsEqoFXcqICIS0Fsq1ZrpZUuam3t8qPfr3sVCrjUKqJVROtXKCguiBCQHWRHSNiXhITsyf37YyYYMZBJMsksuV/XlStzzjznnPvhkE9OzpzzHFFVjDHGBD+HvwswxhjjGxboxhgTIizQjTEmRFigG2NMiLBAN8aYEGGBbowxIaLGQBeROSJyWEQ2nuV9EZGnRGSHiKwXkQG+L9MYY0xNvDlCfxEYdY73rwCSPV9TgOfqX5YxxpjaqjHQVfUT4Pg5mlwDvKxuXwDNRaSdrwo0xhjjnTAfrKMDsK/KdJZn3oEzG4rIFNxH8URFRaV27NixThusqKjA4QiN0//Wl8ATKv0A60ugqk9ftm3bdlRVW1f3ni8C3WuqOhOYCZCWlqarVq2q03qWLl1KRkaGDyvzH+tL4AmVfoD1JVDVpy8i8vXZ3vPFr7tsoOqhdpJnnjHGmEbki0BfANziudplMJCrqt853WKMMaZh1XjKRUReAzKAViKSBTwMuABU9XngfeBKYAdQANzeUMUaY4w5uxoDXVXH1/C+Anf5rCJjTNAoLS0lKyuLoqKiBt9WfHw8W7ZsafDtNAZv+hIZGUlSUhIul8vr9Tbqh6LGmNCSlZVFbGwsXbp0QUQadFt5eXnExsY26DYaS019UVWOHTtGVlYWXbt29Xq9oXENkDHGL4qKikhISGjwMG9qRISEhIRa/+VjgW6MqRcL84ZRl39XC3RjjAkRdg7dGBPUnE4nvXv3Pj09f/58unTp4r+C/MgC3RgT1KKioli7dm2176kqqhoyQwbUpGn00hjTZOzZs4eUlBRuueUWevXqxb59+5g+fTrp6en06dOHhx9++HTbxx57jO7duzNs2DDGjx/P3/72NwAyMjKoHJrk6NGjp4/4y8vLmTp16ul1zZgxA/jmVv4xY8bQo0cPJkyYgPuKbsjMzOSiiy6ib9++DBw4kLy8PEaNGvWtX0LDhg1j3bp19e67HaEbY3zi9+9uYvP+kz5dZ8/2cTx89YXnbFNYWEi/fv0A6Nq1K48//jjbt2/npZdeYvDgwSxatIjt27ezcuVKVJXRo0fzySef0KxZM+bOncvatWspKytjwIABpKamnnNbs2fPJj4+nszMTIqLixk6dCgjRowAYM2aNWzatIn27dszdOhQli1bxsCBAxk7diyvv/466enpnDx5kqioKCZOnMiLL77IE088wbZt2ygqKqJv3771/veyQDfGBLUzT7ns2bOHzp07M3jwYAAWLVrEokWL6N+/PwD5+fls376dvLw8rrvuOqKjowEYPXp0jdtatGgR69ev58033wQgNzeX7du3Ex4ezsCBA0lKSgKgX79+7Nmzh/j4eNq1a0d6ejoAcXFxAFx33XUMHTqU6dOnM2fOHG677Taf/FtYoBtjfKKmI+nG1KxZs9OvVZVp06Zx5513fqvNE088cdblw8LCqKioAPjWteCqytNPP83IkSO/1X7p0qVEREScnnY6nZSVlZ11/dHR0Vx++eW88847zJs3j9WrV3vVr5rYOXRjTEgbOXIkc+bMIT8/H4Ds7GwOHz7MxRdfzPz58yksLCQvL49333339DJdunQ5HbKVR+OV63ruuecoLS0FYNu2bZw6deqs205JSeHAgQNkZmYC7jtEK4N+8uTJ/OxnPyM9PZ0WLVr4pK92hG6MCWkjRoxgy5YtDBkyBICYmBheeeUVBgwYwNixY+nbty9t2rQ5fVoE4IEHHuDGG29k5syZ/OAHPzg9f/LkyezZs4cBAwagqrRu3Zr58+efddvh4eG8/vrr3HPPPRQWFhIVFcXixYsBSE1NJS4ujttv9+F4hpWX9TT2V2pqqtbVRx99VOdlA431JfCESj9UG74vmzdvbtD1V3Xy5MkGXf/DDz+s06dPb9BtVDp58qRmZ2drcnKylpeXn7Vddf++wCo9S67aKRdjjGlkr776KoMGDeKxxx7z6TXydsrFGGOARx55pNG2ddNNN33nQ1pfsCN0Y0y9qOcGGuNbdfl3tUA3xtRZZGQkx44ds1D3MfWMhx4ZGVmr5eyUizGmzpKSksjKyuLIkSMNvq2ioqJaB1yg8qYvlU8sqg0LdGNMnblcrlo9Uac+li5devpuz2DXUH2xUy7GGBMiLNCNMSZEWKAbY0yIsEA3xpgQYYFujDEhwgLdGGNChAW6McaECAt0Y4wJERboxhgTIizQjTEmRARdoB8+WcTHWaU2GJAxxpwh6AL9lRV7eWFjCZNfWsWRvGJ/l2OMMQEj6AL9vu8nc1OPcD7dcZRRT3zC4s2H/F2SMcYEhKALdIdDGNHFxX/uGUZiXCSTX17Fw+9spKi03N+lGWOMXwVdoFdKTozl7bsuYtKwrry0/GuufXYZOw7n+7ssY4zxm6ANdICIMCe/u6onL9yezuG8YkY/8xnvrM32d1nGGOMXXgW6iIwSka0iskNEHqzm/U4i8pGIrBGR9SJype9LPbvhKW34z8+GcWH7OO6du5Zfv72B4jI7BWOMaVpqDHQRcQLPAlcAPYHxItLzjGa/Beapan9gHPB3Xxdak3bxUbx2x2B+fMn5vLpiLzc+v5zsnMLGLsMYY/zGmyP0gcAOVd2lqiXAXOCaM9ooEOd5HQ/s912J3gtzOnjwih7MmJjKriOnuOqpT/l8x1F/lGKMMY1OarpBR0TGAKNUdbJneiIwSFXvrtKmHbAIaAE0Ay5T1dXVrGsKMAUgMTExde7cuXUqOj8/n5iYmHO2OXiqgqfWFHHwlDI2JZwRncMQkTptryF505dgESp9CZV+gPUlUNWnL8OHD1+tqmnVvqmq5/wCxgCzqkxPBJ45o839wC88r4cAmwHHudabmpqqdfXRRx951S6vqFTveClTO//qPb3/9bVaVFpW5202FG/7EgxCpS+h0g9V60ugqk9fgFV6llz15pRLNtCxynSSZ15Vk4B5nl8Qy4FIoJUX625QMRFhPH9zKvddlsy/v8xiwj9WcDTf7i41xoQmbwI9E0gWka4iEo77Q88FZ7TZC3wfQEQuwB3oR3xZaF05HMJ9l3Xn2ZsGsHF/Ltc8s4ytB/P8XZYxxvhcjYGuqmXA3cBCYAvuq1k2icijIjLa0+wXwB0isg54DbjN86dBwPhBn3bMu3MIpeUVjHnucz7dHhC/b4wxxme8ug5dVd9X1e6qer6qPuaZ95CqLvC83qyqQ1W1r6r2U9VFDVl0XfVJas78u4bSoUUUt7+QyeuZe/1dkjHG+ExQ3ylaF+2bR/HGj4dwUbdW/OrfG3hi8TYbitcYExKaXKADxEa6mH1rGmNSk3hi8XamvbWBsvIKf5dljDH1EubvAvzF5XQwfUwf2sdH8tSHOziaX8wzNw0g0uX0d2nGGFMnTfIIvZKIcP+IFP5wzYUs+eowt8xeSW5hqb/LMsaYOmnSgV5p4pAuPDWuP2v2nWDsjOX2JCRjTFCyQPe4um97Zt+azp5jpxg7wwb2MsYEHwv0Ki7u3pp/ThrEkbxibnx+OXuOnvJ3ScYY4zUL9DOkd2nJa1MGU1hazo0zlttTkIwxQcMCvRq9OsQzd8pgKhTGzVzOVwdP+rskY4ypkQX6WXRPjGXenYMJczgYN/MLNu3P9XdJxhhzThbo53Be6xhev3Mw0S4nE2atYGO2hboxJnBZoNegc0Iz5k4ZYqFujAl4Fuhe6JQQzdwpQ2gW7uTm2SvYcsDOqRtjAo8Fupc6JUTz2pTBRIY5uXnWCrYfsjHVjTGBxQK9FjonNOPVOwbhdAjj/7GCXUfskkZjTOCwQK+l81rH8Oodg1BVJsxawb7jBf4uyRhjAAv0OunWJpZ/ThpEQUk5E2at4GBukb9LMsYYC/S66tk+jpd+NJDjp0q4efYKjp8q8XdJxpgmzgK9Hvp1bM6sW9PYd7yAW+esJK/Iht41xviPBXo9DT4vgeduHsCWAyeZ/NIqikrL/V2SMaaJskD3gUt7JPI/N/Zl5Z7j3P3qGnucnTHGLyzQfeSafh145OoLWbzlENPe2mAPnjbGNLom+0zRhnDrRV04dqqEp5ZsJyEmggev6OHvkowxTYgFuo/9/LJkjp8q5vmPd9ImNoIfDevq75KMMU2EBbqPiQi/H92Lo3kl/OE/m2kdG8HVfdv7uyxjTBNg59AbgNMhPDGuH+mdW3L/vLV8vuOov0syxjQBFugNJNLl5B+3pNG1VTPu/Odqe+qRMabBWaA3oPhoFy/cPpDoCCe3v5DJgdxCf5dkjAlhFugNrEPzKF64bSB5RWXc/kKm3U1qjGkwFuiNoGf7OJ67eQA7Dufz0399SandeGSMaQAW6I3ke8mt+dN1vfl0+1F++/ZGu/HIGONzdtliI7oxvSP7ThTw9Ic76JQQzV3Du/m7JGNMCLFAb2T3X96dvccLmL5wK50Toonxd0HGmJBhp1wamYjw1xv6kNa5BffPW8eOEzY6ozHGN7wKdBEZJSJbRWSHiDx4ljY3ishmEdkkIq/6tszQEulyMvOWNNrFR/LkmiJ7jJ0xxidqDHQRcQLPAlcAPYHxItLzjDbJwDRgqKpeCNzn+1JDS8tm4cy5LZ3yCpj0UiYn7XJGY0w9eXOEPhDYoaq7VLUEmAtcc0abO4BnVfUEgKoe9m2Zoen81jHc3T+SXUdO2Tjqxph6k5ounxORMcAoVZ3smZ4IDFLVu6u0mQ9sA4YCTuARVf2gmnVNAaYAJCYmps6dO7dORefn5xMTExofJ+bn57P6RAQvbCrhsk5h3Nwzwt8l1Vmo7JdQ6QdYXwJVffoyfPjw1aqaVt17vrrKJQxIBjKAJOATEemtqjlVG6nqTGAmQFpammZkZNRpY0uXLqWuywaapUuX8vBVGTjf28ysz3aTMaAHNw/u7O+y6iRU9kuo9AOsL4GqofrizSmXbKBjlekkz7yqsoAFqlqqqrtxH60n+6bEpmHalRcwPKU1Dy/YZKMzGmPqxJtAzwSSRaSriIQD44AFZ7SZj/voHBFpBXQHdvmuzNDndAhPje/P+a2b8eNXVrP76Cl/l2SMCTI1BrqqlgF3AwuBLcA8Vd0kIo+KyGhPs4XAMRHZDHwETFXVYw1VdKiKjXQx65Z0nA5h0kuZ5BbalS/GGO95dR26qr6vqt1V9XxVfcwz7yFVXeB5rap6v6r2VNXeqlq3TzsNnRKiee7mVPYeK+Ce1+zKF2OM9+xO0QA0+LwE/nBtLz7ZdoQ//99X/i7HGBMkbCyXADV+YCe2Hsxj9me76dE2lh+mdax5IWNMk2ZH6AHstz+4gKHdEvjN2xtZ/fUJf5djjAlwFugBLMzp4NmbBtCueSR3/nO1PcLOGHNOFugBrnl0OLNuSaOotJwpL6+mqNRGZzTGVM8CPQgkJ8byxNh+bNyfy6/+vd6edmSMqZYFepC4rGciD4xI4Z21+5nxid2zZYz5Lgv0IPLTjPP5QZ92/PWDr1i61Qa0NMZ8mwV6EBERpo/pQ4+2cdzz2hp2Hcn3d0nGmABigR5kosPDmDkxFZfTwR0vr7IHYxhjTrNAD0IdW0bz7E0D2HOsgPtfX0tFhX1IaoyxQA9aQ85P4KGrerJ4y2EeX7zN3+UYYwKA3fofxG4Z0plN+3N5+sMdXNAujit7t/N3ScYYP7Ij9CAmIvzh2l7079ScB95Yx1cHT/q7JGOMH1mgB7mIMCfP35xKTEQYd7y8ipyCEn+XZIzxEwv0EJAYF8nzE1M5lFvM3a/aGOrGNFUW6CFiQKcW/PHaXny24yh/sTHUjWmS7EPREHJjekc27s9l1me7ubBDHNf1T/J3ScaYRmRH6CHmd1f1ZFDXljz47w2sz8rxdznGmEZkgR5iXE4Hf58wgFYxEdz5z9UcySv2d0nGmEZigR6CEmIimDExlRMFJfz0X6spKbMPSY1pCizQQ1SvDvH89YY+ZO45wSPvbvJ3OcaYRmAfioawa/p1YPOBk8z4eBcXto9jwqDO/i7JGNOA7Ag9xP1yZA8u6d6ah9/ZxMrdx/1djjGmAVmghzinQ3hqfH86tozmJ6+sJjvHHjRtTKiyQG8C4qNc/OOWNErKKpjy8ioKS+xB08aEIgv0JqJbmxieHN+PzQdOMvXNdfagaWNCkAV6E3Jpj0SmjkzhvfUH+PvSnf4uxxjjYxboTcxPLjmfq/u252+LtrJkyyF/l2OM8SEL9CZGRPh/N/ThwvZx3Dt3LdsP5fm7JGOMj1igN0FR4U5mTkwj0uVkso2hbkzIsEBvoto3j2LGxFQO5BTx0399SamNoW5M0LNAb8JSO7fgT9f35vOdx/jDe5v9XY4xpp7s1v8mbkxqEtsO5THzk10kJ8YycbAND2BMsLIjdMOvRvVgeEprHlmwic93HPV3OcaYOvIq0EVklIhsFZEdIvLgOdrdICIqImm+K9E0tMrhAc5r1Yyf/OtLdh895e+SjDF1UGOgi4gTeBa4AugJjBeRntW0iwXuBVb4ukjT8GIjXcy+NR2nQ5j0Yia5BaX+LskYU0veHKEPBHao6i5VLQHmAtdU0+4PwF+BIh/WZxpRp4Ronr85lX0nCrjrVbvyxZhgIzWN6SEiY4BRqjrZMz0RGKSqd1dpMwD4jareICJLgQdUdVU165oCTAFITExMnTt3bp2Kzs/PJyYmpk7LBppA7MunWaXM3lhCRscwbu0Zjoh4tVwg9qUuQqUfYH0JVPXpy/Dhw1erarWntet9lYuIOID/BW6rqa2qzgRmAqSlpWlGRkadtrl06VLqumygCcS+ZADhH3zFc0t38r2+3Zk0rKtXywViX+oiVPoB1pdA1VB98eaUSzbQscp0kmdepVigF7BURPYAg4EF9sFocJs6IoVRF7blj//ZzOLNNuaLMcHAm0DPBJJFpKuIhAPjgAWVb6pqrqq2UtUuqtoF+AIYXd0pFxM8HA7h8bH96NU+np/NXcPG7Fx/l2SMqUGNga6qZcDdwEJgCzBPVTeJyKMiMrqhCzT+ExXuZPataTSPcjHppUwO5NrTjowJZF5dh66q76tqd1U9X1Uf88x7SFUXVNM2w47OQ0ebuEjm3J7OqeJybn8hk7wiu5zRmEBld4qaGvVoG8ezEwaw/XA+d726xi5nNCZAWaAbr1zSvTWPXduLT7Yd4aF3Ntoj7IwJQDY4l/HauIGd2HeigGc/2klSi2juGt7N3yUZY6qwQDe18ovLU8g6Ucj0hVtpFx/J9QOS/F2SMcbDAt3UisMhTB/TlyN5xfzyzfW0iY1kWHIrf5dljMHOoZs6CA9z8PzEVLq1ieHHr6xm0/66X6M+f002Q//yIV0f/A9D//Ih89dk17yQCRi2/wKLBbqpk7hIFy/ePpC4yDBueyGTfccLar2O+WuymfbWBrJzClEgO6eQaW9tsFAIErb/Ao8FuqmztvGRvDxpICVlFdwyZyUnS2p35cv0hVspLC3/1rzC0nKmL9zqyzJNA7H9F3gs0E29dGsTy5zb0jiQW8jjq4rILy7zetn9OdXfeXq2+Saw2P4LPBbopt5SO7fk7xMG8HVeBVNeXkVxWXnNCwHtm0fVar4JLLb/Ao8FuvGJS3skMqlXOJ/vPMZ9c9dSXlHz6ZepI1OIcjm/NS/K5WTqyJSGKtP4kO2/wGOBbnxmaAcXv7uqJ/+38SDT3lpf492k1/bvwJ+v702H5lEI0KF5FH++vjfX9u/QOAWberH9F3jsOnTjU5OGdSW3sJSnlmwnLtLFb35wwTmfeHRt/w4WAEHM9l9gsUA3Pvfzy5I5WVjKrM92Exvp4t7Lkv1dkjFNggW68TkR4aGrepJfXMbji7cRHe7kjovP83dZxoQ8C3TTIBwO4a839KGwtJzH3t9CpMvBxCFd/F2WMSHNAt00GKdDePzGfhSVlPO7dzYRHuZgbHonf5dlTMiyq1xMgwoPc/DshAFc3L01D761gX+vzvJ3ScaELAt00+AiXU5mTkzlovMTmPrmOt5Za2N9GNMQLNBNo4h0OZl1SzoDu7bk56+vtQGcjGkAFuim0USFO5lzWzqDuiZw/7y1vL3GTr8Y40sW6KZRRYeHMee2dAafl8D989Yxb9U+f5dkTMiwQDeNLircyexb0xnWrRW/fHM9r3zxtb9LMiYkWKAbv4gKd/KPW9L4fo82/Hb+RmZ/ttvfJRkT9CzQjd9Eupw8d3MqV/Rqyx/e28yTi7fXOKCXMebsLNCNX4WHOXh6fH9uGJDE44u38dh/tlioG1NHdqeo8bswp4PpY/oQGxnGrM92k1tYyp+v702Y0443jKkNC3QTEBwO4eGrexIf5eLJJds5UVDKMzf1J/KMBygYY87ODoFMwBARfn55d34/+kKWfHWIibNXkFNQ4u+yjAkaFugm4Nx6UReeGtefdftyGfP8crJOFPi7JGOCggW6CUhX923PSz8ayKGTRVz/98/ZmJ3r75KMCXgW6CZgDTk/gTd/fBFOh3DjjOV8+NUhf5dkTECzQDcBLaVtLPPvGkrXVs2Y/NIqXvp8j79LMiZgWaCbgJcYF8m8O4dwaY82PLxgEw+9s5Gy8gp/l2VMwPEq0EVklIhsFZEdIvJgNe/fLyKbRWS9iCwRkc6+L9U0Zc0iwpgxMY0pF5/Hy8u/5rYXMsktKPV3WcYElBoDXUScwLPAFUBPYLyI9Dyj2RogTVX7AG8C/8/XhRrjdAi/vvICpo/pw4rdxxj97GdsO5Tn77KMCRjeHKEPBHao6i5VLQHmAtdUbaCqH6lq5bVlXwBJvi3TmG/8MK0jc6cMpqCknOueXcYHGw/4uyRjAoLUNG6GiIwBRqnqZM/0RGCQqt59lvbPAAdV9Y/VvDcFmAKQmJiYOnfu3DoVnZ+fT0xMTJ2WDTTWl7o7UVTB02uK2ZVbwZVdXdyQ7MLpkHqv1/ZJYLK+uA0fPny1qqZV+6aqnvMLGAPMqjI9EXjmLG1vxn2EHlHTelNTU7WuPvroozovG2isL/VTVFqm095ar51/9Z6Om7Fcj+QV1Xudtk8Ck/XFDVilZ8lVb065ZAMdq0wneeZ9i4hcBvwGGK2qxd7+tjGmPiLCnPzput787Yd9+XLvCa588lOW7zzm77KM8QtvAj0TSBaRriISDowDFlRtICL9gRm4w/yw78s05tzGpCYx/66hxESEMWHWFzy1ZDvlFTYMr2laagx0VS0D7gYWAluAeaq6SUQeFZHRnmbTgRjgDRFZKyILzrI6YxrMBe3iWHDPMEb3bc///ncbE2Z9wYHcQn+XZUyj8Wr4XFV9H3j/jHkPVXl9mY/rMqZOFm8+xMrdxwFYses43/+fjxmX3pGFmw6xP6eQ9s2jmDoyhWv7d/D5tuevyWb6wq0Nvh1v/Hb+Bl5bsY/7epUyadr7jB/UkT9e29svtZjGY+Ohm5Axf002097aQGFpOQAKFJaUM2fZntNtsnMKmfbWBgCfhu2Z226o7Xjjt/M38MoXe09Pl6uenrZQD212678JGdMXbj0dqJWqO4teWFrO9IVbG3zbDbEdb7y2Yl+t5pvQYYFuQsb+HO/Pl2fXom19tl2bmnyl/Cz3lpxtvgkdFugmZLRvHuV1W6dD+Gz70Qbfdm1q8hWnVH9z1dnmm9BhgW5CxtSRKUSd8QxSl0NwOb8dZOFOBy2bhXPz7BU88MY6jp+q/2Puqtt2lMvJ1JEp9V53bY0f1LFW803osA9FTcio/PDxzCtNqps3qldbnlqynZmf7GLJlkP8+soLGJOahNTxKPZs2/bHVS6VH3xWnjN3ithVLk2EBboJKdf271BtiFY375ejejC6X3t++/ZGpr65nnmr9vH70b18vm1/+OO1vfnjtb1ZunQpOydk+Lsc00jslItp0nq0jWPenUP46w292XnkFFc9/Sn/3FxMTkH9T8MY09gs0E2T53AIY9M78dEvMrh5cGc+3FvGJdOX8sKy3ZTak5FMELFAN8YjPtrFo9f04tGhUfTqEMfv393MyMc/4YONBytHEzUmoFmgG3OGjrEOXpk0iFm3pOFwCD9+ZTVjnl/Oil02iqMJbBboxlRDRLisZyIf3Ps9/nx9b/YdL2DszC+4Zc5K1mfl+Ls8Y6plgW7MOYQ5HYwf2ImPpw7n11f2YH1WDqOfWcakFzMt2E3AsUA3xgtR4U6mXHw+n/5yOA+M6M6qr08w+pll3PbCSjL3HPd3ecYAFujG1EpspIu7L03ms18NZ+rIFDZk5fLD55dz4/PLWbz5EBX2UA3jRxboxtRBbKSLu4Z347NfXcpDV/UkO6eQyS+v4vLHP+bVFXspLCmveSXG+JgFujH1EBXu5EfDurJ0agZPjutHpMvJr9/ewJC/LOEv//cV+44X+LtE04TYrf/G+IDL6eCafh0Y3bc9mXtO8MKy3cz8ZCczPtnJpSltmDC4E5d0b4PTYSMemoZjgW6MD4kIA7u2ZGDXluzPKeS1lXt5beU+lry4ivbxkfwwrSNjUpPo2DLa36WaEGSBbkwDad88il+MSOGeS5NZsuUQr67cy1MfbufJJdsZcl4CN6QmMapXW2Ii7MfQ+Ib9TzKmgYWHObiidzuu6N2OrBMFvPVlNm+uzuKBN9bx2/kbuLxnW0b3bc/F3VsREeaseYXGnIUFujGNKKlFND/7fjL3XNqNL/ee4O012by3/gDvrttPbGQYI3q25YpebRmW3IpIl4W7qR0LdGP8QERI7dyS1M4tefjqC1m24yjvrjvAfzcf5N9fZhETEcYlKa0Z0TORjJQ2xEe5/F2yCQIW6Mb4mcvpICOlDRkpbSgp683nO4/ywcaDLN5ymP+sP4DTIaR3acGlPdxtktvE1PnJSia0WaAbE0DCw74J94oKZc2+HD786hBLthzmT+9/xZ/e/4p28ZF8L7kVw5Jbc9H5CbSKifB32SZAWKAbE6AcDiG1cwtSO7dg6sge7M8p5JNtR/h42xE+2HiQeauyAOjRNpbB5yUw+LwE0ru0IMECvsmyQDcmSLRvHsW4gZ0YN7AT5RXKhuxclu04yvKdx5ibuZcXP98DQLc2MaR5fhGUn6pAVe0UTRNhgW5MEHI6hH4dm9OvY3PuGt6N4rJyNmTlsnLPcTJ3H+f9DQeYm7kPgL+s/i99k5rTt2Nz+ibF0zspnjaxkX7ugWkIFujGhICIMCdpXVqS1qUlZEBFhbLzSD6vLvqCwuhE1u7L4ZkPt1M5GGSb2Ah6dYinZ7s4eraP44J2cXRqGW1DEwQ5C3RjQpDDISQnxnJJRxcZGX0AKCgpY9P+k6zPymXT/lw2Zufy8bYjlHtSPtLlILlNLN0TY+meGENyYgzdWsfSoUWUBX2QsEA3pomIDg8jvUtL0ru0PD2vqLScHYfz2XzgJFsP5rH1YB6fbj/Cv7/MOt0mIsxB11bN6NqqGV1aNaNrQjM6J0TTpVUzWsdE4LCwDxgW6MY0YZEuJ706xNOrQ/y35ucUlLDjcD47j+Sz43A+u4+eYuuhPP67+RBlVR7iERHmoGPLaJJaRNGxRTQdWkTRoXnU6e+tYiLs6L4RWaAbY76jeXT4N+fkqygrr2B/ThG7j51i7/EC9h0v4Otjp8g6UciavTnkFpZ+q32YQ0iMi6RtfCRt4yJJjIskMS6CNnERtImNpE1sBK1iImge7bIrcXzAAt0Y47Uwp4NOCdF0Sqh++N+8olKycwrJPlHI/twiDuQUcjC3iIMni9hy4CQfbT1MQTVPc3I5hYRmESTEhJMQE0FCs3BaVvlqER3O18fLaXcwjxbRLuKiXDbWTTUs0I0xPhMb6aJHWxc92sadtU1+cRmHThZx+GQxh/OKOJpfwtH8Yo7kFXP8lPv1riP5HD9V8p3w//PKT06/jnQ5iI9yERfpcn+PchEXGUZspItYz/eYyDBiI8JoFhFGjOerWYSTmIgwoiPCiHY5Q+ozAK8CXURGAU8CTmCWqv7ljPcjgJeBVOAYMFZV9/i2VGNC1/w12UxfuJX9OYW0bx7F1JEpvLFqL8t2Hj/dZuj5LflhWqfvtAO+M2/V18d5bcU+7utVyqRp7zN+UEf+eG1vr7Zb3fqu7d/B67ort12uilPkO9uOiQgjpnUMG7Jya9z2I1cn873urTh+qoSPl6+iU/IFnCgo5Yudx/h42xEOnSz2nOaJprC0nB2HyzhZVEpeUdnpq3dqEuVyEh3uJCq88nsYUS4H0eFhRLmcRLgcRLmcRLqcRLocRIZ98zoizP1+RJjndZiDCJeDcKeT8DDHN1/Ob767nIJqwzxMvMZAFxEn8CxwOZAFZIrIAlXdXKXZJOCEqnYTkXHAX4GxDVGwMaFm/ppspr21gcJS99Fodk4h972+9jvtlu08/q2Az84pZOqb60Ch1BNe2TmF3P/6WiqqLFeuyitf7AX4VrBWt92pb6wDgdLyb9Y37a0NAN8J9eqW9/W2H16wiT9f35tr+3fgSIKTjD7tmb8mmw+/Onx62aLSCrJOFJ5uB6CqFJVWkFdcSn5RGfnFnq+iMgpKyjlVUsap4jJOFZdzqriMgtJyCkvKKSgpo7C0gsKSMo7kFVPomV9UWk5hqfu7l78nzmliz3CG13813+HNEfpAYIeq7gIQkbnANUDVQL8GeMTz+k3gGRERbahfQ8aEkOkLt54Op9qqDL+qKqppB/Dain3fCtXqtltaTVoVlpYzfeHW7wR6dcs3xrarW/bMdiJClOeou03sWYqqA1WltFwpKiunuLSCotJySsorKC6toLisnJKyCorLKigpq3DPLyuntEwpLnfPKy2voLSsgtj8vb4rqgqpKXNFZAwwSlUne6YnAoNU9e4qbTZ62mR5pnd62hw9Y11TgCmeyRRgax3rbgUcrbFVcLC+BJ5G7Ud4226pDbXu8oJcnNHfXJJYcnDH6rput+qy9V2+jsu2Ao6ea9kzawxg9fk/1llVW1f3RqN+KKqqM4GZ9V2PiKxS1TQflOR31pfAEyr9AHdfynIPh0xfQmm/NERfHF60yQY6VplO8syrto2IhAHxuD8cNcYY00i8CfRMIFlEuopIODAOWHBGmwXArZ7XY4AP7fy5McY0rhpPuahqmYjcDSzEfdniHFXdJCKPAqtUdQEwG/iniOwAjuMO/YZU79M2AcT6EnhCpR9gfQlUDdKXGj8UNcYYExy8OeVijDEmCFigG2NMiAj4QBeRSBFZKSLrRGSTiPzeM7+riKwQkR0i8rrnA9uAJyJOEVkjIu95poO1H3tEZIOIrBWRVZ55LUXkvyKy3fO9hb/r9IaINBeRN0XkKxHZIiJDgrEvIpLi2R+VXydF5L4g7cvPPT/vG0XkNU8OBOvPyr2efmwSkfs88xpknwR8oAPFwKWq2hfoB4wSkcG4hxd4XFW7ASdwDz8QDO4FtlSZDtZ+AAxX1X5Vrqd9EFiiqsnAEs90MHgS+EBVewB9ce+foOuLqm717I9+uMdVKgDeJsj6IiIdgJ8BaaraC/fFGJVDigTVz4qI9ALuwH3HfV/gKhHpRkPtE1UNmi8gGvgSGIT7Lqswz/whwEJ/1+dF/UmenXcp8B4gwdgPT617gFZnzNsKtPO8bgds9XedXvQjHtiN5wKBYO7LGfWPAJYFY1+ADsA+oCXuK/HeA0YG488K8ENgdpXp3wG/bKh9EgxH6JWnKdYCh4H/AjuBHFUt8zTJwv2fINA9gXtnVg55kUBw9gNAgUUistozpANAoqoe8Lw+CCT6p7Ra6QocAV7wnAqbJSLNCM6+VDUOeM3zOqj6oqrZwN+AvcABIBdYTXD+rGwEviciCSISDVyJ+ybMBtknQRHoqlqu7j8jk3D/6dLDvxXVnohcBRxW1WAZa6Imw1R1AHAFcJeIXFz1TXUfegTDNbFhwADgOVXtD5zijD9/g6gvAHjOLY8G3jjzvWDoi+d88jW4f9m2B5oBo/xaVB2p6hbcp4oWAR8Aa4HyM9r4bJ8ERaBXUtUc4CPcf2419wwzANUPRxBohgKjRWQPMBf3aZcnCb5+AKePolDVw7jP0w4EDolIOwDP98P+q9BrWUCWqq7wTL+JO+CDsS+VrgC+VNVDnulg68tlwG5VPaKqpcBbuH9+gvVnZbaqpqrqxbjP/W+jgfZJwAe6iLQWkeae11G4x2XfgjvYx3ia3Qq845cCvaSq01Q1SVW74P5z+ENVnUCQ9QNARJqJSGzla9znazfy7SEggqIvqnoQ2CciKZ5Z38c9NHTQ9aWK8XxzugWCry97gcEiEi0iwjf7JOh+VgBEpI3neyfgeuBVGmifBPydoiLSB3gJ9yfdDmCeqj4qIufhPtJtCawBblbVYv9V6j0RyQAeUNWrgrEfnprf9kyGAa+q6mMikgDMAzoBXwM3qurxs6wmYIhIP2AWEA7sAm7H83+N4OtLM9yBeJ6q5nrmBd1+8VyePBYow/1zMRn3OfOg+lkBEJFPcX9eVgrcr6pLGmqfBHygG2OM8U7An3IxxhjjHQt0Y4wJERboxhgTIizQjTEmRFigG2NMiGjUh0Qb4y3PZV1LPJNtcd9dd8QzPVBVS/xSWDU8l6GWqOrnfi7FNHEW6CYgqeox3KNrIiKPAPmq+jd/1SMiYVXGETlTBpAPeB3oNazPmDqxUy4maIhIqoh87BkQbGGVW6eXisjjIrLKM555uoi85Rlr+o+eNl08453/y9PmTc9gSTWt9wlxj/d+r4hc7RmPe42ILBaRRBHpAvwY+LlnDPLviciLIjKmSt35nu8ZIvKpiCwANnsGnZsuIpkisl5E7mzUf1ATcizQTbAQ4GlgjKqmAnOAx6q8X6Lucdmfx30b9V1AL+A2z+kbgBTg76p6AXAS+KmIuGpYb7iqpqnq/wCfAYM9g3jNBX6pqns823xc3WORf1pDPwYA96pqd9zjeeeqajqQDtwhIl1r/09jjJudcjHBIgJ3QP/XPbwHTtxDq1Za4Pm+AdhUOTSpiOzCPVxpDrBPVZd52r2C+yEKH9Sw3tervE4CXvccwYfjHke9tlaqauVyI4A+VY7m44HkOq7XGAt0EzQEd1APOcv7lWN6VFR5XTld+f/8zHEu1Iv1nqry+mngf1V1geeD0EfOskwZnr9+RcSBO/yrW58A96jqwrOsx5hasVMuJlgUA61FZAiAiLhE5MJarqNT5fLATbhPoWytxXrj+WbI1lurzM8DYqtM78H9CDhwj0vuOsv6FgI/8Zz2QUS6ewbXMqZOLNBNsKjAPXTqX0VkHe4HBVxUy3Vsxf0wji1AC9wPtSipxXofAd4QkdW4H4dW6V3gusoPRYF/AJd41jeEbx+VVzUL97CwX4rIRmAG9lezqQcbbdE0CZ6rUd5T90OHjQlJdoRujDEhwo7QjTEmRNgRujHGhAgLdGOMCREW6MYYEyIs0I0xJkRYoBtjTIj4/yb9e+wuqsYUAAAAAElFTkSuQmCC\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "# Estimation of the probability of O-ring malfunction\n", "data_pred = pd.DataFrame({'Temperature': np.linspace(start=30, stop=90, num=121), 'Intercept': 1})\n", "data_pred['Frequency'] = logmodel.predict(data_pred[['Intercept','Temperature']])\n", "data_pred.plot(x=\"Temperature\",y=\"Frequency\",kind=\"line\",ylim=[0,1])\n", "plt.scatter(x=data[\"Temperature\"],y=data[\"Frequency\"])\n", "plt.grid(True)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The final prediction of O-ring failure depending on the temperature is completely different than before! Considering all the flights for the analysis, successful and failed, the estimated probability of failure around 30 °F is slightly higher than 80%. This result would have suggested that the Challenger flight was very unadvisable in the morning of January 28, 1986." ] } ], "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.4" }, "toc": { "base_numbering": 1, "nav_menu": {}, "number_sections": false, "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 }