Skip to content

M
mooc-rr
Commit db834875 authored by Samuel Loury's avatar Samuel Loury
Browse files
Options

Add confidence intervals and analyse pressure

parent a00792ce
Hide whitespace changes
Inline Side-by-side
Showing with 669 additions and 181 deletions
module2/exo5/exo5_python_fr.html 0 → 100644
View file @ db834875
<?xml version="1.0" encoding="utf-8"?>
<!DOCTYPE html PUBLIC "-//W3C//DTD XHTML 1.0 Strict//EN"
"http://www.w3.org/TR/xhtml1/DTD/xhtml1-strict.dtd">
<html xmlns="http://www.w3.org/1999/xhtml" lang="fr" xml:lang="fr">
<head>
<!-- 2020-07-20 Mon 11:26 -->
<meta http-equiv="Content-Type" content="text/html;charset=utf-8" />
<meta name="viewport" content="width=device-width, initial-scale=1" />
<title>Analyse du risque de défaillance des joints toriques de la navette Challenger</title>
<meta name="generator" content="Org mode" />
<meta name="author" content="Arnaud Legrand" />
<style type="text/css">
<!--/*--><![CDATA[/*><!--*/
.title { text-align: center;
margin-bottom: .2em; }
.subtitle { text-align: center;
font-size: medium;
font-weight: bold;
margin-top:0; }
.todo { font-family: monospace; color: red; }
.done { font-family: monospace; color: green; }
.priority { font-family: monospace; color: orange; }
.tag { background-color: #eee; font-family: monospace;
padding: 2px; font-size: 80%; font-weight: normal; }
.timestamp { color: #bebebe; }
.timestamp-kwd { color: #5f9ea0; }
.org-right { margin-left: auto; margin-right: 0px; text-align: right; }
.org-left { margin-left: 0px; margin-right: auto; text-align: left; }
.org-center { margin-left: auto; margin-right: auto; text-align: center; }
.underline { text-decoration: underline; }
#postamble p, #preamble p { font-size: 90%; margin: .2em; }
p.verse { margin-left: 3%; }
pre {
border: 1px solid #ccc;
box-shadow: 3px 3px 3px #eee;
padding: 8pt;
font-family: monospace;
overflow: auto;
margin: 1.2em;
}
pre.src {
position: relative;
overflow: visible;
padding-top: 1.2em;
}
pre.src:before {
display: none;
position: absolute;
background-color: white;
top: -10px;
right: 10px;
padding: 3px;
border: 1px solid black;
}
pre.src:hover:before { display: inline;}
/* Languages per Org manual */
pre.src-asymptote:before { content: 'Asymptote'; }
pre.src-awk:before { content: 'Awk'; }
pre.src-C:before { content: 'C'; }
/* pre.src-C++ doesn't work in CSS */
pre.src-clojure:before { content: 'Clojure'; }
pre.src-css:before { content: 'CSS'; }
pre.src-D:before { content: 'D'; }
pre.src-ditaa:before { content: 'ditaa'; }
pre.src-dot:before { content: 'Graphviz'; }
pre.src-calc:before { content: 'Emacs Calc'; }
pre.src-emacs-lisp:before { content: 'Emacs Lisp'; }
pre.src-fortran:before { content: 'Fortran'; }
pre.src-gnuplot:before { content: 'gnuplot'; }
pre.src-haskell:before { content: 'Haskell'; }
pre.src-hledger:before { content: 'hledger'; }
pre.src-java:before { content: 'Java'; }
pre.src-js:before { content: 'Javascript'; }
pre.src-latex:before { content: 'LaTeX'; }
pre.src-ledger:before { content: 'Ledger'; }
pre.src-lisp:before { content: 'Lisp'; }
pre.src-lilypond:before { content: 'Lilypond'; }
pre.src-lua:before { content: 'Lua'; }
pre.src-matlab:before { content: 'MATLAB'; }
pre.src-mscgen:before { content: 'Mscgen'; }
pre.src-ocaml:before { content: 'Objective Caml'; }
pre.src-octave:before { content: 'Octave'; }
pre.src-org:before { content: 'Org mode'; }
pre.src-oz:before { content: 'OZ'; }
pre.src-plantuml:before { content: 'Plantuml'; }
pre.src-processing:before { content: 'Processing.js'; }
pre.src-python:before { content: 'Python'; }
pre.src-R:before { content: 'R'; }
pre.src-ruby:before { content: 'Ruby'; }
pre.src-sass:before { content: 'Sass'; }
pre.src-scheme:before { content: 'Scheme'; }
pre.src-screen:before { content: 'Gnu Screen'; }
pre.src-sed:before { content: 'Sed'; }
pre.src-sh:before { content: 'shell'; }
pre.src-sql:before { content: 'SQL'; }
pre.src-sqlite:before { content: 'SQLite'; }
/* additional languages in org.el's org-babel-load-languages alist */
pre.src-forth:before { content: 'Forth'; }
pre.src-io:before { content: 'IO'; }
pre.src-J:before { content: 'J'; }
pre.src-makefile:before { content: 'Makefile'; }
pre.src-maxima:before { content: 'Maxima'; }
pre.src-perl:before { content: 'Perl'; }
pre.src-picolisp:before { content: 'Pico Lisp'; }
pre.src-scala:before { content: 'Scala'; }
pre.src-shell:before { content: 'Shell Script'; }
pre.src-ebnf2ps:before { content: 'ebfn2ps'; }
/* additional language identifiers per "defun org-babel-execute"
in ob-*.el */
pre.src-cpp:before { content: 'C++'; }
pre.src-abc:before { content: 'ABC'; }
pre.src-coq:before { content: 'Coq'; }
pre.src-groovy:before { content: 'Groovy'; }
/* additional language identifiers from org-babel-shell-names in
ob-shell.el: ob-shell is the only babel language using a lambda to put
the execution function name together. */
pre.src-bash:before { content: 'bash'; }
pre.src-csh:before { content: 'csh'; }
pre.src-ash:before { content: 'ash'; }
pre.src-dash:before { content: 'dash'; }
pre.src-ksh:before { content: 'ksh'; }
pre.src-mksh:before { content: 'mksh'; }
pre.src-posh:before { content: 'posh'; }
/* Additional Emacs modes also supported by the LaTeX listings package */
pre.src-ada:before { content: 'Ada'; }
pre.src-asm:before { content: 'Assembler'; }
pre.src-caml:before { content: 'Caml'; }
pre.src-delphi:before { content: 'Delphi'; }
pre.src-html:before { content: 'HTML'; }
pre.src-idl:before { content: 'IDL'; }
pre.src-mercury:before { content: 'Mercury'; }
pre.src-metapost:before { content: 'MetaPost'; }
pre.src-modula-2:before { content: 'Modula-2'; }
pre.src-pascal:before { content: 'Pascal'; }
pre.src-ps:before { content: 'PostScript'; }
pre.src-prolog:before { content: 'Prolog'; }
pre.src-simula:before { content: 'Simula'; }
pre.src-tcl:before { content: 'tcl'; }
pre.src-tex:before { content: 'TeX'; }
pre.src-plain-tex:before { content: 'Plain TeX'; }
pre.src-verilog:before { content: 'Verilog'; }
pre.src-vhdl:before { content: 'VHDL'; }
pre.src-xml:before { content: 'XML'; }
pre.src-nxml:before { content: 'XML'; }
/* add a generic configuration mode; LaTeX export needs an additional
(add-to-list 'org-latex-listings-langs '(conf " ")) in .emacs */
pre.src-conf:before { content: 'Configuration File'; }
table { border-collapse:collapse; }
caption.t-above { caption-side: top; }
caption.t-bottom { caption-side: bottom; }
td, th { vertical-align:top; }
th.org-right { text-align: center; }
th.org-left { text-align: center; }
th.org-center { text-align: center; }
td.org-right { text-align: right; }
td.org-left { text-align: left; }
td.org-center { text-align: center; }
dt { font-weight: bold; }
.footpara { display: inline; }
.footdef { margin-bottom: 1em; }
.figure { padding: 1em; }
.figure p { text-align: center; }
.equation-container {
display: table;
text-align: center;
width: 100%;
}
.equation {
vertical-align: middle;
}
.equation-label {
display: table-cell;
text-align: right;
vertical-align: middle;
}
.inlinetask {
padding: 10px;
border: 2px solid gray;
margin: 10px;
background: #ffffcc;
}
#org-div-home-and-up
{ text-align: right; font-size: 70%; white-space: nowrap; }
textarea { overflow-x: auto; }
.linenr { font-size: smaller }
.code-highlighted { background-color: #ffff00; }
.org-info-js_info-navigation { border-style: none; }
#org-info-js_console-label
{ font-size: 10px; font-weight: bold; white-space: nowrap; }
.org-info-js_search-highlight
{ background-color: #ffff00; color: #000000; font-weight: bold; }
.org-svg { width: 90%; }
/*]]>*/-->
</style>
<link rel="stylesheet" type="text/css" href="http://www.pirilampo.org/styles/readtheorg/css/htmlize.css"/>
<link rel="stylesheet" type="text/css" href="http://www.pirilampo.org/styles/readtheorg/css/readtheorg.css"/>
<script src="https://ajax.googleapis.com/ajax/libs/jquery/2.1.3/jquery.min.js"></script>
<script src="https://maxcdn.bootstrapcdn.com/bootstrap/3.3.4/js/bootstrap.min.js"></script>
<script type="text/javascript" src="http://www.pirilampo.org/styles/lib/js/jquery.stickytableheaders.js"></script>
<script type="text/javascript" src="http://www.pirilampo.org/styles/readtheorg/js/readtheorg.js"></script>
<script type="text/javascript">
/*
@licstart The following is the entire license notice for the
JavaScript code in this tag.
Copyright (C) 2012-2020 Free Software Foundation, Inc.
The JavaScript code in this tag is free software: you can
redistribute it and/or modify it under the terms of the GNU
General Public License (GNU GPL) as published by the Free Software
Foundation, either version 3 of the License, or (at your option)
any later version. The code is distributed WITHOUT ANY WARRANTY;
without even the implied warranty of MERCHANTABILITY or FITNESS
FOR A PARTICULAR PURPOSE. See the GNU GPL for more details.
As additional permission under GNU GPL version 3 section 7, you
may distribute non-source (e.g., minimized or compacted) forms of
that code without the copy of the GNU GPL normally required by
section 4, provided you include this license notice and a URL
through which recipients can access the Corresponding Source.
@licend The above is the entire license notice
for the JavaScript code in this tag.
*/
<!--/*--><![CDATA[/*><!--*/
function CodeHighlightOn(elem, id)
{
var target = document.getElementById(id);
if(null != target) {
elem.cacheClassElem = elem.className;
elem.cacheClassTarget = target.className;
target.className = "code-highlighted";
elem.className = "code-highlighted";
}
}
function CodeHighlightOff(elem, id)
{
var target = document.getElementById(id);
if(elem.cacheClassElem)
elem.className = elem.cacheClassElem;
if(elem.cacheClassTarget)
target.className = elem.cacheClassTarget;
}
/*]]>*///-->
</script>
<script type="text/x-mathjax-config">
MathJax.Hub.Config({
displayAlign: "center",
displayIndent: "0em",
"HTML-CSS": { scale: 100,
linebreaks: { automatic: "false" },
webFont: "TeX"
},
SVG: {scale: 100,
linebreaks: { automatic: "false" },
font: "TeX"},
NativeMML: {scale: 100},
TeX: { equationNumbers: {autoNumber: "AMS"},
MultLineWidth: "85%",
TagSide: "right",
TagIndent: ".8em"
}
});
</script>
<script type="text/javascript"
src="https://cdnjs.cloudflare.com/ajax/libs/mathjax/2.7.0/MathJax.js?config=TeX-AMS_HTML"></script>
</head>
<body>
<div id="content">
<h1 class="title">Analyse du risque de défaillance des joints toriques de la navette Challenger</h1>
<div id="table-of-contents">
<h2>Table des matières</h2>
<div id="text-table-of-contents">
<ul>
<li><a href="#org496f632">1. Chargement des données</a></li>
<li><a href="#org28650b8">2. Inspection graphique des données</a></li>
<li><a href="#orgf095e70">3. Estimation de l'influence de la température</a></li>
<li><a href="#orgb01bc65">4. Estimation de la probabilité de dysfonctionnant des joints toriques</a></li>
</ul>
</div>
</div>
<p>
Le 27 Janvier 1986, veille du décollage de la navette <i>Challenger</i>, eu
lieu une télé-conférence de trois heures entre les ingénieurs de la
Morton Thiokol (constructeur d'un des moteurs) et de la NASA. La
discussion portait principalement sur les conséquences de la
température prévue au moment du décollage de 31°F (juste en dessous de
0°C) sur le succès du vol et en particulier sur la performance des
joints toriques utilisés dans les moteurs. En effet, aucun test
n'avait été effectué à cette température.
</p>
<p>
L'étude qui suit reprend donc une partie des analyses effectuées cette
nuit là et dont l'objectif était d'évaluer l'influence potentielle de
la température et de la pression à laquelle sont soumis les joints
toriques sur leur probabilité de dysfonctionnement. Pour cela, nous
disposons des résultats des expériences réalisées par les ingénieurs
de la NASA durant les 6 années précédant le lancement de la navette
Challenger.
</p>
<div id="outline-container-org496f632" class="outline-2">
<h2 id="org496f632"><span class="section-number-2">1</span> Chargement des données</h2>
<div class="outline-text-2" id="text-1">
<p>
Nous commençons donc par charger ces données:
</p>
<div class="org-src-container">
<pre class="src src-python"><span style="color: #F0DFAF; font-weight: bold;">import</span> numpy <span style="color: #F0DFAF; font-weight: bold;">as</span> np
<span style="color: #F0DFAF; font-weight: bold;">import</span> pandas <span style="color: #F0DFAF; font-weight: bold;">as</span> pd
<span style="color: #DFAF8F;">data</span> = pd.read_csv(<span style="color: #CC9393;">"shuttle.csv"</span>)
data
</pre>
</div>
<pre class="example">
Date Count Temperature Pressure Malfunction
0 4/12/81 6 66 50 0
1 11/12/81 6 70 50 1
2 3/22/82 6 69 50 0
3 11/11/82 6 68 50 0
4 4/04/83 6 67 50 0
5 6/18/82 6 72 50 0
6 8/30/83 6 73 100 0
7 11/28/83 6 70 100 0
8 2/03/84 6 57 200 1
9 4/06/84 6 63 200 1
10 8/30/84 6 70 200 1
11 10/05/84 6 78 200 0
12 11/08/84 6 67 200 0
13 1/24/85 6 53 200 2
14 4/12/85 6 67 200 0
15 4/29/85 6 75 200 0
16 6/17/85 6 70 200 0
17 7/29/85 6 81 200 0
18 8/27/85 6 76 200 0
19 10/03/85 6 79 200 0
20 10/30/85 6 75 200 2
21 11/26/85 6 76 200 0
22 1/12/86 6 58 200 1
</pre>
<p>
Le jeu de données nous indique la date de l'essai, le nombre de joints
toriques mesurés (il y en a 6 sur le lanceur principal), la
température (en Fahrenheit) et la pression (en psi), et enfin le
nombre de dysfonctionnements relevés.
</p>
</div>
</div>
<div id="outline-container-org28650b8" class="outline-2">
<h2 id="org28650b8"><span class="section-number-2">2</span> Inspection graphique des données</h2>
<div class="outline-text-2" id="text-2">
<p>
Comment la fréquence d'échecs varie-t-elle avec la température ?
</p>
<div class="org-src-container">
<pre class="src src-python"><span style="color: #F0DFAF; font-weight: bold;">import</span> matplotlib.pyplot <span style="color: #F0DFAF; font-weight: bold;">as</span> plt
plt.clf()
<span style="color: #DFAF8F;">data</span>[<span style="color: #CC9393;">"Frequency"</span>]=data.Malfunction/data.Count
data.plot(x=<span style="color: #CC9393;">"Temperature"</span>,y=<span style="color: #CC9393;">"Frequency"</span>,kind=<span style="color: #CC9393;">"scatter"</span>,ylim=[0,1])
plt.grid(<span style="color: #BFEBBF;">True</span>)
plt.savefig(matplot_lib_filename)
<span style="color: #F0DFAF; font-weight: bold;">print</span>(matplot_lib_filename)
</pre>
</div>
<div class="figure">
<p><img src="freq_temp_python.png" alt="freq_temp_python.png" />
</p>
</div>
<p>
Comment la fréquence d'échecs varie-t-elle avec la pression ?
</p>
<div class="org-src-container">
<pre class="src src-python"><span style="color: #F0DFAF; font-weight: bold;">import</span> matplotlib.pyplot <span style="color: #F0DFAF; font-weight: bold;">as</span> plt
plt.clf()
<span style="color: #DFAF8F;">data</span>[<span style="color: #CC9393;">"Frequency"</span>]=data.Malfunction/data.Count
data.plot(x=<span style="color: #CC9393;">"Pressure"</span>,y=<span style="color: #CC9393;">"Frequency"</span>,kind=<span style="color: #CC9393;">"scatter"</span>,ylim=[0,1])
plt.grid(<span style="color: #BFEBBF;">True</span>)
plt.savefig(matplot_lib_filename)
<span style="color: #F0DFAF; font-weight: bold;">print</span>(matplot_lib_filename)
</pre>
</div>
<div class="figure">
<p><img src="freq_pressure_python.png" alt="freq_pressure_python.png" />
</p>
</div>
<p>
Avec si peu de points de pression, on ne peut pas étudier grand chose. De
plus, visuellement, aucune relation n'a l'air de sauter aux yeux.
</p>
<p>
À première vue, ce n'est pas flagrant mais bon, essayons quand même
d'estimer l'impact de la température \(t\) sur la probabilité de
dysfonctionnements d'un joint.
</p>
</div>
</div>
<div id="outline-container-orgf095e70" class="outline-2">
<h2 id="orgf095e70"><span class="section-number-2">3</span> Estimation de l'influence de la température</h2>
<div class="outline-text-2" id="text-3">
<p>
Supposons que chacun des 6 joints toriques est endommagé avec la même
probabilité et indépendamment des autres et que cette probabilité ne
dépend que de la température. Si on note \(p(t)\) cette probabilité, le
nombre de joints \(D\) dysfonctionnant lorsque l'on effectue le vol à
température \(t\) suit une loi binomiale de paramètre \(n=6\) et
\(p=p(t)\). Pour relier \(p(t)\) à \(t\), on va donc effectuer une
régression logistique.
</p>
<div class="org-src-container">
<pre class="src src-python"><span style="color: #F0DFAF; font-weight: bold;">import</span> statsmodels.api <span style="color: #F0DFAF; font-weight: bold;">as</span> sm
<span style="color: #DFAF8F;">data</span>[<span style="color: #CC9393;">"Intercept"</span>]=1
<span style="color: #DFAF8F;">logmodel</span>=sm.GLM(data[<span style="color: #CC9393;">'Frequency'</span>], data[ [<span style="color: #CC9393;">'Intercept'</span>, <span style="color: #CC9393;">'Temperature'</span>] ], family=sm.families.Binomial(sm.families.links.logit)).fit()
logmodel.summary()
</pre>
</div>
<pre class="example">
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: Mon, 20 Jul 2020 Deviance: 3.0144
Time: 11:26:47 Pearson chi2: 5.00
No. Iterations: 6
Covariance Type: nonrobust
===============================================================================
coef std err z P&gt;|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
===============================================================================
</pre>
<p>
L'estimateur le plus probable du paramètre de température est 0.0014
et l'erreur standard de cet estimateur est de 0.122, autrement dit on
ne peut pas distinguer d'impact particulier et il faut prendre nos
estimations avec des pincettes.
</p>
</div>
</div>
<div id="outline-container-orgb01bc65" class="outline-2">
<h2 id="orgb01bc65"><span class="section-number-2">4</span> Estimation de la probabilité de dysfonctionnant des joints toriques</h2>
<div class="outline-text-2" id="text-4">
<p>
La température prévue le jour du décollage est de 31°F. Essayons
d'estimer la probabilité de dysfonctionnement des joints toriques à
cette température à partir du modèle que nous venons de construire:
</p>
<div class="org-src-container">
<pre class="src src-python"><span style="color: #F0DFAF; font-weight: bold;">import</span> matplotlib.pyplot <span style="color: #F0DFAF; font-weight: bold;">as</span> plt
<span style="color: #DFAF8F;">data_pred</span> = pd.DataFrame({<span style="color: #CC9393;">'Temperature'</span>: np.linspace(start=30, stop=90, num=121), <span style="color: #CC9393;">'Intercept'</span>: 1})
<span style="color: #DFAF8F;">prediction</span> = logmodel.get_prediction(data_pred[ [<span style="color: #CC9393;">'Intercept'</span>,<span style="color: #CC9393;">'Temperature'</span>] ]).summary_frame()
<span style="color: #DFAF8F;">data_pred</span>[<span style="color: #CC9393;">'Frequency'</span>] = prediction[<span style="color: #CC9393;">"mean"</span>]
data_pred.plot(x=<span style="color: #CC9393;">"Temperature"</span>,y=<span style="color: #CC9393;">"Frequency"</span>,kind=<span style="color: #CC9393;">"line"</span>,ylim=[0,1])
plt.scatter(x=data[<span style="color: #CC9393;">"Temperature"</span>], y=data[<span style="color: #CC9393;">"Frequency"</span>])
plt.fill_between(
<span style="background-color: #4d4d4d;"> </span>data_pred[<span style="color: #CC9393;">"Temperature"</span>],
<span style="background-color: #4d4d4d;"> </span>prediction[<span style="color: #CC9393;">"mean_ci_lower"</span>],
<span style="background-color: #4d4d4d;"> </span>prediction[<span style="color: #CC9393;">"mean_ci_upper"</span>],
<span style="background-color: #4d4d4d;"> </span>alpha=0.2,
)
plt.grid(<span style="color: #BFEBBF;">True</span>)
plt.savefig(matplot_lib_filename)
<span style="color: #F0DFAF; font-weight: bold;">print</span>(matplot_lib_filename)
</pre>
</div>
<div class="figure">
<p><img src="proba_estimate_python.png" alt="proba_estimate_python.png" />
</p>
</div>
<p>
Hmm, l'intervalle de confiance semble très grand. Peut-on déduire quoi que ce
soit de ces données ?
</p>
</div>
</div>
</div>
<div id="postamble" class="status">
<p class="author">Auteur: Arnaud Legrand</p>
<p class="date">Created: 2020-07-20 Mon 11:26</p>
<p class="validation"><a href="http://validator.w3.org/check?uri=referer">Validate</a></p>
</div>
</body>
</html>
module2/exo5/exo5_python_fr.org
View file @ db834875
...@@ -32,191 +32,161 @@ de la NASA durant les 6 années précédant le lancement de la navette ...@@ -32,191 +32,161 @@ de la NASA durant les 6 années précédant le lancement de la navette
Challenger. Challenger.
* Chargement des données * Chargement des données
Nous commençons donc par charger ces données: Nous commençons donc par charger ces données:
#+begin_src python :results value :session *python* :exports both #+begin_src python :results value :session *python* :exports both
import numpy as np import numpy as np
import pandas as pd import pandas as pd
data = pd.read_csv("shuttle.csv") data = pd.read_csv("shuttle.csv")
data data
#+end_src #+end_src
#+RESULTS: #+RESULTS:
#+begin_example #+begin_example
Date Count Temperature Pressure Malfunction Date Count Temperature Pressure Malfunction
0 4/12/81 6 66 50 0 0 4/12/81 6 66 50 0
1 11/12/81 6 70 50 1 1 11/12/81 6 70 50 1
2 3/22/82 6 69 50 0 2 3/22/82 6 69 50 0
3 11/11/82 6 68 50 0 3 11/11/82 6 68 50 0
4 4/04/83 6 67 50 0 4 4/04/83 6 67 50 0
5 6/18/82 6 72 50 0 5 6/18/82 6 72 50 0
6 8/30/83 6 73 100 0 6 8/30/83 6 73 100 0
7 11/28/83 6 70 100 0 7 11/28/83 6 70 100 0
8 2/03/84 6 57 200 1 8 2/03/84 6 57 200 1
9 4/06/84 6 63 200 1 9 4/06/84 6 63 200 1
10 8/30/84 6 70 200 1 10 8/30/84 6 70 200 1
11 10/05/84 6 78 200 0 11 10/05/84 6 78 200 0
12 11/08/84 6 67 200 0 12 11/08/84 6 67 200 0
13 1/24/85 6 53 200 2 13 1/24/85 6 53 200 2
14 4/12/85 6 67 200 0 14 4/12/85 6 67 200 0
15 4/29/85 6 75 200 0 15 4/29/85 6 75 200 0
16 6/17/85 6 70 200 0 16 6/17/85 6 70 200 0
17 7/29/85 6 81 200 0 17 7/29/85 6 81 200 0
18 8/27/85 6 76 200 0 18 8/27/85 6 76 200 0
19 10/03/85 6 79 200 0 19 10/03/85 6 79 200 0
20 10/30/85 6 75 200 2 20 10/30/85 6 75 200 2
21 11/26/85 6 76 200 0 21 11/26/85 6 76 200 0
22 1/12/86 6 58 200 1 22 1/12/86 6 58 200 1
#+end_example #+end_example
Le jeu de données nous indique la date de l'essai, le nombre de joints Le jeu de données nous indique la date de l'essai, le nombre de joints
toriques mesurés (il y en a 6 sur le lançeur principal), la toriques mesurés (il y en a 6 sur le lanceur principal), la
température (en Fahrenheit) et la pression (en psi), et enfin le température (en Fahrenheit) et la pression (en psi), et enfin le
nombre de dysfonctionnements relevés. nombre de dysfonctionnements relevés.
* Inspection graphique des données * Inspection graphique des données
Les vols où aucun incident n'est relevé n'apportant aucune information
sur l'influence de la température ou de la pression sur les Comment la fréquence d'échecs varie-t-elle avec la température ?
dysfonctionnements, nous nous concentrons sur les expériences où au #+begin_src python :results output file :var matplot_lib_filename="freq_temp_python.png" :exports both :session *python*
moins un joint a été défectueux. import matplotlib.pyplot as plt
#+begin_src python :results value :session *python* :exports both plt.clf()
data = data[data.Malfunction>0] data["Frequency"]=data.Malfunction/data.Count
data data.plot(x="Temperature",y="Frequency",kind="scatter",ylim=[0,1])
#+end_src plt.grid(True)
#+RESULTS: plt.savefig(matplot_lib_filename)
: Date Count Temperature Pressure Malfunction print(matplot_lib_filename)
: 1 11/12/81 6 70 50 1 #+end_src
: 8 2/03/84 6 57 200 1
: 9 4/06/84 6 63 200 1 #+RESULTS:
: 10 8/30/84 6 70 200 1 [[file:freq_temp_python.png]]
: 13 1/24/85 6 53 200 2
: 20 10/30/85 6 75 200 2 Comment la fréquence d'échecs varie-t-elle avec la pression ?
: 22 1/12/86 6 58 200 1 #+begin_src python :results output file :var matplot_lib_filename="freq_pressure_python.png" :exports both :session *python*
import matplotlib.pyplot as plt
Très bien, nous avons une variabilité de température importante mais
la pression est quasiment toujours égale à 200, ce qui devrait plt.clf()
simplifier l'analyse. data["Frequency"]=data.Malfunction/data.Count
data.plot(x="Pressure",y="Frequency",kind="scatter",ylim=[0,1])
Comment la fréquence d'échecs varie-t-elle avec la température ? plt.grid(True)
#+begin_src python :results output file :var matplot_lib_filename="freq_temp_python.png" :exports both :session *python*
import matplotlib.pyplot as plt plt.savefig(matplot_lib_filename)
print(matplot_lib_filename)
plt.clf() #+end_src
data["Frequency"]=data.Malfunction/data.Count
data.plot(x="Temperature",y="Frequency",kind="scatter",ylim=[0,1]) #+RESULTS:
plt.grid(True) [[file:freq_pressure_python.png]]
plt.savefig(matplot_lib_filename) Avec si peu de points de pression, on ne peut pas étudier grand chose. De
print(matplot_lib_filename) plus, visuellement, aucune relation n'a l'air de sauter aux yeux.
#+end_src
À première vue, ce n'est pas flagrant mais bon, essayons quand même
#+RESULTS: d'estimer l'impact de la température $t$ sur la probabilité de
[[file:freq_temp_python.png]] dysfonctionnements d'un joint.
À première vue, ce n'est pas flagrant mais bon, essayons quand même
d'estimer l'impact de la température $t$ sur la probabilité de
dysfonctionnements d'un joint.
* Estimation de l'influence de la température * Estimation de l'influence de la température
Supposons que chacun des 6 joints toriques est endommagé avec la même Supposons que chacun des 6 joints toriques est endommagé avec la même
probabilité et indépendamment des autres et que cette probabilité ne probabilité et indépendamment des autres et que cette probabilité ne
dépend que de la température. Si on note $p(t)$ cette probabilité, le dépend que de la température. Si on note $p(t)$ cette probabilité, le
nombre de joints $D$ dysfonctionnant lorsque l'on effectue le vol à nombre de joints $D$ dysfonctionnant lorsque l'on effectue le vol à
température $t$ suit une loi binomiale de paramètre $n=6$ et température $t$ suit une loi binomiale de paramètre $n=6$ et
$p=p(t)$. Pour relier $p(t)$ à $t$, on va donc effectuer une $p=p(t)$. Pour relier $p(t)$ à $t$, on va donc effectuer une
régression logistique. régression logistique.
#+begin_src python :results value :session *python* :exports both #+begin_src python :results value :session *python* :exports both
import statsmodels.api as sm import statsmodels.api as sm
data["Success"]=data.Count-data.Malfunction data["Intercept"]=1
data["Intercept"]=1
logmodel=sm.GLM(data['Frequency'], data[ ['Intercept', 'Temperature'] ], family=sm.families.Binomial(sm.families.links.logit)).fit()
# logit_model=sm.Logit(data["Frequency"],data[["Intercept","Temperature"]]).fit() logmodel.summary()
logmodel=sm.GLM(data['Frequency'], data[['Intercept','Temperature']], family=sm.families.Binomial(sm.families.links.logit)).fit() #+end_src
logmodel.summary() #+RESULTS:
#+end_src #+begin_example
Generalized Linear Model Regression Results
#+RESULTS: ==============================================================================
#+begin_example Dep. Variable: Frequency No. Observations: 23
Generalized Linear Model Regression Results Model: GLM Df Residuals: 21
============================================================================== Model Family: Binomial Df Model: 1
Dep. Variable: Frequency No. Observations: 7 Link Function: logit Scale: 1.0000
Model: GLM Df Residuals: 5 Method: IRLS Log-Likelihood: -3.9210
Model Family: Binomial Df Model: 1 Date: Mon, 20 Jul 2020 Deviance: 3.0144
Link Function: logit Scale: 1.0 Time: 11:22:26 Pearson chi2: 5.00
Method: IRLS Log-Likelihood: -3.6370 No. Iterations: 6
Date: Fri, 20 Jul 2018 Deviance: 3.3763 Covariance Type: nonrobust
Time: 16:56:08 Pearson chi2: 0.236 ===============================================================================
No. Iterations: 5 coef std err z P>|z| [0.025 0.975]
=============================================================================== -------------------------------------------------------------------------------
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
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 #+end_example
===============================================================================
#+end_example L'estimateur le plus probable du paramètre de température est 0.0014
et l'erreur standard de cet estimateur est de 0.122, autrement dit on
L'estimateur le plus probable du paramètre de température est 0.0014 ne peut pas distinguer d'impact particulier et il faut prendre nos
et l'erreur standard de cet estimateur est de 0.122, autrement dit on estimations avec des pincettes.
ne peut pas distinguer d'impact particulier et il faut prendre nos
estimations avec des pincettes.
* Estimation de la probabilité de dysfonctionnant des joints toriques * Estimation de la probabilité de dysfonctionnant des joints toriques
La température prévue le jour du décollage est de 31°F. Essayons La température prévue le jour du décollage est de 31°F. Essayons
d'estimer la probabilité de dysfonctionnement des joints toriques à d'estimer la probabilité de dysfonctionnement des joints toriques à
cette température à partir du modèle que nous venons de construire: cette température à partir du modèle que nous venons de construire:
#+begin_src python :results output file :var matplot_lib_filename="proba_estimate_python.png" :exports both :session *python* #+begin_src python :results output file :var matplot_lib_filename="proba_estimate_python.png" :exports both :session *python*
import matplotlib.pyplot as plt import matplotlib.pyplot as plt
data_pred = pd.DataFrame({'Temperature': np.linspace(start=30, stop=90, num=121), 'Intercept': 1}) data_pred = pd.DataFrame({'Temperature': np.linspace(start=30, stop=90, num=121), 'Intercept': 1})
data_pred['Frequency'] = logmodel.predict(data_pred[['Intercept','Temperature']]) prediction = logmodel.get_prediction(data_pred[ ['Intercept','Temperature'] ]).summary_frame()
data_pred.plot(x="Temperature",y="Frequency",kind="line",ylim=[0,1]) data_pred['Frequency'] = prediction["mean"]
plt.scatter(x=data["Temperature"],y=data["Frequency"]) data_pred.plot(x="Temperature",y="Frequency",kind="line",ylim=[0,1])
plt.grid(True) plt.scatter(x=data["Temperature"], y=data["Frequency"])
plt.fill_between(
plt.savefig(matplot_lib_filename) data_pred["Temperature"],
print(matplot_lib_filename) prediction["mean_ci_lower"],
#+end_src prediction["mean_ci_upper"],
alpha=0.2,
#+RESULTS: )
[[file:proba_estimate_python.png]] plt.grid(True)
Comme on pouvait s'attendre au vu des données initiales, la plt.savefig(matplot_lib_filename)
température n'a pas d'impact notable sur la probabilité d'échec des print(matplot_lib_filename)
joints toriques. Elle sera d'environ 0.2, comme dans les essais #+end_src
précédents où nous il y a eu défaillance d'au moins un joint. Revenons
à l'ensemble des données initiales pour estimer la probabilité de #+RESULTS:
défaillance d'un joint: [[file:proba_estimate_python.png]]
#+begin_src python :results output :session *python* :exports both Hmm, l'intervalle de confiance semble très grand. Peut-on déduire quoi que ce
data = pd.read_csv("shuttle.csv") soit de ces données ?
print(np.sum(data.Malfunction)/np.sum(data.Count))
#+end_src
#+RESULTS:
: 0.06521739130434782
Cette probabilité est donc d'environ $p=0.065$, sachant qu'il existe
un joint primaire un joint secondaire sur chacune des trois parties du
lançeur, la probabilité de défaillance des deux joints d'un lançeur
est de $p^2 \approx 0.00425$. La probabilité de défaillance d'un des
lançeur est donc de $1-(1-p^2)^3 \approx 1.2%$. Ça serait vraiment
pas de chance... Tout est sous contrôle, le décollage peut donc avoir
lieu demain comme prévu.
Seulement, le lendemain, la navette Challenger explosera et emportera
avec elle ses sept membres d'équipages. L'opinion publique est
fortement touchée et lors de l'enquête qui suivra, la fiabilité des
joints toriques sera directement mise en cause. Au delà des problèmes
de communication interne à la NASA qui sont pour beaucoup dans ce
fiasco, l'analyse précédente comporte (au moins) un petit
problème... Saurez-vous le trouver ? Vous êtes libre de modifier cette
analyse et de regarder ce jeu de données sous tous les angles afin
d'expliquer ce qui ne va pas.
module2/exo5/freq_temp_python.png
View replaced file @ a00792ce View file @ db834875
module2/exo5/freq_temp_python.png

12.3 KB | W: | H:

module2/exo5/freq_temp_python.png

86.3 KB | W: | H:

module2/exo5/freq_temp_python.png
module2/exo5/freq_temp_python.png
module2/exo5/freq_temp_python.png
module2/exo5/freq_temp_python.png
  • 2-up
  • Swipe
  • Onion skin
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment