From 5bd15e6a368fa759a5cb5c2d432eebbef28b53bd Mon Sep 17 00:00:00 2001
From: 5785ce2009281db763af8064eb1752c5
<5785ce2009281db763af8064eb1752c5@app-learninglab.inria.fr>
Date: Sun, 14 Apr 2024 16:44:05 +0000
Subject: [PATCH] Replace exo5_python_fr.org
---
module2/exo5/exo5_python_fr.org | 300 +++++++++++++++++++++++---------
1 file changed, 219 insertions(+), 81 deletions(-)
diff --git a/module2/exo5/exo5_python_fr.org b/module2/exo5/exo5_python_fr.org
index afff5e0..b230237 100644
--- a/module2/exo5/exo5_python_fr.org
+++ b/module2/exo5/exo5_python_fr.org
@@ -1,7 +1,16 @@
#+TITLE: Analyse du risque de défaillance des joints toriques de la navette Challenger
+:preambule:
#+AUTHOR: Arnaud Legrand
#+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:
#+HTML_HEAD:
#+HTML_HEAD:
@@ -11,6 +20,38 @@
#+LATEX_HEADER: \usepackage{a4}
#+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
@@ -33,39 +74,39 @@ Challenger.
* Chargement des 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 pandas as pd
-data = pd.read_csv("shuttle.csv")
+data = pd.read_csv("module2_exo5_shuttle.csv")
data
#+end_src
#+RESULTS:
#+begin_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
+ 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
#+end_example
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
dysfonctionnements, nous nous concentrons sur les expériences où au
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
#+end_src
#+RESULTS:
-: Date Count Temperature Pressure Malfunction
-: 1 11/12/81 6 70 50 1
-: 8 2/03/84 6 57 200 1
-: 9 4/06/84 6 63 200 1
-: 10 8/30/84 6 70 200 1
-: 13 1/24/85 6 53 200 2
-: 20 10/30/85 6 75 200 2
-: 22 1/12/86 6 58 200 1
+: Date Count Temperature Pressure Malfunction Frequency
+: 1 11/12/81 6 70 50 1 0.166667
+: 8 2/03/84 6 57 200 1 0.166667
+: 9 4/06/84 6 63 200 1 0.166667
+: 10 8/30/84 6 70 200 1 0.166667
+: 13 1/24/85 6 53 200 2 0.333333
+: 20 10/30/85 6 75 200 2 0.333333
+: 22 1/12/86 6 58 200 1 0.166667
Très bien, nous avons une variabilité de température importante mais
la pression est quasiment toujours égale à 200, ce qui devrait
simplifier l'analyse.
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
#+RESULTS:
-[[file:freq_temp_python.png]]
+[[./.ob-jupyter/3d4ddf0ba39f4d97099fa73e230ec4ef702db5f1.png]]
+
+
À 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
@@ -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
régression logistique.
-#+begin_src python :results value :session *python* :exports both
-import statsmodels.api as sm
+#+begin_src jupyter-python :results value :exports both
+ import pandas as pd
+ import statsmodels.api as sm
-data["Success"]=data.Count-data.Malfunction
-data["Intercept"]=1
+ data["Success"]=data.Count-data.Malfunction
+ data["Intercept"]=1
-
-# 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()
-
-logmodel.summary()
+ logit_model=sm.Logit(data["Frequency"],data[["Intercept","Temperature"]]).fit()
+
+ logit_model.summary()
#+end_src
#+RESULTS:
-#+begin_example
- 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.0
-Method: IRLS Log-Likelihood: -3.6370
-Date: Fri, 20 Jul 2018 Deviance: 3.3763
-Time: 16:56:08 Pearson chi2: 0.236
-No. Iterations: 5
-===============================================================================
- 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_example
+:RESULTS:
+: Optimization terminated successfully.
+: Current function value: 0.393555
+: Iterations 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
+: return 1 - self.llf/self.llnull
+#+begin_export html
+
+ Logit Regression Results
+ Logit Regression Results
+
+ 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: 20:21:52 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
+
+#+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
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
d'estimer la probabilité de dysfonctionnement des joints toriques à
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
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])
plt.scatter(x=data["Temperature"],y=data["Frequency"])
plt.grid(True)
@@ -185,38 +265,96 @@ print(matplot_lib_filename)
#+end_src
#+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
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
à l'ensemble des données initiales pour estimer la probabilité de
défaillance d'un joint:
-#+begin_src python :results output :session *python* :exports both
-data = pd.read_csv("shuttle.csv")
+#+begin_src jupyter-python :results output :exports both
+data = pd.read_csv("module2_exo5_shuttle.csv")
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
+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
+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.
+analyse et de regarder ce jeu de données sous tous les angles afin 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]]
--
2.18.1