{ "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
131/24/856532002
82/03/846572001
221/12/866582001
94/06/846632001
04/12/81666500
144/12/856672000
1211/08/846672000
44/04/83667500
311/11/82668500
23/22/82669500
111/12/81670501
166/17/856702000
711/28/836701000
108/30/846702001
56/18/82672500
68/30/836731000
154/29/856752000
2010/30/856752002
188/27/856762000
2111/26/856762000
1110/05/846782000
1910/03/856792000
177/29/856812000
\n", "
" ], "text/plain": [ " Date Count Temperature Pressure Malfunction\n", "13 1/24/85 6 53 200 2\n", "8 2/03/84 6 57 200 1\n", "22 1/12/86 6 58 200 1\n", "9 4/06/84 6 63 200 1\n", "0 4/12/81 6 66 50 0\n", "14 4/12/85 6 67 200 0\n", "12 11/08/84 6 67 200 0\n", "4 4/04/83 6 67 50 0\n", "3 11/11/82 6 68 50 0\n", "2 3/22/82 6 69 50 0\n", "1 11/12/81 6 70 50 1\n", "16 6/17/85 6 70 200 0\n", "7 11/28/83 6 70 100 0\n", "10 8/30/84 6 70 200 1\n", "5 6/18/82 6 72 50 0\n", "6 8/30/83 6 73 100 0\n", "15 4/29/85 6 75 200 0\n", "20 10/30/85 6 75 200 2\n", "18 8/27/85 6 76 200 0\n", "21 11/26/85 6 76 200 0\n", "11 10/05/84 6 78 200 0\n", "19 10/03/85 6 79 200 0\n", "17 7/29/85 6 81 200 0" ] }, "execution_count": 1, "metadata": {}, "output_type": "execute_result" } ], "source": [ "import numpy as np\n", "import pandas as pd\n", "data = pd.read_csv(\"shuttle.csv\")\n", "data.sort_values('Temperature')" ] }, { "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": "iVBORw0KGgoAAAANSUhEUgAAAYUAAAEKCAYAAAD9xUlFAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMi4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvIxREBQAAFYRJREFUeJzt3XuQpXV95/H3Zy7AIBMhsJm4MxBBCFlKAXG4GEx2IokLbgmxiBHcDS5ZMqGE3TK7m8BariHGVEWM2WiJjiOLCqmERFEgu+MiJNUaExCQTIaLgcwiQjMGBFFoHObW3/3jnHlyprun5/TQzzlM9/tV1TXnufa3vz6cj8/l/E6qCkmSABYMuwBJ0kuHoSBJahgKkqSGoSBJahgKkqSGoSBJarQWCkmuSfJkkvt2szxJPppkY5INSU5qqxZJUn/aPFP4DHDmNMvPAo7p/qwGPtFiLZKkPrQWClX1VeB706xyDnBtddwBHJzkFW3VI0nas0VD/N3Lgcd6pke7874zccUkq+mcTbBkyZLXHX744QMp8MUaHx9nwQJv2/SyJ5PZk6nZl8leTE8eeuihp6rqX+xpvWGGQqaYN+WYG1W1FlgLsHLlyrr77rvbrGvWjIyMsGrVqmGX8ZJiTyazJ1OzL5O9mJ4k+XY/6w0zhkeB3v/LvwLYNKRaJEkMNxRuBi7oPoV0GvCDqpp06UiSNDitXT5K8qfAKuCwJKPAbwOLAapqDbAOeDOwEfghcGFbtUiS+tNaKFTV+XtYXsAlbf1+SdLMeWtfktQwFCRJDUNBktQwFCRJDUNBktQwFCRJDUNBktQwFCRJDUNBktQwFCRJDUNBktQwFCRJDUNBktQwFCRJDUNBktQwFCRJDUNBktQwFCRJDUNBktQwFCRJDUNBktQwFCRJDUNBktQwFCRJDUNBktQwFCRJDUNBktQwFCRJDUNBktQwFCRJDUNBktQwFCRJDUNBktQwFCRJDUNBktRoNRSSnJnkwSQbk1w+xfKXJ/mLJH+f5P4kF7ZZjyRpeq2FQpKFwFXAWcBxwPlJjpuw2iXAA1V1ArAK+HCS/dqqSZI0vTbPFE4BNlbVw1W1FbgeOGfCOgUsTRLgIOB7wPYWa5IkTWNRi/teDjzWMz0KnDphnY8BNwObgKXA26tqfOKOkqwGVgMsW7aMkZGRNuqddWNjY/tMrYNiTyazJ1OzL5MNoidthkKmmFcTpv8NsB54I/Aq4NYkf11Vz+6yUdVaYC3AypUra9WqVbNfbQtGRkbYV2odFHsymT2Zmn2ZbBA9afPy0ShweM/0CjpnBL0uBL5QHRuBbwE/1WJNkqRptBkKdwHHJDmye/P4PDqXino9CpwBkGQZcCzwcIs1SZKm0drlo6ranuRS4BZgIXBNVd2f5OLu8jXA7wKfSXIvnctNl1XVU23VJEmaXpv3FKiqdcC6CfPW9LzeBLypzRokSf3zE82SpIahIElqGAqSpIahIElqGAqSpIahIElqGAqSpIahIElqGAqSpIahIElqGAqSpIahIElqGAqSpIahIElqGAqSpIahIElqGAqSpIahIElqGAqSpIahIElqGAqSpIahIElqGAqSpIahIElqGAqSpIahIElqGAqSpIahIElqGAqSpIahIElqGAqSpIahIElqGAqSpIahIElqtBoKSc5M8mCSjUku3806q5KsT3J/kq+0WY8kaXqL+lkpyaur6r6Z7DjJQuAq4BeAUeCuJDdX1QM96xwMfBw4s6oeTfJjM/kdkqTZ1e+ZwpokdyZ5V/eNvB+nABur6uGq2gpcD5wzYZ13AF+oqkcBqurJPvctSWpBX2cKVfWGJMcAvwrcneRO4NNVdes0my0HHuuZHgVOnbDOTwKLk4wAS4GPVNW1E3eUZDWwGmDZsmWMjIz0U/bQjY2N7TO1Doo9mcyeTM2+TDaInvQVCgBV9Y9J3gvcDXwUeG2SAO+pqi9MsUmm2s0Uv/91wBnAEuD2JHdU1UMTfvdaYC3AypUra9WqVf2WPVQjIyPsK7UOij2ZzJ5Mzb5MNoie9HtP4XjgQuDfArcCb6mqe5L8S+B2YKpQGAUO75leAWyaYp2nqup54PkkXwVOAB5CkjRw/d5T+BhwD3BCVV1SVfcAVNUm4L272eYu4JgkRybZDzgPuHnCOjcBP5NkUZID6Vxe+uZM/whJ0uzo9/LRm4HNVbUDIMkC4ICq+mFVXTfVBlW1PcmlwC3AQuCaqro/ycXd5Wuq6ptJ/i+wARgHrp7pU06SpNnTbyjcBvw8MNadPhD4MvDT021UVeuAdRPmrZkw/SHgQ33WIUlqUb+Xjw6oqp2BQPf1ge2UJEkaln5D4fkkJ+2cSPI6YHM7JUmShqXfy0fvBj6XZOfTQ68A3t5OSZKkYen3w2t3Jfkp4Fg6nz/4h6ra1mplkqSB6/vDa8DJwCu727w2CVN9+liStO/q98Nr1wGvAtYDO7qzCzAUJGkO6fdMYSVwXFVNHKZCkjSH9Pv00X3Aj7dZiCRp+Po9UzgMeKA7OuqWnTOr6uxWqpIkDUW/oXBFm0VIkl4a+n0k9StJfgI4pqpu6w5et7Dd0iRJg9bXPYUkvwZ8Hvhkd9Zy4Ma2ipIkDUe/N5ovAU4HnoXOF+4Afp+yJM0x/YbClu73LAOQZBGTv0VNkrSP6zcUvpLkPcCSJL8AfA74i/bKkiQNQ7+hcDnwXeBe4NfpfEfC7r5xTZK0j+r36aNx4FPdH0nSHNXv2EffYop7CFV11KxXJEkampmMfbTTAcDbgB+d/XIkScPU1z2Fqnq65+fxqvoj4I0t1yZJGrB+Lx+d1DO5gM6Zw9JWKpIkDU2/l48+3PN6O/AI8MuzXo0kaaj6ffro59ouRJI0fP1ePvov0y2vqj+cnXIkScM0k6ePTgZu7k6/Bfgq8FgbRUmShmMmX7JzUlU9B5DkCuBzVXVRW4VJkgav32EujgC29kxvBV4569VIkoaq3zOF64A7k3yRzieb3wpc21pVkqSh6Pfpo99L8iXgZ7qzLqyqv2uvLEnSMPR7+QjgQODZqvoIMJrkyJZqkiQNSb9fx/nbwGXAf+/OWgz8cVtFSZKGo98zhbcCZwPPA1TVJhzmQpLmnH5DYWtVFd3hs5O8rL2SJEnD0m8o/HmSTwIHJ/k14Db8wh1JmnP6ffroD7rfzfwscCzwvqq6tdXKJEkDt8czhSQLk9xWVbdW1W9W1X/rNxCSnJnkwSQbk1w+zXonJ9mR5JdmUrwkaXbtMRSqagfwwyQvn8mOkywErgLOAo4Dzk9y3G7W+yBwy0z2L0maff1+ovkF4N4kt9J9Agmgqv7zNNucAmysqocBklwPnAM8MGG9/wTcQGfAPUnSEPUbCv+n+zMTy9l1FNVR4NTeFZIsp/O46xuZJhSSrAZWAyxbtoyRkZEZljIcY2Nj+0ytg2JPJrMnU7Mvkw2iJ9OGQpIjqurRqvrsXuw7U8yrCdN/BFxWVTuSqVbvblS1FlgLsHLlylq1atVelDN4IyMj7Cu1Doo9mcyeTM2+TDaInuzpnsKNO18kuWGG+x4FDu+ZXgFsmrDOSuD6JI8AvwR8PMkvzvD3SJJmyZ4uH/X+3/ejZrjvu4BjumMkPQ6cB7yjd4WqasZPSvIZ4H9X1Y1IkoZiT6FQu3m9R1W1PcmldJ4qWghcU1X3J7m4u3zNjCqVJLVuT6FwQpJn6ZwxLOm+pjtdVfUj021cVeuAdRPmTRkGVfUf+qpYktSaaUOhqhYOqhBJ0vDN5PsUJElznKEgSWoYCpKkhqEgSWrMq1B4emwLf//Y93l6bMuwS5GkGXl6bAubt+1o/f1r3oTCTesf5/QP/hX//uqvc/oH/4qb1z8+7JIkqS8737++9d3nW3//mheh8PTYFi67YQMvbBvnuS3beWHbOL91wwbPGCS95PW+f+2oav39a16Ewugzm1m8YNc/dfGCBYw+s3lIFUlSfwb9/jUvQmHFIUvYNj6+y7xt4+OsOGTJkCqSpP4M+v1rXoTCoQftz5XnHs8BixewdP9FHLB4AVeeezyHHrT/sEuTpGn1vn8tTFp//+r3S3b2eWefuJzTjz6M0Wc2s+KQJQaCpH3GzvevO2//Gn9z9htaff+aN6EAncQ1DCTtiw49aH+WLF7Y+nvYvLh8JEnqj6EgSWoYCpKkhqEgSWoYCpKkhqEgSWoYCpKkhqEgSWoYCpKkhqEgSWoYCpKkhqEgSWoYCpKkhqEgSWoYCpKkhqEgSWoYCpKkhqEgSWoYCpKkhqEgSWoYCpKkRquhkOTMJA8m2Zjk8imW/7skG7o/f5vkhDbrkSRNr7VQSLIQuAo4CzgOOD/JcRNW+xbwr6vqeOB3gbVt1SNJ2rM2zxROATZW1cNVtRW4Hjind4Wq+tuqeqY7eQewosV6JEl7sKjFfS8HHuuZHgVOnWb9/wh8aaoFSVYDqwGWLVvGyMjILJXYrrGxsX2m1kGxJ5PZk6nZl8kG0ZM2QyFTzKspV0x+jk4ovGGq5VW1lu6lpZUrV9aqVatmqcR2jYyMsK/UOij2ZDJ7MjX7MtkgetJmKIwCh/dMrwA2TVwpyfHA1cBZVfV0i/VIkvagzXsKdwHHJDkyyX7AecDNvSskOQL4AvArVfVQi7VIkvrQ2plCVW1PcilwC7AQuKaq7k9ycXf5GuB9wKHAx5MAbK+qlW3VJEmaXpuXj6iqdcC6CfPW9Ly+CLiozRrmi6fHtjD6zGZWHLKEQw/av/Xt5jJ7Mnwbn3iOZ364jY1PPMfRy5YOu5x5pdVQ0GDctP5xLrthA4sXLGDb+DhXnns8Z5+4vLXt5jJ7Mnzvu/Ferr3jUf7ra7bzG//zq1zw+iN4/zmvGXZZ84bDXOzjnh7bwmU3bOCFbeM8t2U7L2wb57du2MDTY1ta2W4usyfDt/GJ57j2jkd3mXft7Y+y8YnnhlTR/GMo7ONGn9nM4gW7/s+4eMECRp/Z3Mp2c5k9Gb71j31/RvM1+wyFfdyKQ5awbXx8l3nbxsdZcciSVraby+zJ8J14+MEzmq/ZZyjs4w49aH+uPPd4Dli8gKX7L+KAxQu48tzj93iDdG+3m8vsyfAdvWwpF7z+iF3mXfD6I7zZPEDeaJ4Dzj5xOacffdiMn5jZ2+3mMnsyfO8/5zVccNorufcbd3Dbb5xmIAyYoTBHHHrQ/nv1Bra3281l9mT4jl62lNEDFxsIQ+DlI0lSw1CQJDUMBUlSw1CQJDUMBUlSw1CQJDUMBUlSw1CQJDUMBUlSw1CQJDUMBUlSw1CQJDUMBUlSw1CQJDUMBUlSw1CQJDUMBUlSw1CQJDUMBUlSw1CQJDUMBUlSw1CQJDUMBUlSw1CQJDUMBUlSw1CQJDUMBUlSw1CQJDVaDYUkZyZ5MMnGJJdPsTxJPtpdviHJSW3WI0maXmuhkGQhcBVwFnAccH6S4yasdhZwTPdnNfCJtuqRJO1Zm2cKpwAbq+rhqtoKXA+cM2Gdc4Brq+MO4OAkr2ixJknSNBa1uO/lwGM906PAqX2ssxz4Tu9KSVbTOZMAGEvy4OyW2prDgKeGXcRLjD2ZzJ5Mzb5M9mJ68hP9rNRmKGSKebUX61BVa4G1s1HUICW5u6pWDruOlxJ7Mpk9mZp9mWwQPWnz8tEocHjP9Apg016sI0kakDZD4S7gmCRHJtkPOA+4ecI6NwMXdJ9COg34QVV9Z+KOJEmD0drlo6ranuRS4BZgIXBNVd2f5OLu8jXAOuDNwEbgh8CFbdUzJPvcJa8BsCeT2ZOp2ZfJWu9JqiZdwpckzVN+olmS1DAUJEkNQ2EWJXkkyb1J1ie5uzvviiSPd+etT/LmYdc5SEkOTvL5JP+Q5JtJXp/kR5PcmuQfu/8eMuw6B2k3PZm3x0mSY3v+7vVJnk3y7vl8nEzTk9aPE+8pzKIkjwArq+qpnnlXAGNV9QfDqmuYknwW+Ouqurr7FNqBwHuA71XV73fHxDqkqi4baqEDtJuevJt5fJzs1B0e53E6H3S9hHl8nOw0oScX0vJx4pmCWpPkR4CfBf4XQFVtrarv0xne5LPd1T4L/OJwKhy8aXqijjOA/1dV32YeHycT9PakdYbC7Crgy0m+0R2aY6dLu6PAXjOfToGBo4DvAp9O8ndJrk7yMmDZzs+jdP/9sWEWOWC76wnM3+Ok13nAn3Zfz+fjpFdvT6Dl48RQmF2nV9VJdEZ/vSTJz9IZ+fVVwIl0xnT68BDrG7RFwEnAJ6rqtcDzwKQh1OeZ3fVkPh8nAHQvpZ0NfG7YtbxUTNGT1o8TQ2EWVdWm7r9PAl8ETqmqJ6pqR1WNA5+iM3rsfDEKjFbV17vTn6fzhvjEztFwu/8+OaT6hmHKnszz42Sns4B7quqJ7vR8Pk522qUngzhODIVZkuRlSZbufA28CbhvwlDgbwXuG0Z9w1BV/wQ8luTY7qwzgAfoDG/yzu68dwI3DaG8odhdT+bzcdLjfHa9TDJvj5Meu/RkEMeJTx/NkiRH0Tk7gM4lgj+pqt9Lch2dU70CHgF+fT6N75TkROBqYD/gYTpPTywA/hw4AngUeFtVfW9oRQ7YbnryUeb3cXIgnWH0j6qqH3TnHcr8Pk6m6knr7yeGgiSp4eUjSVLDUJAkNQwFSVLDUJAkNQwFSVKjtW9ekwat+wjjX3YnfxzYQWdICeh8kHDrUAqbRpJfBdZ1P78gDZ2PpGpOeimNTptkYVXt2M2yrwGXVtX6GexvUVVtn7UCpR5ePtK8kOSdSe7sjkH/8SQLkixK8v0kH0pyT5Jbkpya5CtJHt45Vn2Si5J8sbv8wSTv7XO/H0hyJ3BKkt9JcleS+5KsScfb6XwQ6c+62++XZDTJwd19n5bktu7rDyT5ZJJb6QymtyjJH3Z/94YkFw2+q5qLDAXNeUleTWdIgJ+uqhPpXDY9r7v45cCXuwMZbgWuoDP0xNuA9/fs5pTuNicB70hyYh/7vaeqTqmq24GPVNXJwGu6y86sqj8D1gNvr6oT+7i89VrgLVX1K8Bq4MmqOgU4mc4AjEfsTX+kXt5T0Hzw83TeOO9OArCEzvABAJur6tbu63uBH1TV9iT3Aq/s2cctVfUMQJIbgTfQ+e9nd/vdyj8PewJwRpLfBA4ADgO+AXxphn/HTVX1Qvf1m4B/laQ3hI6hMxyEtNcMBc0HAa6pqv+xy8xkEZ03753GgS09r3v/+5h48632sN/N1b1h1x3D5mN0RkN9PMkH6ITDVLbzz2fwE9d5fsLf9K6q+kukWeTlI80HtwG/nOQw6DyltBeXWt6UzncrH0jnG8H+Zgb7XUInZJ7qjqR7bs+y54ClPdOPAK/rvu5db6JbgHd1A2jnd/oumeHfJE3imYLmvKq6N8nvALclWQBsAy4GNs1gN18D/oTOF5xct/NpoX72W1VPp/O9zPcB3wa+3rP408DVSTbTuW9xBfCpJP8E3DlNPZ+kM3ro+u6lqyfphJX0ovhIqrQH3Sd7Xl1V7x52LVLbvHwkSWp4piBJanimIElqGAqSpIahIElqGAqSpIahIElq/H/IxmFZztFAcQAAAABJRU5ErkJggg==\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", "
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: Thu, 05 Nov 2020 Deviance: 0.22231
Time: 02:29:10 Pearson chi2: 0.236
No. Iterations: 4 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: Thu, 05 Nov 2020 Deviance: 0.22231\n", "Time: 02:29:10 Pearson chi2: 0.236\n", "No. Iterations: 4 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(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": "iVBORw0KGgoAAAANSUhEUgAAAXcAAAEKCAYAAADpfBXhAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMi4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvIxREBQAAGzdJREFUeJzt3X+UVOWd5/H3t6tBGhohoGGAJoHM4cA6UX41jUriNkYBc+KvWQ2io4k7LHEnJJPdIxs5J7OaWT1n57S7h0zWiIwyTOLR1nEVNcsG1E3HiauxQRAEhh9riDadBDGj0Noo3f3dP+6t6qrqbrq6qO6qevy8zulD3VvPfe7z7aI+dfupW7fM3RERkbBUFHsAIiJSeAp3EZEAKdxFRAKkcBcRCZDCXUQkQAp3EZEA9RvuZrbBzI6a2Rt93G9m9rdmdsjMdpnZ3MIPU0REBiKXI/eNwNLT3H8FMD3+WQncf+bDEhGRM9FvuLv7i8AfTtPkauDHHnkFGGtmEws1QBERGbjKAvQxGXg7bbklXvfb7IZmtpLo6J6qqqp5U6ZMyWuHXV1dVFSE8XaBailNodQSSh2gWpIOHDhwzN3P7a9dIcLdelnX6zUN3H09sB6gtrbWt23bltcOm5qaqK+vz2vbUqNaSlMotYRSB6iWJDP7TS7tCvEy2AKkH4LXAK0F6FdERPJUiHB/BrglPmvmQuB9d+8xJSMiIkOn32kZM3sUqAfOMbMW4E5gGIC7rwM2A18GDgEfArcO1mBFRCQ3/Ya7uy/v534HvlmwEYlIWTh16hQtLS2cPHlySPY3ZswY9u3bNyT7Gmy51DJixAhqamoYNmxYXvsoxBuqIvIJ1NLSwujRo5k6dSpmvZ1XUVgnTpxg9OjRg76fodBfLe7Ou+++S0tLC9OmTctrH2GcVyQiQ+7kyZOMHz9+SIL9k8bMGD9+/Bn9VaRwF5G8KdgHz5n+bhXuIiIB0py7iJStRCLB+eefn1retGkTU6dOLd6ASojCXUTKVlVVFTt37uzz/o6ODiorP5kxp2kZEQnKxo0buf7667nyyitZvHgxAA0NDcyfP58LLriAO++8M9X2nnvuYcaMGVx22WUsX76ce++9F4D6+nqSl0c5duxY6q+Bzs5OVq9enerrgQceALovJ3Ddddcxc+ZMbrrpJqKzxKG5uZmLL76YWbNmUVdXx4kTJ1iyZEnGi9LChQvZtWtXQX8Pn8yXNBEpqO8/u4e9rccL2ud5k87mziv/5LRt2tvbmT17NgDTpk3jqaeeAuDll19m165djBs3jq1bt3Lw4EFeffVV3J2rrrqKF198kVGjRtHY2MiOHTvo6Ohg7ty5zJs377T7e+ihhxgzZgzNzc189NFHLFy4MPUCsmPHDvbs2cOkSZNYuHAhL730EnV1dSxbtozHHnuM+fPnc/z4caqqqrjlllvYuHEja9eu5cCBA3z00UdccMEFBfitdVO4i0jZ6mta5vLLL2fcuHEAbN26la1btzJnzhwA2traOHjwICdOnODaa69l5MiRAFx11VX97m/r1q3s2rWLJ554AoD333+fgwcPMnz4cOrq6qipqQFg9uzZHD58mDFjxjBx4kTmz58PwNlnnw3Atddey8KFC2loaGDDhg18/etfP7NfRC8U7iJyxvo7wh5qo0aNSt12d9asWcM3vvGNjDZr167t83TDyspKurq6ADLONXd3fvjDH7JkyZKM9k1NTZx11lmp5UQiQUdHB+7e6z5GjhzJ5ZdfztNPP83jjz9OvlfIPR3NuYtI0JYsWcKGDRtoa2sD4MiRIxw9epRLLrmEp556ivb2dk6cOMGzzz6b2mbq1Kls374dIHWUnuzr/vvv59SpUwAcOHCADz74oM99z5w5k9bWVpqbm4Hok6kdHR0ArFixgm9/+9vMnz8/9VdGIenIXUSCtnjxYvbt28dFF10EQHV1NQ8//DBz585l2bJlzJ49m89+9rN88YtfTG1z++2389WvfpWf/OQnXHrppan1K1as4PDhw8ydOxd359xzz2XTpk197nv48OE89thjfOtb36K9vZ2qqiqef/55AObNm8fZZ5/NrbcO0rUW3b0oP/PmzfN8/fznP89721KjWkpTKLUMZh179+4dtL57c/z48UHt/8477/SGhoZB3UfS8ePH/ciRIz59+nTv7Ozss11vv2Ngm+eQsZqWEREZYo888ggLFizgnnvuGbSvDtS0jIgIcNdddw3Zvm688cYeb/AWmo7cRSRv7r1+XbIUwJn+bhXuIpKXESNG8O677yrgB4HH13MfMWJE3n1oWkZE8lJTU0NLSwvvvPPOkOzv5MmTZxR2pSSXWpLfxJQvhbuI5GXYsGF5f0tQPpqamlKfMi13Q1GLpmVERAKkcBcRCZDCXUQkQAp3EZEAKdxFRAKkcBcRCZDCXUQkQAp3EZEAKdxFRAKkcBcRCZDCXUQkQAp3EZEAKdxFRAKkcBcRCZDCXUQkQAp3EZEA5RTuZrbUzPab2SEzu6OX+8eY2bNm9rqZ7TGzWws/VBERyVW/4W5mCeA+4ArgPGC5mZ2X1eybwF53nwXUA//NzIYXeKwiIpKjXI7c64BD7v6mu38MNAJXZ7VxYLSZGVAN/AHoKOhIRUQkZ9bfN5eb2XXAUndfES/fDCxw91VpbUYDzwAzgdHAMnf/X730tRJYCTBhwoR5jY2NeQ26ra2N6urqvLYtNaqlNIVSSyh1gGpJWrRo0XZ3r+2vXS5fkG29rMt+RVgC7AQuBf4YeM7M/sndj2ds5L4eWA9QW1vr9fX1Oey+p6amJvLdttSoltIUSi2h1AGqZaBymZZpAaakLdcArVltbgWe9Mgh4NdER/EiIlIEuYR7MzDdzKbFb5LeQDQFk+4t4EsAZjYBmAG8WciBiohI7vqdlnH3DjNbBWwBEsAGd99jZrfF968D/guw0cx2E03jfNfdjw3iuEVE5DRymXPH3TcDm7PWrUu73QosLuzQREQkX/qEqohIgBTuIiIBUriLiARI4S4iEiCFu4hIgBTuIiIBUriLiARI4S4iEiCFu4hIgBTuIiIBUriLiARI4S4iEiCFu4hIgBTuIiIBUriLiARI4S4iEiCFu4hIgBTuIiIBUriLiARI4S4iEiCFu4hIgBTuIiIBUriLiARI4S4iEiCFu4hIgBTuIiIBUriLiARI4S4iEiCFu4hIgBTuIiIBUriLiARI4S4iEiCFu4hIgBTuIiIByinczWypme03s0NmdkcfberNbKeZ7TGzXxR2mCIiMhCV/TUwswRwH3A50AI0m9kz7r43rc1Y4EfAUnd/y8w+PVgDFhGR/uVy5F4HHHL3N939Y6ARuDqrzY3Ak+7+FoC7Hy3sMEVEZCDM3U/fwOw6oiPyFfHyzcACd1+V1mYtMAz4E2A08AN3/3Evfa0EVgJMmDBhXmNjY16Dbmtro7q6Oq9tS41qKU2h1BJKHaBakhYtWrTd3Wv7a9fvtAxgvazLfkWoBOYBXwKqgJfN7BV3P5Cxkft6YD1AbW2t19fX57D7npqamsh321KjWkpTKLWEUgeoloHKJdxbgClpyzVAay9tjrn7B8AHZvYiMAs4gIiIDLlc5tybgelmNs3MhgM3AM9ktXka+KKZVZrZSGABsK+wQxURkVz1e+Tu7h1mtgrYAiSADe6+x8xui+9f5+77zOxnwC6gC3jQ3d8YzIGLiEjfcpmWwd03A5uz1q3LWm4AGgo3NBERyZc+oSoiEiCFu4hIgBTuIiIBUriLiARI4S4iEiCFu4hIgBTuIiIBUriLiARI4S4iEiCFu4hIgBTuIiIBUriLiARI4S4iEiCFu4hIgBTuIiIBUriLiARI4S4iEiCFu4hIgBTuIiIBUriLiARI4S4iEiCFu4hIgBTuIiIBUriLiARI4S4iEiCFu4hIgBTuIiIBUriLiARI4S4iEiCFu4hIgBTuIiIBUriLiARI4S4iEiCFu4hIgBTuIiIByinczWypme03s0Nmdsdp2s03s04zu65wQxQRkYHqN9zNLAHcB1wBnAcsN7Pz+mj3N8CWQg9SREQGJpcj9zrgkLu/6e4fA43A1b20+xbwP4GjBRyfiIjkwdz99A2iKZal7r4iXr4ZWODuq9LaTAYeAS4FHgJ+6u5P9NLXSmAlwIQJE+Y1NjbmNei2tjaqq6vz2rbUqJbSFEotodQBqiVp0aJF2929tr92lTn0Zb2sy35FWAt81907zXprHm/kvh5YD1BbW+v19fU57L6npqYm8t221KiW0hRKLaHUAaploHIJ9xZgStpyDdCa1aYWaIyD/Rzgy2bW4e6bCjJKEREZkFzCvRmYbmbTgCPADcCN6Q3cfVrytpltJJqWUbCLiBRJv+Hu7h1mtoroLJgEsMHd95jZbfH96wZ5jCIiMkC5HLnj7puBzVnreg11d//6mQ9LRETOhD6hKiISIIW7iEiAFO4iIgFSuIuIBEjhLiISoJzOlhEZLJt2HKFhy35a32tn0tgqVi+ZwTVzJhd7WJIjPX6lS+EuRbNpxxHWPLmb9lOdABx5r501T+4GUECUAT1+pU3TMlI0DVv2p4Ihqf1UJw1b9hdpRDIQevxKm8Jdiqb1vfYBrZfSosevtCncpWgmja0a0HopLXr8SpvCXYpm9ZIZVA1LZKyrGpZg9ZIZRRqRDIQev9KmN1SlaJJvuulsi/Kkx6+0KdylqK6ZM1lhUMb0+JUuTcuIiARI4S4iEiCFu4hIgBTuIiIBUriLiARI4S4iEiCFu4hIgBTuIiIBUriLiARI4S4iEiCFu4hIgBTuIiIBUriLiARI4S4iEiCFu4hIgBTuIiIBUriLiARI4S4iEiCFu4hIgPQdqiKSoavL6XSns8txhy6Plru6onVdyXVdTpdHbTq7utt0xcvpbaKf7vVdXaS173ubzi5Sbfa9dYq3X/lNarm7PXGfcT/ueHx/97hJG0e0nN0muW16ne6Zv4vM9k6nkzEeT+/D02qNl5P7rJ8E9fWD+zgq3KVkZD8hkqHicRBET55kG9KehJlP8vQnWXKb04ZH6knXHQxvHDnFse0tWWPKDLfsfWYHTPc2PYOlOxy6t+m+Py18soIlo413B1SX9+yz050PPzzJWS+/EPeT2Wdv++3yYv8v6MfeN3JqVmFQYUZFhZEwi5YrjES8bPG6RIVRYfH6CsMMEvFytH3cT2pdtDy8siJubyTifsBIVGT22b0tGcuf+uh3g/t7IsdwN7OlwA+ABPCgu//XrPtvAr4bL7YB/97dXy/kQAdD8kmf+YSnO0RSr7ZpbdJDp7cnd3/bpx21uDu7ftfB8ddb+36SZh0Z9NwHWUcSWUHQ1d1PakwZYdMzJLt6CZLMcfUMlU532to+ZETzz+Px9AzAjPBLO5JJLpec3fn9F06GQPTETwuF7Cd7HDLpwVIRbxMFUWawGFEflRUVnFVpcX+k2ifDKn2/R3//eyZNPKc76LL67G2/6X2mB1oqsNLWRTWQ6is7BLtDMnus3WNJtkn1k94mDthEhfHKyy/zhYULM+usyBpDfNvMCvt/ocCamo4N+j76DXczSwD3AZcDLUCzmT3j7nvTmv0a+Nfu/i9mdgWwHlgwGAP+xYF3+N4vP6TqtV9kBFnGEVp6qPU4ius+qvJSyZOdO864i+wndnqYdB8xdD+Zsp8MmUcZWdvHbXo7WrG0J+CxinYm/tHYVD+ZT0LSjpjiJ3bqqCrtyCrtiZ1ZU+YTurfwSMT9ZuwjOzxS9VtakKSHV7T9tldf5eKLLszpd5PcPjn2UtLU1ER9/axiD6MgPjWignNHn1XsYZSNXI7c64BD7v4mgJk1AlcDqXB39/+b1v4VoKaQg0xXfVaCCaMqmPDp6h5BlgqKrCd1X0FD9vanOUrp7jvz6Cc9ELqDI3P/GUdyWcG0fds2Llwwv8f22WPv3nd32Gb3XWxRkMwp9jAK4u1RFUwZN7LYwxDJm3k/h69mdh2w1N1XxMs3AwvcfVUf7W8HZibbZ923ElgJMGHChHmNjY15DbqtrY3q6uq8ti01qqU0hVJLKHWAaklatGjRdnev7a9dLkfuvR0S9vqKYGaLgD8HvtDb/e6+nmjKhtraWq/P8+3i6Agxv21LjWopTaHUEkodoFoGKpdwbwGmpC3XAK3ZjczsAuBB4Ap3f7cwwxMRkXzk8iGmZmC6mU0zs+HADcAz6Q3M7DPAk8DN7n6g8MMUEZGB6PfI3d07zGwVsIXoVMgN7r7HzG6L718H/GdgPPCj+I29jlzmhEREZHDkdJ67u28GNmetW5d2ewXQ4w1UkaG2accRGrbsp/W9diaNrWL1khkAPdZdM2fykOx7MPaTi+9t2s2jv3qb73z+FH++ZjPLF0zh7mvOL8pYpDj0CVUJxqYdR1jz5G7aT3UCcOS9dlb/4+tgcKrTU+vWPLkboKDB29u+B2M/ufjept08/MpbqeVO99SyAv6TQxcOk2A0bNmfCtekU12eCvak9lOdNGzZP+j7Hoz95OLRX709oPUSJoW7BKP1vfZBaXsm/RV6P7no7OOzK32tlzAp3CUYk8ZWDUrbM+mv0PvJRaKPTyv3tV7CpHCXYKxeMoOqYYmMdcMqjGGJzFCrGpZIvdE6mPsejP3kYvmCKQNaL2HSG6oSjOQbl8U4W6avfRfjbJnkm6bJOfaEmc6W+QRSuEtQrpkzuddAHYqQ7WvfxXD3Nedz9zXn09TUxP+7qb7Yw5Ei0LSMiEiAFO4iIgFSuIuIBEjhLiISIIW7iEiAFO4iIgFSuIuIBEjhLiISIIW7iEiAFO4iIgFSuIuIBEjhLiISIIW7iEiAFO4iIgFSuIuIBEjhLiISIIW7iEiAFO4iIgFSuIuIBEjhLiISIIW7iEiAFO4iIgFSuIuIBEjhLiISIIW7iEiAFO4iIgFSuIuIBEjhLiISoJzC3cyWmtl+MztkZnf0cr+Z2d/G9+8ys7mFH6qIiOSq33A3swRwH3AFcB6w3MzOy2p2BTA9/lkJ3F/gcYqIyADkcuReBxxy9zfd/WOgEbg6q83VwI898gow1swmFnisIiKSo8oc2kwG3k5bbgEW5NBmMvDb9EZmtpLoyB6gzcz2D2i03c4BjuW5balRLaUplFpCqQNUS9Jnc2mUS7hbL+s8jza4+3pgfQ77PP2AzLa5e+2Z9lMKVEtpCqWWUOoA1TJQuUzLtABT0pZrgNY82oiIyBDJJdybgelmNs3MhgM3AM9ktXkGuCU+a+ZC4H13/212RyIiMjT6nZZx9w4zWwVsARLABnffY2a3xfevAzYDXwYOAR8Ctw7ekIECTO2UENVSmkKpJZQ6QLUMiLn3mBoXEZEyp0+oiogESOEuIhKgkg93MxthZq+a2etmtsfMvh+vH2dmz5nZwfjfTxV7rLkws4SZ7TCzn8bL5VrHYTPbbWY7zWxbvK5caxlrZk+Y2T+b2T4zu6gcazGzGfHjkfw5bmbfKdNa/kP8fH/DzB6Nc6Ds6gAws7+M69hjZt+J1w16LSUf7sBHwKXuPguYDSyNz8i5A3jB3acDL8TL5eAvgX1py+VaB8Aid5+ddr5uudbyA+Bn7j4TmEX0+JRdLe6+P348ZgPziE5ueIoyq8XMJgPfBmrd/fNEJ3LcQJnVAWBmnwf+HdEn/WcBXzGz6QxFLe5eNj/ASOA1ok/I7gcmxusnAvuLPb4cxl8TP5CXAj+N15VdHfFYDwPnZK0ru1qAs4FfE59cUM61ZI1/MfBSOdZC9yfexxGd0ffTuJ6yqiMe5/XAg2nLfwX8p6GopRyO3JNTGTuBo8Bz7v4rYILH59LH/366mGPM0VqiB7YrbV051gHRJ5C3mtn2+LISUJ61fA54B/j7eLrsQTMbRXnWku4G4NH4dlnV4u5HgHuBt4guYfK+u2+lzOqIvQFcYmbjzWwk0SnjUxiCWsoi3N2906M/NWuAuvhPnbJiZl8Bjrr79mKPpUAWuvtcoiuCftPMLin2gPJUCcwF7nf3OcAHlMGf+6cTf9jwKuAfiz2WfMTzz1cD04BJwCgz+7Pijio/7r4P+BvgOeBnwOtAx1DsuyzCPcnd3wOagKXA75NXnoz/PVrEoeViIXCVmR0murLmpWb2MOVXBwDu3hr/e5RoXreO8qylBWiJ/xoEeIIo7MuxlqQrgNfc/ffxcrnVchnwa3d/x91PAU8CF1N+dQDg7g+5+1x3vwT4A3CQIail5MPdzM41s7Hx7SqiB/6fiS558LW42deAp4szwty4+xp3r3H3qUR/Mv8fd/8zyqwOADMbZWajk7eJ5kPfoAxrcfffAW+b2Yx41ZeAvZRhLWmW0z0lA+VXy1vAhWY20syM6DHZR/nVAYCZfTr+9zPAnxI9NoNeS8l/QtXMLgD+gegd8wrgcXf/azMbDzwOfIboP8P17v6H4o00d2ZWD9zu7l8pxzrM7HNER+sQTWs84u73lGMtAGY2G3gQGA68SXT5jArKs5aRRG9Gfs7d34/Xld3jEp/yvIxoCmMHsAKopszqADCzfwLGA6eA/+juLwzFY1Ly4S4iIgNX8tMyIiIycAp3EZEAKdxFRAKkcBcRCZDCXUQkQLl8QbbIkIpPE3shXvwjoJPoEgEAde7+cVEGdhpm9m+BzfF58yJFp1MhpaSZ2V1Am7vfWwJjSbh7Zx/3/RJY5e47B9BfpbsPyUfR5ZNH0zJSVszsaxZd33+nmf3IzCrMrNLM3jOzBjN7zcy2mNkCM/uFmb1pZl+Ot11hZk/F9+83s+/l2O/dZvYq0XWNvm9mzfH1uddZZBnR5agfi7cfbmYtaZ+svtDMno9v321mD5jZc0QXK6s0s/8e73uXma0Y+t+qhEjhLmUjvmDctcDF8YXkKoku5QAwBtgaX8zsY+Auoo+tXw/8dVo3dfE2c4EbzWx2Dv2+5u517v4y8AN3nw+cH9+31N0fA3YCyzy6nnp/00ZzgCvd/WZgJdEF5eqA+UQXYftMPr8fkXSac5dychlRAG6LLjlCFdFH7QHa3f25+PZuosvEdpjZbmBqWh9b3P1fAMxsE/AFoudBX/1+TPelFgC+ZGargRHAOcB24H8PsI6n3f1kfHsx8K/MLP3FZDrRR9JF8qZwl3JiwAZ3/6uMlWaVRCGc1EX0DV7J2+n/z7PfZPJ++m33+I2p+Lot/wOY6+5HzOxuopDvTQfdfxlnt/kgq6a/cPcXECkgTctIOXke+KqZnQPRWTV5TGEstug7U0cSXTP8pQH0W0X0YnEsvirmv0m77wQwOm35MNFX3ZHVLtsW4C/iF5Lk96BWDbAmkR505C5lw913x1cLfN7MKoiusncb0DqAbn4JPAL8MfCT5NktufTr7u+a2T8QXd74N8Cv0u7+e+BBM2snmte/C/g7M/sd8OppxvMA0ZUBd8ZTQkeJXnREzohOhZRPjPhMlM+7+3eKPRaRwaZpGRGRAOnIXUQkQDpyFxEJkMJdRCRACncRkQAp3EVEAqRwFxEJ0P8HfLcy7/zjy3oAAAAASUVORK5CYII=\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(\"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": [ "### New Analysis\n", "

Instead of ignoring the cases where malfunction = 0, we will use all the data

" ] }, { "cell_type": "code", "execution_count": 7, "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYUAAAEKCAYAAAD9xUlFAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMi4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvIxREBQAAF9JJREFUeJzt3X2UXXV97/H3d5IACYlAg02VQAFJsVyBCOFJtDfx6Qa7JPUCBbyCl940ZUlul9y2htvVa6m1a1V8qHpFY+SiQldNVRBom14e1Ii0IASM4UHBuYBhEhogBshASGYy3/vH2bN7Mkxmzhlmz5lzeL/WmpWz9/mdne939pz5zN5nn9+JzESSJICuVhcgSZo8DAVJUslQkCSVDAVJUslQkCSVDAVJUqmyUIiIqyPiqYh4YC/3R0R8PiK6I2JDRJxQVS2SpMZUeaTwNWDxCPefAcwrvpYBX6qwFklSAyoLhcy8HfjlCEOWANdkzV3AgRHxuqrqkSSNbmoL/+9DgCfqlnuKdU8OHRgRy6gdTTB9+vQTDz300AkpsFEDAwN0dXXmyzOd2pt9tZ9O7W2i+nrkkUeeyczXjjaulaEQw6wbds6NzFwFrAJYsGBBrlu3rsq6mrZ27VoWLlzY6jIq0am92Vf76dTeJqqviPhFI+NaGbs9QP2f/HOBzS2qRZJEa0PhJuDC4iqkU4HnMvNlp44kSROnstNHEfENYCFwcET0AH8OTAPIzJXAGuA9QDfwInBRVbVIkhpTWShk5vmj3J/AJVX9/5Kk5nXeS/mSpDEzFCRJJUNBklQyFCRJJUNBklQyFCRJJUNBklQyFCRJJUNBklQyFCRJJUNBklQyFCRJJUNBklQyFCRJJUNBklQyFCRJJUNBklQyFCRJJUNBklQyFCRJJUNBklQyFCRJJUNBklQyFCRJJUNBklQyFCRJJUNBklQyFCRJJUNBklQyFCRJJUNBklQyFCRJJUNBklQyFCRJJUNBklSqNBQiYnFEPBwR3RFx2TD3HxAR/xARP4mIByPioirrkSSNrLJQiIgpwJXAGcAxwPkRccyQYZcAD2Xm8cBC4NMRsU9VNUmSRlblkcLJQHdmPpqZu4DVwJIhYxKYFREBzAR+CfRXWJMkaQSRmdVsOOJsYHFmLi2WLwBOyczldWNmATcBbwRmAedm5j8Ns61lwDKAOXPmnLh69epKah6r3t5eZs6c2eoyKtGpvdlX++nU3iaqr0WLFt2bmQtGGze1whpimHVDE+g/AeuBtwNvAG6NiB9m5vN7PChzFbAKYMGCBblw4cLxr/YVWLt2LZOtpvHSqb3ZV/vp1N4mW19Vnj7qAQ6tW54LbB4y5iLg+qzpBh6jdtQgSWqBKkPhHmBeRBxRvHh8HrVTRfU2Au8AiIg5wNHAoxXWJEkaQWWnjzKzPyKWAzcDU4CrM/PBiLi4uH8l8JfA1yLifmqnm1Zk5jNV1SRJGlmVrymQmWuANUPWray7vRl4d5U1SJIa5zuaJUklQ0GSVDIUJEklQ0GSVDIUJEklQ0GSVDIUJEklQ0GSVDIUJEklQ0GSVDIUJEklQ0GSVDIUJEklQ0GSVDIUJEklQ0GSVDIUJEklQ0GSVDIUJEklQ0GSVDIUJEklQ0GSVDIUJEklQ0GSVDIUJEklQ0GSVDIUJEklQ0GSVDIUJEklQ0GSVDIUJEklQ0GSVDIUJEklQ0GSVKo0FCJicUQ8HBHdEXHZXsYsjIj1EfFgRPygynokSSOb2sigiHhTZj7QzIYjYgpwJfAuoAe4JyJuysyH6sYcCHwRWJyZGyPiV5v5PyRJ46vRI4WVEXF3RHyo+EXeiJOB7sx8NDN3AauBJUPGvB+4PjM3AmTmUw1uW5JUgcjMxgZGzAN+DzgHuBv4ambeOsL4s6kdASwtli8ATsnM5XVjPgtMA/4DMAv4XGZeM8y2lgHLAObMmXPi6tWrG+tugvT29jJz5sxWl1GJTu3NvtpPp/Y2UX0tWrTo3sxcMOrAzGz4C5gCnAVsAn4K/Az4z3sZew5wVd3yBcD/HjLmC8BdwP7AwcDPgd8YqYYTTzwxJ5vvf//7rS6hMp3am321n07tbaL6AtZlA7/nG31N4TjgIuC3gVuB92bmfRHxeuBO4PphHtYDHFq3PBfYPMyYZzLzBeCFiLgdOB54pJG6JEnjq9HXFL4A3Accn5mXZOZ9AJm5GfizvTzmHmBeRBwREfsA5wE3DRlzI/C2iJgaETOAU6gdgUiSWqChIwXgPcCOzNwNEBFdwH6Z+WJmXjvcAzKzPyKWAzdTO+10dWY+GBEXF/evzMyfRsT/BTYAA9RONzV1lZMkafw0Ggq3Ae8EeovlGcAtwFtGelBmrgHWDFm3csjyJ4FPNliHJKlCjZ4+2i8zBwOB4vaMakqSJLVKo6HwQkScMLgQEScCO6opSZLUKo2ePvow8K2IGLx66HXAudWUJElqlYZCITPviYg3AkcDAfwsM/sqrUySNOEaPVIAOAk4vHjMmyOCHObdx5Kk9tXom9euBd4ArAd2F6sTMBQkqYM0eqSwADimeKu0JKlDNXr10QPAr1VZiCSp9Ro9UjgYeCgi7gZ2Dq7MzDMrqUqS1BKNhsLlVRYhSZocGr0k9QcR8evAvMy8rZi8bkq1pUmSJlpDrylExO8D3wa+XKw6BLihqqIkSa3R6AvNlwCnA88DZObPAT9PWZI6TKOhsDNrn7MMQERMpfY+BUlSB2k0FH4QEX8KTI+IdwHfAv6hurIkSa3QaChcBjwN3A/8AbXPSNjbJ65JktpUo1cfDQBfKb4kSR2q0bmPHmOY1xAy88hxr0iS1DLNzH00aD/gHOBXxr8cSVIrNfSaQmZurfvalJmfBd5ecW2SpAnW6OmjE+oWu6gdOcyqpCJJUss0evro03W3+4HHgd8d92okSS3V6NVHi6ouRJLUeo2ePvofI92fmZ8Zn3IkSa3UzNVHJwE3FcvvBW4HnqiiKElSazTzITsnZOZ2gIi4HPhWZi6tqjBJ0sRrdJqLw4Bddcu7gMPHvRpJUks1eqRwLXB3RHyH2jub3wdcU1lVkqSWaPTqo7+KiH8G3lasuigzf1xdWZKkVmj09BHADOD5zPwc0BMRR1RUkySpRRr9OM4/B1YA/7NYNQ3426qKkiS1RqNHCu8DzgReAMjMzTjNhSR1nEZDYVdmJsX02RGxf3UlSZJapdFQ+GZEfBk4MCJ+H7gNP3BHkjpOo1cffar4bObngaOBj2bmrZVWJkmacKMeKUTElIi4LTNvzcw/ycw/bjQQImJxRDwcEd0RcdkI406KiN0RcXYzxUuSxteooZCZu4EXI+KAZjYcEVOAK4EzgGOA8yPimL2M+wRwczPblySNv0bf0fwScH9E3EpxBRJAZv7hCI85GejOzEcBImI1sAR4aMi4/w5cR23CPUlSCzUaCv9UfDXjEPacRbUHOKV+QEQcQu1y17czQihExDJgGcCcOXNYu3Ztk6VUq7e3d9LVNF46tTf7aj+d2ttk62vEUIiIwzJzY2Z+fQzbjmHW5ZDlzwIrMnN3xHDDiwdlrgJWASxYsCAXLlw4hnKqs3btWiZbTeOlU3uzr/bTqb1Ntr5Ge03hhsEbEXFdk9vuAQ6tW54LbB4yZgGwOiIeB84GvhgRv9Pk/yNJGiejnT6q//P9yCa3fQ8wr5gjaRNwHvD++gGZWc6fFBFfA/4xM29AktQSo4VC7uX2qDKzPyKWU7uqaApwdWY+GBEXF/evbKpSSVLlRguF4yPieWpHDNOL2xTLmZmvGenBmbkGWDNk3bBhkJn/taGKJUmVGTEUMnPKRBUiSWq9Zj5PQZLU4QwFSVLJUJAklQwFSVLpVRMKW3t38pMnnmVr785WlyKpSVt7d7Kjb7fP3wnwqgiFG9dv4vRPfI8PXPUjTv/E97hp/aZWlySpQYPP38eefsHn7wTo+FDY2ruTFddt4KW+Abbv7OelvgE+ct0G/+KQ2kD983d3ps/fCdDxodCzbQfTuvZsc1pXFz3bdrSoIkmN8vk78To+FOYeNJ2+gYE91vUNDDD3oOktqkhSo3z+TryOD4XZM/flirOOY79pXczadyr7TeviirOOY/bMfVtdmqRR1D9/p0T4/J0AjX7ITls7c/4hnH7UwfRs28Hcg6b7AyW1kcHn79133sG/nPlWn78Ve1WEAtT+4vCHSWpPs2fuy/RpU3wOT4COP30kSWqcoSBJKhkKkqSSoSBJKhkKkqSSoSBJKhkKkqSSoSBJKhkKkqSSoSBJKhkKkqSSoSBJKhkKkqSSoSBJKhkKkqSSoSBJKhkKkqSSoSBJKhkKkqSSoSBJKhkKkqRSpaEQEYsj4uGI6I6Iy4a5/79ExIbi618j4vgq65EkjayyUIiIKcCVwBnAMcD5EXHMkGGPAf8xM48D/hJYVVU9kqTRVXmkcDLQnZmPZuYuYDWwpH5AZv5rZm4rFu8C5lZYjyRpFJGZ1Ww44mxgcWYuLZYvAE7JzOV7Gf/HwBsHxw+5bxmwDGDOnDknrl69upKax6q3t5eZM2e2uoxKdGpv9tV+OrW3iepr0aJF92bmgtHGTa2whhhm3bAJFBGLgP8GvHW4+zNzFcWppQULFuTChQvHqcTxsXbtWiZbTeOlU3uzr/bTqb1Ntr6qDIUe4NC65bnA5qGDIuI44CrgjMzcWmE9kqRRVPmawj3AvIg4IiL2Ac4DbqofEBGHAdcDF2TmIxXWIklqQGVHCpnZHxHLgZuBKcDVmflgRFxc3L8S+CgwG/hiRAD0N3LOS5JUjSpPH5GZa4A1Q9atrLu9FHjZC8uCrb076dm2g7kHTWf2zH3HbWw76dS+qtK9ZTvbXuyje8t2jpozq9XlqE1VGgoamxvXb2LFdRuY1tVF38AAV5x1HGfOP+QVj20nndpXVT56w/1cc9dG/ujYfi79m9u58LTD+NiSY1tdltqQ01xMMlt7d7Liug281DfA9p39vNQ3wEeu28DW3p2vaGw76dS+qtK9ZTvX3LVxj3XX3LmR7i3bW1SR2pmhMMn0bNvBtK49d8u0ri56tu14RWPbSaf2VZX1Tzzb1HppJIbCJDP3oOn0DQzssa5vYIC5B01/RWPbSaf2VZX5hx7Y1HppJIbCJDN75r5ccdZx7Deti1n7TmW/aV1ccdZxw77Q2szYdtKpfVXlqDmzuPC0w/ZYd+Fph/lis8bEF5onoTPnH8LpRx3c0JU3zYxtJ53aV1U+tuRYLjz1cO6/9y5uu/RUA0FjZihMUrNn7tvwL8JmxraTTu2rKkfNmUXPjGkGgl4RTx9JkkqGgiSpZChIkkqGgiSpZChIkkqGgiSpZChIkkqGgiSpZChIkkqGgiSpZChIkkqGgiSpZChIkkqGgiSpZChIkkqGgiSpZChIkkqGgiSpZChIkkqGgiSpZChIkkqGgiSpZChIkkqGgiSpZChIkkqGgiSpZChIkkqGgiSpVGkoRMTiiHg4Iroj4rJh7o+I+Hxx/4aIOKHKeqRmbe3dyU+eeJatvTtHHbvusa185paHWffY1nHbZjNju7dsZ9uLfXRv2T7q2GZUVW+zNezo2z3qdru3bOfb657o2O9BFdsdampVG46IKcCVwLuAHuCeiLgpMx+qG3YGMK/4OgX4UvGv1HI3rt/Eius2MK2ri76BAa446zjOnH/IsGM/cNVd3NFdC4PPf6+btx01m2uXnvqKttnM2I/ecD/X3LWRPzq2n0v/5nYuPO0wPrbk2DF2Xn29Y6nhD3+zj0s/8b29bnfwezCoE78H473d4VR5pHAy0J2Zj2bmLmA1sGTImCXANVlzF3BgRLyuwpqkhmzt3cmK6zbwUt8A23f281LfAB+5bsOwf6Wte2xrGQiDfti99WVHDM1ss5mx3Vu27/HLEOCaOze+4r+Wq6p3rDXsztzrdl8t34Px3O7eRGZWs+GIs4HFmbm0WL4AOCUzl9eN+UfgrzPzjmL5u8CKzFw3ZFvLgGXF4tHAw5UUPXYHA8+0uoiKdGpvI/YV06bPmHrQ634jurqmDK7LgYHd/duefCT7drxYP3bKrINfP2X/A1/2x8zuF559cvf2ZzaPZZvNjO2accDsqa957eEAu198jikzDgCg//mnHx948bmRz2WNoKp6x1rDYG/Dbbf+e1CvTb4H4/azOIpfz8zXjjaostNHQAyzbmgCNTKGzFwFrBqPoqoQEesyc0Gr66hCp/bWyX31P/dUx/UFndvbZPtZrPL0UQ9waN3yXGDzGMZIkiZIlaFwDzAvIo6IiH2A84Cbhoy5CbiwuArpVOC5zHyywpokSSOo7PRRZvZHxHLgZmAKcHVmPhgRFxf3rwTWAO8BuoEXgYuqqqdik/bU1jjo1N7sq/10am+Tqq/KXmiWJLUf39EsSSoZCpKkkqEwBhHxeETcHxHrI2Jdse7yiNhUrFsfEe9pdZ3NiogDI+LbEfGziPhpRJwWEb8SEbdGxM+Lfw9qdZ3N2ktfnbC/jq6rf31EPB8RH273fTZCX52wzy6NiAcj4oGI+EZE7DfZ9pevKYxBRDwOLMjMZ+rWXQ70ZuanWlXXKxURXwd+mJlXFVeMzQD+FPhlZv51MX/VQZm5oqWFNmkvfX2YNt9f9YppZTZRmybmEtp8nw0a0tdFtPE+i4hDgDuAYzJzR0R8k9rFNscwifaXRwoCICJeA/wW8H8AMnNXZj5LbSqSrxfDvg78TmsqHJsR+uo07wD+X2b+gjbfZ0PU99UJpgLTI2IqtT9ONjPJ9pehMDYJ3BIR9xZTcAxaXsz2enWrDwHH4EjgaeCrEfHjiLgqIvYH5gy+d6T491dbWeQY7K0vaO/9NdR5wDeK2+2+z+rV9wVtvM8ycxPwKWAj8CS192XdwiTbX4bC2JyemSdQm+X1koj4LWozvL4BmE9th3+6hfWNxVTgBOBLmflm4AXgZdOdt6G99dXu+6tUnBI7E/hWq2sZT8P01db7rAixJcARwOuB/SPiA62t6uUMhTHIzM3Fv08B3wFOzswtmbk7MweAr1CbJbad9AA9mfmjYvnb1H6Zbhmcubb496kW1TdWw/bVAfur3hnAfZm5pVhu9302aI++OmCfvRN4LDOfzsw+4HrgLUyy/WUoNCki9o+IWYO3gXcDD8SeU36/D3igFfWNVWb+G/BERBxdrHoH8BC1qUg+WKz7IHBjC8obs7311e77a4jz2fMUS1vvszp79NUB+2wjcGpEzIiIoPaz+FMm2f7y6qMmRcSR1I4OoHZq4u8y868i4lpqh7UJPA78QbvN4xQR84GrgH2AR6ld7dEFfBM4jNoP9TmZ+cuWFTkGe+nr87T5/gKIiBnAE8CRmflcsW427b/PhuurE55jfwGcC/QDPwaWAjOZRPvLUJAklTx9JEkqGQqSpJKhIEkqGQqSpJKhIEkqVfbJa9JEKy7F/G6x+GvAbmpTXEDtDYa7WlLYCCLi94A1xfsppJbzklR1pMk0a21ETMnM3Xu57w5geWaub2J7UzOzf9wKlOp4+kivChHxwYi4u5iH/4sR0RURUyPi2Yj4ZETcFxE3R8QpEfGDiHh0cL7+iFgaEd8p7n84Iv6swe1+PCLuBk6OiL+IiHuKefRXRs251N6M9ffF4/eJiJ6IOLDY9qkRcVtx++MR8eWIuJXa5H5TI+Izxf+9ISKWTvx3VZ3IUFDHi4g3UZsW4S2ZOZ/aadPzirsPAG4pJjjcBVxObfqBc4CP1W3m5OIxJwDvj4j5DWz3vsw8OTPvBD6XmScBxxb3Lc7MvwfWA+dm5vwGTm+9GXhvZl4ALAOeysyTgZOoTcx42Fi+P1I9X1PQq8E7qf3iXFebcobp1KZQANiRmbcWt++nNp1xf0TcDxxet42bM3MbQETcALyV2vNnb9vdxb9PhwLwjoj4E2A/4GDgXuCfm+zjxsx8qbj9buA3I6I+hOZRmyZBGjNDQa8GAVydmf9rj5W1Dzqp/+t8ANhZd7v++TH0xbccZbs7snjBrpjH5wvUZmfdFBEfpxYOw+nn34/gh455YUhPH8rM7yKNI08f6dXgNuB3I+JgqF2lNIZTLe+O2mc9z6A2J/6/NLHd6dRC5pliht2z6u7bDsyqW34cOLG4XT9uqJuBDxUBNPi5xtOb7El6GY8U1PEy8/5idsrbIqIL6AMupvZRiI26A/g7ah/ycu3g1UKNbDczt0btc6IfAH4B/Kju7q8CV0XEDmqvW1wOfCUi/g24e4R6vkxtVs31xamrp6iFlfSKeEmqNIriyp43ZeaHW12LVDVPH0mSSh4pSJJKHilIkkqGgiSpZChIkkqGgiSpZChIkkr/HzHofwgP0tIHAAAAAElFTkSuQmCC\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "data = pd.read_csv(\"shuttle.csv\")\n", "%matplotlib inline\n", "pd.set_option('mode.chained_assignment',None) # this removes a useless warning from pandas\n", "import matplotlib.pyplot as plt\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": [ "

As we can see using visual inspection there is some tendency the tests with temperatures between 65 and 90 have less failure. Of course we can't really conclude using only this observation since we don't have enough data for the lower temperatures. However, we can see clearly that there were 2 failures occured on the lowest temperature which is 53, and there is no successful experiment below 65. This visual inspection should have raised suspicion before launching the challenger shuttle.

" ] }, { "cell_type": "code", "execution_count": 8, "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXcAAAEKCAYAAADpfBXhAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMi4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvIxREBQAAIABJREFUeJzt3Xl4VOXd//H3dyb7QmLYISA7yA5hEXEBrYK2KiriinVBpHWp7SNVn199tE+16oNt1VZxQ3GpgisupYJa44JbQBBkX8UEkJ0kkD33748ZMGAgQzLJLPm8rivXzDlzn3O+dwY+c3LmnPuYcw4REYkunlAXICIiwadwFxGJQgp3EZEopHAXEYlCCncRkSikcBcRiUI1hruZPW1mW83s28O8bmb2sJmtMbPFZjYw+GWKiMjRCGTPfTow+givnwl09f9MBKbWvSwREamLGsPdOfcxsPMITc4FnnM+XwDpZtY6WAWKiMjRiwnCOtoC31eZzvXP23xoQzObiG/vnsTExKx27drVaoOVlZV4PNHxdYH6Ep6ipS/R0g9QX/ZbtWrVdudc85raBSPcrZp51Y5p4Jx7AngCYNCgQW7+/Pm12mB2djYjRoyo1bLhRn0JT9HSl2jpB6gv+5nZd4G0C8bHYC5QdRc8E9gUhPWKiEgtBSPc3wKu8J81czywxzn3k0MyIiLScGo8LGNmLwEjgGZmlgvcCcQCOOceA2YDZwFrgH3AVfVVrIiIBKbGcHfOXVLD6w64PmgViUhEKCsrIzc3l+Li4gbZXlpaGsuXL2+QbdW3QPqSkJBAZmYmsbGxtdpGML5QFZFGKDc3l9TUVDp06IBZdedVBFdBQQGpqan1vp2GUFNfnHPs2LGD3NxcOnbsWKttRMd5RSLS4IqLi2natGmDBHtjY2Y0bdq0Tn8VKdxFpNYU7PWnrr9bhbuISBTSMXcRiVher5c+ffocmJ41axYdOnQIXUFhROEuIhErMTGRRYsWHfb18vJyYmIaZ8zpsIyIRJXp06dz4YUXcvbZZ3PGGWcAMGXKFAYPHkzfvn258847D7S955576N69Oz/72c+45JJLeOCBBwAYMWIE+4dH2b59+4G/BioqKpg8efKBdT3++OPAj8MJjB07lh49enDZZZfhO0sccnJyOOGEE+jXrx9DhgyhoKCAUaNGHfShNHz4cBYvXhzU30Pj/EgTkaD649tLWbYpP6jr7NmmCXee3euIbYqKiujfvz8AHTt25I033gDg888/Z/HixWRkZDB37lxWr17NV199hXOOc845h48//pjk5GRmzJjBwoULKS8vZ+DAgWRlZR1xe9OmTSMtLY2cnBxKSkoYPnz4gQ+QhQsXsnTpUtq0acPw4cOZN28eQ4YM4aKLLmLmzJkMHjyY/Px8EhMTueKKK5g+fToPPvggq1atoqSkhL59+wbht/YjhbuIRKzDHZY5/fTTycjIAGDu3LnMnTuXAQMGAFBYWMjq1aspKCjgvPPOIykpCYBzzjmnxu3NnTuXxYsX8+qrrwKwZ88eVq9eTVxcHEOGDCEzMxOA/v37s2HDBtLS0mjdujWDBw8GoEmTJgCcd955DB8+nClTpvD0009z5ZVX1u0XUQ2Fu4jUWU172A0tOTn5wHPnHLfffjvXXXfdQW0efPDBw55uGBMTQ2VlJcBB55o75/j73//OqFGjDmqfnZ1NfHz8gWmv10t5eTnOuWq3kZSUxOmnn86bb77Jyy+/TG1HyD0SHXMXkag2atQonn76aQoLCwHIy8tj69atnHzyybzxxhsUFRVRUFDA22+/fWCZDh06sGDBAoADe+n71zV16lTKysoAWLVqFXv37j3stnv06MGmTZvIyckBfFemlpeXAzBhwgRuuukmBg8efOCvjGDSnruIRLUzzjiD5cuXM2zYMABSUlJ44YUXGDhwIBdddBH9+/fn2GOP5aSTTjqwzC233MK4ceN4/vnnOfXUUw/MnzBhAhs2bGDgwIE452jevDmzZs067Lbj4uKYOXMmN954I0VFRSQmJvL+++8DkJWVRZMmTbjqqnoaa9E5F5KfrKwsV1sffvhhrZcNN+pLeIqWvtRnP5YtW1Zv665Ofn5+va7/zjvvdFOmTKnXbeyXn5/v8vLyXNeuXV1FRcVh21X3OwbmuwAyVodlREQa2IsvvsjQoUO555576u3WgTosIyIC3HXXXQ22rUsvvfQnX/AGm/bcRaTWnKv2dskSBHX93SrcRaRWEhIS2LFjhwK+Hjj/eO4JCQm1XocOy4hIrWRmZpKbm8u2bdsaZHvFxcV1CrtwEkhf9t+JqbYU7iJSK7GxsbW+S1BtZGdnH7jKNNI1RF90WEZEJAop3EVEopDCXUQkCincRUSikMJdRCQKKdxFRKKQwl1EJAop3EVEopDCXUQkCincRUSiUMSF+77Sct7bUEZ5RWWoSxERCVsRF+7vLN7MP1eUMu7xz/lux+HvXSgi0phFXLiPG9SOSX3jWbO1kDMf+oSZORs15KiIyCEiLtwBjm8Tw7s3n0z/dunc+toSbnhxIXuKykJdlohI2IjIcAdok57IC9cM5dbRPZizdAtnPfQJX2/cFeqyRETCQsSGO4DHY/xqRGde/dUJeDww7rHPefLjdTpMIyKNXkDhbmajzWylma0xs9uqeT3NzN42s2/MbKmZXRX8Ug+vf7t03rnxJE47rgX3zF7OxOcX6DCNiDRqNYa7mXmBR4AzgZ7AJWbW85Bm1wPLnHP9gBHAX8wsLsi1HlFaYiyPXZ7F//yiJx+u2Mq5//iUFVvyG7IEEZGwEcie+xBgjXNunXOuFJgBnHtIGwekmpkBKcBOoDyolQbAzLj6xI7MmHg8+0orGPPIPN76ZlNDlyEiEnJW0/FpMxsLjHbOTfBPjweGOuduqNImFXgL6AGkAhc55/5VzbomAhMBWrZsmTVjxoxaFV1YWEhKSsoR2+wuqeTRRSWs2lXJWR1jGdstFo9ZrbZXnwLpS6RQX8JPtPQD1Jf9Ro4cucA5N6jGhs65I/4AFwJPVZkeD/z9kDZjgb8BBnQB1gNNjrTerKwsV1sffvhhQO1Kyirc/3tjsTv21nfcFdO+dHuKSmu9zfoSaF8igfoSfqKlH86pL/sB810Nue2cC+iwTC7Qrsp0JnDosY6rgNf9217jD/ceAay7XsXFeLh7TB/uPb8P89Zs5/xHP9NVrSLSKAQS7jlAVzPr6P+S9GJ8h2Cq2gicBmBmLYHuwLpgFloXlwxpz/PXDGV7YQnnPjKPL9btCHVJIiL1qsZwd86VAzcAc4DlwMvOuaVmNsnMJvmb/Qk4wcyWAB8AtzrnttdX0bUxrHNT3rx+OE2T4xg/7UveWJgb6pJEROpNTCCNnHOzgdmHzHusyvNNwBnBLS34jm2azOu/Gs6kFxbw25nfsHFHETed1gULwy9aRUTqIqKvUK2NtKRYnr16COcPbMvf3l/Fba8t0fDBIhJ1AtpzjzZxMR7+cmE/MtMTefg/a9hWWMI/Lh1AUlyj/HWISBRqdHvu+5kZvzujO/ec15vslVu59Mkv2bW3NNRliYgERaMN9/0uG3osUy/PYtnmfC58/HM27S4KdUkiInXW6MMdYFSvVjx71RC27Clm7NTPWLetMNQliYjUicLdb1jnpsyYeDzF5ZWMe/xzlm/WoGMiErkU7lX0bpvGy9cNI8bj4aLHP2ehbv4hIhFK4X6ILi1SeGXSMNKT4hg/7StyNuwMdUkiIkdN4V6NdhlJvHzdMFo0ieeKaV/x2ZqwuthWRKRGCvfDaJWWwMyJw2ifkcRV03P4eNW2UJckIhIwhfsRNE+N56WJx9OxWTITnpuvgBeRiKFwr0FGchwvXns8nZunMOG5+XykgBeRCKBwD0BGchwvThhK5+YpTHxuPvN0DF5EwpzCPUDHJMfxzwlD6dA0mWuezdGY8CIS1hTuRyEjOY5/XjuUzGOSuHp6Dgu+02mSIhKeFO5HqVlKPC9eO5SWTRK48ukcvs3bE+qSRER+QuFeCy1SE/jnhKE0SYxl/LQvWbmlINQliYgcROFeS23SE3nx2qHExXi4fNqXuvG2iIQVhXsdHNs0mReuGUp5RSWXPfUlW/YUh7okERFA4V5nXVum8uzVQ9i9r4zLp+mGHyISHhTuQdA3M50nrxjExp37uHJ6DntLykNdkog0cgr3IBnWuSn/uGQAS3J3M+mFBZSW66bbIhI6CvcgOqNXK+47vy+frN7Of73yDZWVLtQliUgjFRPqAqLNuMHt2LG3lPvfXUHzlHju+MVxmFmoyxKRRkbhXg8mndKJrQXFPD1vPS2bxHPdKZ1DXZKINDIK93pgZtzx855sKyjh3n+voEWTeM4bkBnqskSkEVG41xOPx/jLuH7sKCzl968upkVqAsO7NAt1WSLSSOgL1XoUH+PlsfFZdGqWwqTnF7B8c36oSxKRRkLhXs/SEmN55qrBJMfHcNUzOWzeUxTqkkSkEVC4N4A26Yk8c9VgCkvKueqZHAqKy0JdkohEOYV7AzmudRMevWwgq7cWcv2LCymr0EVOIlJ/FO4N6ORuzfnzeb35eNU2/ufNpTini5xEpH7obJkGdtHg9ny3Yx+PZq+lY7MkuoW6IBGJStpzD4FbzujOz/u05t5/r2D+Fg0yJiLBF1C4m9loM1tpZmvM7LbDtBlhZovMbKmZfRTcMqPL/nPg+7dL54nFJSzO3R3qkkQkytQY7mbmBR4BzgR6ApeYWc9D2qQDjwLnOOd6ARfWQ61RJSHWyxPjB5EaZ0x4dr5OkRSRoApkz30IsMY5t845VwrMAM49pM2lwOvOuY0AzrmtwS0zOjVPjee3WQnsK63gmunzNQ68iASN1XTGhpmNBUY75yb4p8cDQ51zN1Rp8yAQC/QCUoGHnHPPVbOuicBEgJYtW2bNmDGjVkUXFhaSkpJSq2XDTWFhIeuKEvjbghL6t/By44B4PBE6imS0vS/R0Jdo6QeoL/uNHDlygXNuUE3tAjlbprqkOfQTIQbIAk4DEoHPzewL59yqgxZy7gngCYBBgwa5ESNGBLD5n8rOzqa2y4ab7OxsbvrFCJq0Wc9dby8jp6Q1t47uEeqyaiXa3pdo6Eu09APUl6MVSLjnAu2qTGcCm6pps905txfYa2YfA/2AVUhAfnlCB1ZvLWRq9lq6NE/hgiyNIikitRfIMfccoKuZdTSzOOBi4K1D2rwJnGRmMWaWBAwFlge31OhmZtx1Ti9O6NyU219fwoLvdoa6JBGJYDWGu3OuHLgBmIMvsF92zi01s0lmNsnfZjnwLrAY+Ap4yjn3bf2VHZ1ivR4evWwgrdMTuO75BeTt1hk0IlI7AZ3n7pyb7Zzr5pzr7Jy7xz/vMefcY1XaTHHO9XTO9XbOPVhfBUe79KQ4pv1yECVllUx4VmfQiEjt6ArVMNSlRSoPXzqAlVvy+a+XdaNtETl6CvcwNbJ7C/77rON4d+kWHv7P6lCXIyIRRgOHhbFrTuzIii0FPPj+arq3TOXMPq1DXZKIRAjtuYcxM+Oe83ozsH06v3v5G5Zt0m36RCQwCvcwt/8+rGmJsVz73Hx2FJaEuiQRiQAK9wjQIjWBJ67IYnthCb/+59e6i5OI1EjhHiH6ZqZz/wV9+XL9Tv749tJQlyMiYU5fqEaQMQPasnxzPo9/vI6erdO4dGj7UJckImFKe+4R5veje3BKt+bc+da3zN+gIQpEpHoK9wjj9RgPXzyAzGOSmPTCAjZpiAIRqYbCPQKlJcXy5BVZFJdVct3zCyguqwh1SSISZhTuEapLi1T+dlF/luTt4fbXl1DTTVdEpHFRuEew03u25Hend+ONhXlM+3R9qMsRkTCicI9wN4zswuherfjz7OV8snpbqMsRkTChcI9wHo/xwLh+dGmRwg0vLmTjjn2hLklEwoDCPQqkxMfw5BW+++Ve+5zGgBcRhXvUOLZpMv+4dACrtxZwyyvf6AtWkUZO4R5FTuranNvO7MG/v93CIx+uCXU5IhJCCvcoc+1JnRjTvw1/eW8V7y/7IdTliEiIKNyjjJlx3wV96dWmCTfPXMSarYWhLklEQkDhHoUSYr08Pn4Q8TEeJj43nz1FZaEuSUQamMI9SrVNT2Tq5Vls3LmP38xYSIVusi3SqCjco9iQjhnceU4vsldu44G5K0Ndjog0II3nHuUuH9qeZZvymZq9lp6tm3B2vzahLklEGoD23KOcmfHHc3ox6NhjmPzqN3ybtyfUJYlIA1C4NwJxMR6mXp7FMUlxTHxuPtt1k22RqKdwbySap8bz5BWD2LG3lF+9sIDSct1kWySaKdwbkd5t0/i/sX3J2bCLO99aqiEKRKKYvlBtZM7t35blmwt47KO19GydyvhhHUJdkojUA+25N0KTR3Xn1B4t+OPby/hs7fZQlyMi9UDh3gh5PcZDF/enQ7Nkrv/n1xoDXiQKKdwbqdSEWJ66YhCVDiY8l0NBsYYoEIkmCvdGrEOzZB69bCBrt+3l5hmLNESBSBRRuDdyw7s0486ze/LBiq1MmaMhCkSihc6WEcYffywrt/jOoOnWMoXzB2aGuiQRqaOA9tzNbLSZrTSzNWZ22xHaDTazCjMbG7wSpb6ZGXed04thnZpy22tLWPDdzlCXJCJ1VGO4m5kXeAQ4E+gJXGJmPQ/T7n5gTrCLlPoX6/Uw9fKBtElPYOJzC8jdpTNoRCJZIHvuQ4A1zrl1zrlSYAZwbjXtbgReA7YGsT5pQOlJcUy7cjClFZVMeHa+zqARiWBW0yXo/kMso51zE/zT44GhzrkbqrRpC7wInApMA95xzr1azbomAhMBWrZsmTVjxoxaFV1YWEhKSkqtlg034diXpdsr+MuCYno383LzwHg8ZgEtF459qa1o6Uu09APUl/1Gjhy5wDk3qKZ2gXyhWt3/7EM/ER4EbnXOVdgRgsA59wTwBMCgQYPciBEjAtj8T2VnZ1PbZcNNOPZlBJCW+R1/mPUtnxS24M6zewW0XDj2pbaipS/R0g9QX45WIOGeC7SrMp0JbDqkzSBghj/YmwFnmVm5c25WUKqUBnf58ceybttenp63nk7NkjUGjUiECSTcc4CuZtYRyAMuBi6t2sA513H/czObju+wjII9wv2/nx/Hdzv2cudbS8nMSGJk9xahLklEAlTjF6rOuXLgBnxnwSwHXnbOLTWzSWY2qb4LlNDxeoyHLxnAca2bcMM/v2bZpvxQlyQiAQroPHfn3GznXDfnXGfn3D3+eY855x6rpu2V1X2ZKpEpOT6Gab8cTGpCLNc8m8OWPcWhLklEAqDhB6RGrdISmHblIPKLyrh6eg6FJeWhLklEaqBwl4D0apPGI5cNZOUPBdzw4teUV+g2fSLhTOEuARvRvQV/Orc32Su3cceb3+o2fSJhTAOHyVG5dGh7cnft49HstbRNT+SGU7uGuiQRqYbCXY7a5FHd2bKnmAfmrqJVWiJjszSKpEi4UbjLUTMz7rugL1sLSrjttcU0T43nlG7Na7WuWQvzmDJnJZt2F9EmPZHJo7ozZkDbIFcs9UXvX/jSMXeplbgY3yiS3Vqm8qsXFvDN97uPeh2zFuZx++tLyNtdhAPydhdx++tLmLUwL/gFS9Dp/QtvCneptdSEWKZfPZimKXFcPT2HLXuP7gyaKXNWUlRWcdC8orIK3REqQuj9C28Kd6mTFqkJPHvVEBzwwPxifsgP/CKnTbuLjmq+hBe9f+FN4S511ql5CtOvGkxhqeOKaV+xZ19g48C3SU88qvkSXvT+hTeFuwRF38x0bhqYwPrte7n62Rz2ldZ8FevkUd1JjPUeNC8x1svkUd3rq0wJIr1/4U3hLkHTs6mXBy/uz8KNu5j0wteUlh/5GPyYAW259/w+tE1PxIC26Ynce34fnW0RIfT+hTedCilBdVaf1tx3fl9+/9pifjtzEQ9fMgCv5/A3cBkzoK3CIILp/QtfCncJunGD25FfXMbd/1pOcryX+87vi+cIAS8iwadwl3ox4aROFBSX89AHq0mM9XLXOb040i0YRSS4FO5Sb27+WVf2lZbz5CfrSYjzctvoHgp4kQaicJd6Y2b891nHsa+0gsc/Wkd8jJffnd4t1GWJNAoKd6lXZsafzu1NaXklD3+wmliPceNpGklSpL4p3KXeeTy+gcYqKh1/eW8VMV4PvxrROdRliUQ1hbs0CK/HmHJhP8orHfe/uwKH49cjuoS6LJGopXCXBuP1GH8d1w8z+L93V+IcXD9SAS9SHxTu0qBivB7+Oq4/hm9UwfIKx02nddFZNCJBpnCXBuf1GH8Z1x+vx8Pf3l9FSXkFk0d1V8CLBJHCXULC6zGmjO1LXIyHR7PXUlxWyR2/OE4BLxIkCncJGY/H+PN5vYmP8fD0vPXsKy3nnvP6HHEsGhEJjMJdQsrMuPPsnqQmxPD3/6yhsKScv47rT1yMBiwVqQuFu4ScmfFfZ3QnJT6Ge/+9goLicqZePpCkOP3zFKkt7R5J2LjulM7cd34fPlm9jUuf/JJde0tDXZJIxFK4S1i5eEh7Hr0si2Wb87nw8c/J0/04RWpF4S5hZ3TvVjx71RB+yC/mvEfmsWxTfqhLEok4CncJS8M6N+WVScPwmDHu8c/5ZPW2UJckElEU7hK2erRqwhvXn0DmMYlc9UwOM3M2hrokkYihcJew1jotkVcmDWNY56bc+toS7n93BZWVLtRliYQ9hbuEvdSEWJ6+cjCXDGnP1Oy1/OqfC9hXWh7qskTCWkDhbmajzWylma0xs9uqef0yM1vs//nMzPoFv1RpzGK9Hv58Xm/u+EVP3lv2A2Onfs4mnUkjclg1hruZeYFHgDOBnsAlZtbzkGbrgVOcc32BPwFPBLtQETPjmhM7Mu3KwXy/cx9n//1Tvlq/M9RliYSlQPbchwBrnHPrnHOlwAzg3KoNnHOfOed2+Se/ADKDW6bIj0Z2b8Eb1w8nLTGWS5/8guc/34BzOg4vUpXV9J/CzMYCo51zE/zT44GhzrkbDtP+FqDH/vaHvDYRmAjQsmXLrBkzZtSq6MLCQlJSUmq1bLhRX2pvX5nj8cUlfLOtguFtYriiVxzx3uAMOhYt70u09APUl/1Gjhy5wDk3qMaGzrkj/gAXAk9VmR4P/P0wbUcCy4GmNa03KyvL1daHH35Y62XDjfpSNxUVle5v7610HW57x43620du/bbCoKw3Wt6XaOmHc+rLfsB8V0O+OucCOiyTC7SrMp0JbDq0kZn1BZ4CznXO7QhgvSJ15vEYN/+sG89cOZjNe4o5+++f8q/Fm0NdlkjIBRLuOUBXM+toZnHAxcBbVRuYWXvgdWC8c25V8MsUObIR3Vvwr5tOpHOLFK5/8WvumPUtxWUVoS5LJGRqDHfnXDlwAzAH3yGXl51zS81skplN8jf7H6Ap8KiZLTKz+fVWschhZB6TxMvXDePakzry/BffMeaReaz+oSDUZYmEREADZjvnZgOzD5n3WJXnE4CffIEq0tBmL9nM7CVbAFj1QwFnPfwJY/q3Zd6a7WzeU0yb9EQmj+rOmAFtg77tWQvzmDJnJZt2F9XrdgLxh1lLeOnL77m5dxnX3D6bS4a24+4xfUJSi4SG7oYgUWPWwjxuf30JRf7DMZUOXIXjlQW5B9rk7S7i9teXAAQ1eA/ddn1tJxB/mLWEF774cRyeCucOTCvgGw8NPyBRY8qclQfCdb/qTvQtKqtgypyV9b7t+thOIF768vujmi/RSeEuUeNohiMI9k1ADrftUAyRUHGYa1cON1+ik8Jdokab9MSA23oM3lyUF7QrWw+37aOpKVi8Vv2FXIebL9FJ4S5RY/Ko7iTGeg+aF+sxYg+5ajU+xkPmMUn8ZsYiJjw7n9xd++pl24mxXiaP6l7ndR+tS4a2O6r5Ep0U7hI1xgxoy73n96FteiIGtE1PZMqF/Zgytt9B8+6/oC8f3jKCP/z8OD5bu4PT//oxT3y8lrKKyqBu+97z+4TkbJm7x/Th8uPbH9hT95px+fHt9WVqI6OzZSSqjBnQttpArW7ehJM6Mbp3K+56ayl/nr2C1xbk8cdze3F8p6ZB3XYo3D2mD3eP6UN2djZrLxsR6nIkBLTnLo1a5jFJPPXLwTwxPovCknIufuILbnppITuLa78XLxIOtOcuApzRqxUndW3O1I/W8thHa3m3spLvvKu47pROJMXpv4lEHu25i/glxnn53end+OB3p9C/hZeHPljNyAeymZmzkQrdt1UijMJd5BDtMpL4df8EXpk0jNZpidz62hLOfOhj3lv2g24KIhFD4S5yGIM7ZPDGr0/g0csGUlbhuPa5+Yx59DM+Xb1dIS9hT+EucgRmxll9WvPeb0/m/gv6sC2/mMunfcm4xz9n3hqFvIQvhbtIAGK8Hi4a3J7/3DKC/z23F9/vLOKyp77kgqmf8Z8VOlwj4UfhLnIUEmK9XDGsA9mTR/Cnc3vxQ34JV0+fz1kPf8obC3PrdCGUSDAp3EVqISHWy3h/yD9wYT/KKir57cxvOPn/PuTxj9ayZ19ZqEuURk4n8IrUQazXw9isTM4f0JbsVVt54uN13PvvFTz4/mouyGrLFcM60K1laqjLlEZI4S4SBB6PcWqPlpzaoyVLN+1h+rwNvDw/lxe+2MjQjhlcfvyxnNGrJfEx3ppXJhIECneRIOvVJo0pF/bj9rOO4+X53/PCF99x40sLyUiO44KBbRk3qB1dtTcv9UzhLlJPMpLjmHRKZyae1IlP12znpa828sy8DTz5yXr6ZaYxNiuTn/dtQ0ZyXKhLlSikcBepZx6PcXK35pzcrTnbC0uYtTCPVxfkcsebS/nj28sY0b05Z/drw8+Oa0lyvP5LSnDoX5JIA2qWEs+Ekzox4aROLN+cz6yFeby5aBPvL99KQqyHU3u04MzerRnZowUpCnqpA/3rEQmR41o34bjWTbh1dA/mf7eLt7/ZxL+/3cLsJVuIi/FwUpdmnN6zJacd15LmqfGhLlcijMJdJMQ8HmNIxwyGdMzgrnN6seC7Xcxespn3lv3AByu2YraEvpnpnNajBSO6N6d3mzQ8Ht0PVY5M4S4SRrxVgv7Os3uyfHMBHyz3hfzf3l/FX99bRUZyHCd2acaJXZtxYpc9/m12AAANA0lEQVRmIbkJt4Q/hbtImDIzerZpQs82TbjxtK5sLyzh09XbyV65lU/XbOetbzYB0LFZMsd3yuD4Tk0Z0jGD1mkKe1G4i0SMZinxB+7T6pxj5Q8FfLp6O5+v3cE732zmpa++B6BdRiKDj81g4LHH4Aoqqah0eHUYp9FRuItEIDOjR6sm9GjVhAkndaK8opLlmwv4asNOvlq/g49Xb+P1hXkA3Jczh76Z6fRrl06/zDT6tkunTVoCZgr8aKZwF4kCMV4PfTLT6JOZxjUndsQ5x8ad+3jh3c8oSWnNwo27mfbpOsoqfEMTH5MUS++2afRs7Tvsc1zrJnRslkysV2MJRguFu0gUMjOObZrM8LaxjBjRG4CS8gpWbC5gce5ulm7K59tNe3hm3gZK/cMUx3k9dGqeTPdWqXRtkULXlql0aZFC+4wkhX4EUriLNBLxMV7foZl26QfmlVVUsm7bXpZvzmfFlgJWbMln/oZdvLlo04E2sV6jfUYSHZul0Kl5Msc2TaJj02TaN02idVqijueHKYW7SCMW6/XQvVUq3VsdPJBZYUk5a7YWsnZrIWu2FbJ+217Wb9/Lx6u3UVpeWWV5o216Iu0yksg8JonMYxJpm55I22MSaZOeSMvUeGK01x8SCncR+YmU+Bj6t0unf5W9fIDKSseW/GI2bN/Lhh37+H7XPjbu2Efurn3M3bSFHXtLD2rvMWiRmkCrtARaNfE9tmgST4vUBFqkxtOiSTzNUuLJSIrThVlBpnAXkYB5PEabdN9e+Qldfvr6vtJyNu0uIndXEZv3FLN5dxGb9hTzQ34xa7YVMm/tdgqKy3+ynNdjZCTH0TQ5jmYp8TRNiSMjOY6MpDgyUuI4JimO73ZU0HJzPsckxZGWGEtCrEdn/ByBwl1EgiYpLoYuLVLp0uLw49XvKy1na34JWwtK2FpQzPaCErYXlrKtoIQde0vZXljCxp372LW3lIKSgz8I7s/55MDzOK+HtKRYmiTEkJYYS5PEWJokxJKaEEPqgccYUuJ//En2//iee0mKi4na7wwCCnczGw08BHiBp5xz9x3yuvlfPwvYB1zpnPs6yLWKRK1ZC/OYMmclm3YX0SY9kcmjuvPK/I3MW7vzQJvhnTO4cFD7n7QDfjJv/nc7eenL77m5dxnX3D6bS4a24+4xfQLa7pgBbQ87P5Dl92+7wjm8Zj/ZdlJcDB2axbDo+9019uXOs3tycrfm7NxXyofzcujQrSe79pXx2drtZK/cxraCEgqKy/B6jLIKx4bte8kvLqeguOzAaZ81iY/xkBTnC/qkOC+JcV4SY70kxXlJiPU9T4jzkhDjJSHWQ0Lsj4/xMR7iY/yPsT8+j9v/4/3xebzXS2yMEev14FxgtdVFjeFuZl7gEeB0IBfIMbO3nHPLqjQ7E+jq/xkKTPU/ikgNZi3M4/bXl1BUVgFA3u4ibp656Cft5q3deVDY5+0uYvKr34CDskp3YN7vZi6isspyFc7xwhcbAQ4K2eq2e/vrS5j/3U5eW5D3k/nAQQFf3fJ12fbkV74B40Ao5+0u4o43l3Lv+X0YM6AtW5p6GdGnNbMW5vHB8q0Hli0uq+T7nUUH2u1XXFZBYUk5hcXlFBSXU1hSzt6ScvaWlrO3pIK9JeXsK61gb2k5+0p9z4tKKw487thbemC6pLyC4rJKisoqqKisezCf2TGWkSPrvJojCmTPfQiwxjm3DsDMZgDnAlXD/VzgOef7OPrCzNLNrLVzbnPQKxaJMlPmrDwQVEerur3TymraAbz05fcHBWx12y0qqziw133o/ClzVh4UntUtX5dtl1UTmoFut7p2vj1sL81SgjtccllFJcVlFZSUV1JS7nte6n/ue6ygpKyS0grf9IHH8krKKnw/nl0bg1pTdaymPw/MbCww2jk3wT89HhjqnLuhSpt3gPucc5/6pz8AbnXOzT9kXROBif7J7sDKWtbdDNhey2XDjfoSnhqsL3GtumTV17or9u3Bm5R2YLp0y5oFddluXZYPwrLNgO1HWrbqNsJcXf59Heuca15To0D23Kv7tuHQT4RA2uCcewJ4IoBtHrkgs/nOuUF1XU84UF/CU7T0xczml+/ZGvH9gOh5T6Bh+hLI1QW5QLsq05nAplq0ERGRBhJIuOcAXc2so5nFARcDbx3S5i3gCvM5Htij4+0iIqFT42EZ51y5md0AzMF3KuTTzrmlZjbJ//pjwGx8p0GuwXcq5FX1VzIQhEM7YUR9CU/R0pdo6QeoL0elxi9URUQk8mhEHxGRKKRwFxGJQmEf7maWYGZfmdk3ZrbUzP7on59hZu+Z2Wr/4zGhrjUQZuY1s4X+awMiuR8bzGyJmS0ys/n+eZHal3Qze9XMVpjZcjMbFol9MbPu/vdj/0++md0coX35rf//+7dm9pI/ByKuHwBm9ht/P5aa2c3+efXel7APd6AEONU51w/oD4z2n5FzG/CBc64r8IF/OhL8BlheZTpS+wEw0jnXv8r5upHal4eAd51zPYB++N6fiOuLc26l//3oD2ThO7nhDSKsL2bWFrgJGOSc643vRI6LibB+AJhZb+BafFf69wN+YWZdaYi+OOci5gdIAr7GN27NSqC1f35rYGWo6wug/kz/G3kq8I5/XsT1w1/rBqDZIfMiri9AE2A9/pMLIrkvh9R/BjAvEvsCtAW+BzLwndH3jr8/EdUPf50X4htscf/0HcDvG6IvkbDnvv9QxiJgK/Cec+5LoKXzn0vvf2wRyhoD9CC+N7bqEByR2A/wXYE818wW+IeVgMjsSydgG/CM/3DZU2aWTGT2paqLgZf8zyOqL865POABYCOwGd91M3OJsH74fQucbGZNzSwJ3ynj7WiAvkREuDvnKpzvT81MYIj/T52IYma/ALY65yJl7IuaDHfODcQ3Iuj1ZnZyqAuqpRhgIDDVOTcA2EsE/Ll/JP6LDc8BXgl1LbXhP/58LtARaAMkm9nloa2qdpxzy4H7gfeAd4FvgJ/eraQeRES47+ec2w1kA6OBH8ysNYD/cWsISwvEcOAcM9sAzABONbMXiLx+AOCc2+R/3IrvuO4QIrMvuUCu/69BgFfxhX0k9mW/M4GvnXM/+KcjrS8/A9Y757Y558qA14ETiLx+AOCcm+acG+icOxnYCaymAfoS9uFuZs3NLN3/PBHfG78C35AHv/Q3+yXwZmgqDIxz7nbnXKZzrgO+P5n/45y7nAjrB4CZJZtZ6v7n+I6HfksE9sU5twX43sy6+2edhm8464jrSxWX8OMhGYi8vmwEjjezJDMzfO/JciKvHwCYWQv/Y3vgfHzvTb33JeyvUDWzvsCz+L4x9wAvO+f+18yaAi8D7fH9Y7jQObfz8GsKH2Y2ArjFOfeLSOyHmXXCt7cOvsMaLzrn7onEvgCYWX/gKSAOWIdv+AwPkdmXJHxfRnZyzu3xz4u498V/yvNF+A5hLAQmAClEWD8AzOwToClQBvzOOfdBQ7wnYR/uIiJy9ML+sIyIiBw9hbuISBRSuIuIRCGFu4hIFFK4i4hEoUBukC3SoPyniX3gn2wFVOAbIgBgiHOuNCSFHYGZXQ3M9p83LxJyOhVSwpqZ3QUUOuceCINavM65isO89ilwg3Nu0VGsL8Y51yCXokvjo8MyElHM7JfmG99/kZk9amYeM4sxs91mNsXMvjazOWY21Mw+MrN1ZnaWf9kJZvaG//WVZvaHANd7t5l9hW9coz+aWY5/fO7HzOcifMNRz/QvH2dmuVWurD7ezN73P7/bzB43s/fwDVYWY2Z/9W97sZlNaPjfqkQjhbtEDP+AcecBJ/gHkovBN5QDQBow1z+YWSlwF77L1i8E/rfKaob4lxkIXGpm/QNY79fOuSHOuc+Bh5xzg4E+/tdGO+dmAouAi5xvPPWaDhsNAM52zo0HJuIbUG4IMBjfIGzta/P7EalKx9wlkvwMXwDO9w05QiK+S+0Bipxz7/mfL8E3TGy5mS0BOlRZxxzn3C4AM5sFnIjv/8Hh1lvKj0MtAJxmZpOBBKAZsAD491H2403nXLH/+RnAcWZW9cOkK75L0kVqTeEukcSAp51zdxw00ywGXwjvV4nvDl77n1f9d37ol0yuhvUWOf8XU/5xW/4BDHTO5ZnZ3fhCvjrl/PiX8aFt9h7Sp1875z5AJIh0WEYiyfvAODNrBr6zampxCOMM890zNQnfmOHzjmK9ifg+LLb7R8W8oMprBUBqlekN+G51xyHtDjUH+LX/g2T/fVATj7JPIj+hPXeJGM65Jf7RAt83Mw++UfYmAZuOYjWfAi8CnYHn95/dEsh6nXM7zOxZfMMbfwd8WeXlZ4CnzKwI33H9u4AnzWwL8NUR6nkc38iAi/yHhLbi+9ARqROdCimNhv9MlN7OuZtDXYtIfdNhGRGRKKQ9dxGRKKQ9dxGRKKRwFxGJQgp3EZEopHAXEYlCCncRkSj0/wHRUJwHFwSFegAAAABJRU5ErkJggg==\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "data[\"Success\"]=data.Count-data.Malfunction\n", "data[\"Intercept\"]=1\n", "logmodel=sm.GLM(data['Frequency'], data[['Intercept','Temperature']], family=sm.families.Binomial(sm.families.links.logit)).fit()\n", "\n", "%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": {}, "source": [ "

Now on the plot above we reuse the same regression technique, however, this time we include all of dataset. Now we can see clearly that the temperature of the experiment indeed has an influence on the probability of the failure. From the corresponding plot we can see the probability is arround 80 percent which is pretty high.

" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [] } ], "metadata": { "celltoolbar": "Hide code", "kernelspec": { "display_name": "Python 3", "language": "python", "name": "python3" }, "language_info": { "codemirror_mode": { "name": "ipython", "version": 3 }, "file_extension": ".py", "mimetype": "text/x-python", "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", "version": "3.6.4" } }, "nbformat": 4, "nbformat_minor": 2 }