Skip to content

M
mooc-rr
Commit 5bd15e6a authored by 5785ce2009281db763af8064eb1752c5's avatar 5785ce2009281db763af8064eb1752c5
Browse files
Options

Replace exo5_python_fr.org

parent c502c48e
Hide whitespace changes
Inline Side-by-side
Showing with 219 additions and 81 deletions
module2/exo5/exo5_python_fr.org
View file @ 5bd15e6a
#+TITLE: Analyse du risque de défaillance des joints toriques de la navette Challenger #+TITLE: Analyse du risque de défaillance des joints toriques de la navette Challenger
:preambule:
#+AUTHOR: Arnaud Legrand #+AUTHOR: Arnaud Legrand
#+LANGUAGE: fr #+LANGUAGE: fr
#+PROPERTY: header-args :tangle yes
#+PROPERTY: header-args:jupyter-python :exports both
#+PROPERTY: header-args:jupyter-python+ :session /jpy:localhost#8889:jlab
#+PROPERTY: header-args:jupyter-R :exports both
#+PROPERTY: header-args:jupyter-R+ :session /jpy::ir
#+HTML_HEAD: <link rel="stylesheet" type="text/css" href="http://www.pirilampo.org/styles/readtheorg/css/htmlize.css"/> #+HTML_HEAD: <link rel="stylesheet" type="text/css" href="http://www.pirilampo.org/styles/readtheorg/css/htmlize.css"/>
#+HTML_HEAD: <link rel="stylesheet" type="text/css" href="http://www.pirilampo.org/styles/readtheorg/css/readtheorg.css"/> #+HTML_HEAD: <link rel="stylesheet" type="text/css" href="http://www.pirilampo.org/styles/readtheorg/css/readtheorg.css"/>
#+HTML_HEAD: <script src="https://ajax.googleapis.com/ajax/libs/jquery/2.1.3/jquery.min.js"></script> #+HTML_HEAD: <script src="https://ajax.googleapis.com/ajax/libs/jquery/2.1.3/jquery.min.js"></script>
...@@ -11,6 +20,38 @@ ...@@ -11,6 +20,38 @@
#+LATEX_HEADER: \usepackage{a4} #+LATEX_HEADER: \usepackage{a4}
#+LATEX_HEADER: \usepackage[french]{babel} #+LATEX_HEADER: \usepackage[french]{babel}
:end:
#+name: scrollImage
#+begin_src elisp :exports none
(defun org-image-resize (frame)
(setq org-export-allow-bind-keywords t)
(when (derived-mode-p 'org-mode)
(if (< (window-total-width) 80)
(setq org-image-actual-width (* 4 window-pixel-width))
(setq org-image-actual-width (* 70 (window-font-width))))
(org-redisplay-inline-images)))
(add-hook 'window-size-change-functions 'org-image-resize)
(visual-line-mode t)
(variable-pitch-mode t)
(add-hook 'org-babel-after-execute-hook 'org-redisplay-inline-images)
#+end_src
#+RESULTS: scrollImage
| org-redisplay-inline-images |
#+begin_src emacs-lisp
(setq org-image-actual-width '(200))
#+end_src
#+RESULTS:
| 200 |
# #+PROPERTY: header-args :session :exports both # #+PROPERTY: header-args :session :exports both
...@@ -33,39 +74,39 @@ Challenger. ...@@ -33,39 +74,39 @@ 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 jupyter-python :results value :exports both :display text
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("module2_exo5_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
...@@ -79,40 +120,55 @@ sur l'influence de la température ou de la pression sur les ...@@ -79,40 +120,55 @@ sur l'influence de la température ou de la pression sur les
dysfonctionnements, nous nous concentrons sur les expériences où au dysfonctionnements, nous nous concentrons sur les expériences où au
moins un joint a été défectueux. moins un joint a été défectueux.
#+begin_src python :results value :session *python* :exports both #+begin_src jupyter-python :results value :exports both :display text
data = data[data.Malfunction>0] data = data[data.Malfunction>0]
data data
#+end_src #+end_src
#+RESULTS: #+RESULTS:
: Date Count Temperature Pressure Malfunction : Date Count Temperature Pressure Malfunction Frequency
: 1 11/12/81 6 70 50 1 : 1 11/12/81 6 70 50 1 0.166667
: 8 2/03/84 6 57 200 1 : 8 2/03/84 6 57 200 1 0.166667
: 9 4/06/84 6 63 200 1 : 9 4/06/84 6 63 200 1 0.166667
: 10 8/30/84 6 70 200 1 : 10 8/30/84 6 70 200 1 0.166667
: 13 1/24/85 6 53 200 2 : 13 1/24/85 6 53 200 2 0.333333
: 20 10/30/85 6 75 200 2 : 20 10/30/85 6 75 200 2 0.333333
: 22 1/12/86 6 58 200 1 : 22 1/12/86 6 58 200 1 0.166667
Très bien, nous avons une variabilité de température importante mais Très bien, nous avons une variabilité de température importante mais
la pression est quasiment toujours égale à 200, ce qui devrait la pression est quasiment toujours égale à 200, ce qui devrait
simplifier l'analyse. simplifier l'analyse.
Comment la fréquence d'échecs varie-t-elle avec la température ? Comment la fréquence d'échecs varie-t-elle avec la température ?
#+begin_src python :results output file :var matplot_lib_filename="freq_temp_python.png" :exports both :session *python*
import matplotlib.pyplot as plt
plt.clf()
data["Frequency"]=data.Malfunction/data.Count
data.plot(x="Temperature",y="Frequency",kind="scatter",ylim=[0,1])
plt.grid(True)
plt.savefig(matplot_lib_filename)
print(matplot_lib_filename)
#+begin_src jupyter-python :exports both :results raw drawer
import matplotlib.pyplot as plt
cm = 1/2.54
plt.rcParams["figure.facecolor"] = "#A4B9C5"
#plt.rcParams['text.usetex'] = True
c_orange="#FF9A00"
c_darkblue="#5A7F96"
#+end_src
#+RESULTS:
#+begin_src jupyter-python :results raw drawer :exports both
import matplotlib.pyplot as plt
pd.options.mode.copy_on_write = True
#plt.clf()
data["Frequency"]=data["Malfunction"]/data["Count"]
data.plot(x="Temperature",y="Frequency",kind="scatter",ylim=[0,1])
plt.grid(True)
#+end_src #+end_src
#+RESULTS: #+RESULTS:
[[file:freq_temp_python.png]] [[./.ob-jupyter/3d4ddf0ba39f4d97099fa73e230ec4ef702db5f1.png]]
À première vue, ce n'est pas flagrant mais bon, essayons quand même À 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 d'estimer l'impact de la température $t$ sur la probabilité de
...@@ -128,38 +184,62 @@ température $t$ suit une loi binomiale de paramètre $n=6$ et ...@@ -128,38 +184,62 @@ 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 jupyter-python :results value :exports both
import statsmodels.api as sm import pandas as pd
import statsmodels.api as sm
data["Success"]=data.Count-data.Malfunction data["Success"]=data.Count-data.Malfunction
data["Intercept"]=1 data["Intercept"]=1
logit_model=sm.Logit(data["Frequency"],data[["Intercept","Temperature"]]).fit()
# logit_model=sm.Logit(data["Frequency"],data[["Intercept","Temperature"]]).fit()
logmodel=sm.GLM(data['Frequency'], data[['Intercept','Temperature']], family=sm.families.Binomial(sm.families.links.logit)).fit() logit_model.summary()
logmodel.summary()
#+end_src #+end_src
#+RESULTS: #+RESULTS:
#+begin_example :RESULTS:
Generalized Linear Model Regression Results : Optimization terminated successfully.
============================================================================== : Current function value: 0.393555
Dep. Variable: Frequency No. Observations: 7 : Iterations 5
Model: GLM Df Residuals: 5 : /home/vincent/.pyenv/versions/3.10.1/envs/default/lib/python3.10/site-packages/statsmodels/discrete/discrete_model.py:4469: RuntimeWarning: divide by zero encountered in scalar divide
Model Family: Binomial Df Model: 1 : return 1 - self.llf/self.llnull
Link Function: logit Scale: 1.0 #+begin_export html
Method: IRLS Log-Likelihood: -3.6370
Date: Fri, 20 Jul 2018 Deviance: 3.3763 Logit Regression Results
Time: 16:56:08 Pearson chi2: 0.236 Logit Regression Results
No. Iterations: 5
=============================================================================== Dep. Variable: Frequency No. Observations: 7
coef std err z P>|z| [0.025 0.975] Model: Logit Df Residuals: 5
------------------------------------------------------------------------------- Method: MLE Df Model: 1
Intercept -1.3895 7.828 -0.178 0.859 -16.732 13.953 Date: Sat, 13 Apr 2024 Pseudo R-squ.: inf
Temperature 0.0014 0.122 0.012 0.991 -0.238 0.240 Time: 20:21:52 Log-Likelihood: -2.7549
=============================================================================== converged: True LL-Null: 0.0000
#+end_example Covariance Type: nonrobust LLR p-value: 1.000
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
#+end_export
:END:
| Dep. Variable: | Frequency | No. Observations: 7 | |
| Model: | Logit | Df Residuals: | 5 |
| Method: | MLE | Df Model: | 1 |
| Date: | Sat, 13 Apr 2024 | Pseudo R-squ.: | inf |
| Time: | 16:44:31 | Log-Likelihood: | -2.7549 |
| converged: | True | LL-Null: | 0.0000 |
| Covariance Type: | nonrobust | LLR p-value: | 1.000 |
| | 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 |
L'estimateur le plus probable du paramètre de température est 0.0014 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 et l'erreur standard de cet estimateur est de 0.122, autrement dit on
...@@ -171,11 +251,11 @@ La température prévue le jour du décollage est de 31°F. Essayons ...@@ -171,11 +251,11 @@ 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 jupyter-python :results output file :var matplot_lib_filename="proba_estimate_python.png" :exports both
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']]) data_pred['Frequency'] = logit_model.predict(data_pred[['Intercept','Temperature']])
data_pred.plot(x="Temperature",y="Frequency",kind="line",ylim=[0,1]) data_pred.plot(x="Temperature",y="Frequency",kind="line",ylim=[0,1])
plt.scatter(x=data["Temperature"],y=data["Frequency"]) plt.scatter(x=data["Temperature"],y=data["Frequency"])
plt.grid(True) plt.grid(True)
...@@ -185,38 +265,96 @@ print(matplot_lib_filename) ...@@ -185,38 +265,96 @@ print(matplot_lib_filename)
#+end_src #+end_src
#+RESULTS: #+RESULTS:
[[file:proba_estimate_python.png]] :RESULTS:
: proba_estimate_python.png
[[./.ob-jupyter/32327fd28a2604a503ea9a032641035fd0b02e56.png]]
:END:
Comme on pouvait s'attendre au vu des données initiales, la Comme on pouvait s'attendre au vu des données initiales, la
température n'a pas d'impact notable sur la probabilité d'échec des température n'a pas d'impact notable sur la probabilité d'échec des
joints toriques. Elle sera d'environ 0.2, comme dans les essais joints toriques. Elle sera d'environ 0,2, comme dans les essais
précédents où nous il y a eu défaillance d'au moins un joint. Revenons 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 à l'ensemble des données initiales pour estimer la probabilité de
défaillance d'un joint: défaillance d'un joint:
#+begin_src python :results output :session *python* :exports both #+begin_src jupyter-python :results output :exports both
data = pd.read_csv("shuttle.csv") data = pd.read_csv("module2_exo5_shuttle.csv")
print(np.sum(data.Malfunction)/np.sum(data.Count)) print(np.sum(data.Malfunction)/np.sum(data.Count))
#+end_src #+end_src
#+RESULTS: #+RESULTS:
: 0.06521739130434782 : 0.06521739130434782
Cette probabilité est donc d'environ $p=0.065$, sachant qu'il existe 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 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 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 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 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 pas de chance... Tout est sous contrôle, le décollage peut donc avoir
lieu demain comme prévu. lieu demain comme prévu.
Seulement, le lendemain, la navette Challenger explosera et emportera Seulement, le lendemain, la navette Challenger explosera et emportera
avec elle ses sept membres d'équipages. L'opinion publique est 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 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 joints toriques sera directement mise en cause. Au-delà des problèmes
de communication interne à la NASA qui sont pour beaucoup dans ce de communication interne à la NASA qui sont pour beaucoup dans ce
fiasco, l'analyse précédente comporte (au moins) un petit fiasco, l'analyse précédente comporte (au moins) un petit
problème... Saurez-vous le trouver ? Vous êtes libre de modifier cette 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 analyse et de regarder ce jeu de données sous tous les angles afin d'expliquer ce qui ne va pas.
d'expliquer ce qui ne va pas.
* Autre analyse
Visualisation de l'ensemble du nuage avec intervalle de confiance (par
bootstrap).
#+begin_src jupyter-python :exports both :results raw drawer
from __future__ import print_function
import numpy as np
import statsmodels.api as sm
from sklearn.linear_model import LogisticRegression
from sklearn.linear_model import LinearRegression
import matplotlib.pyplot as plt
from scipy import stats
import seaborn as sns
sns.set(style='darkgrid', font_scale=1.2)
%matplotlib inline
#+end_src
#+RESULTS:
#+begin_src jupyter-python :kernel python3 :exports :results raw drawer :display text
import numpy as np
import pandas as pd
data = pd.read_csv("module2_exo5_shuttle.csv")
data["Frequency"]=data["Malfunction"]/data["Count"]
#+end_src
#+RESULTS:
#+begin_src jupyter-python :kernel python3 :exports :results raw drawer
fig, ax = plt.subplots(figsize=(10,6))
plot=sns.regplot(x='Temperature',
y='Frequency',
data=data,
logistic=True,
n_boot=500,
y_jitter=0,
ax=ax)
plt.show()
#+end_src
#+RESULTS:
[[./.ob-jupyter/dd04d70ba0976bf2224db67c19eefcbcf2398b08.png]]
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