From 16b6a365d5ad1a165483c68e9540532004065f2f Mon Sep 17 00:00:00 2001 From: 9b01737421a268a54ae68902c9bdbf26 <9b01737421a268a54ae68902c9bdbf26@app-learninglab.inria.fr> Date: Sat, 3 May 2025 16:04:27 +0000 Subject: [PATCH] Add new file --- module2/exo5/exo5_python_fr.org | 440 ++++++++++++++++++++++++++++++++ 1 file changed, 440 insertions(+) create mode 100644 module2/exo5/exo5_python_fr.org diff --git a/module2/exo5/exo5_python_fr.org b/module2/exo5/exo5_python_fr.org new file mode 100644 index 0000000..99402df --- /dev/null +++ b/module2/exo5/exo5_python_fr.org @@ -0,0 +1,440 @@ +#+TITLE: Analyse du risque de défaillance des joints toriques de la navette Challenger +#+AUTHOR: Arnaud Legrand +#+LANGUAGE: fr + +#+HTML_HEAD: +#+HTML_HEAD: +#+HTML_HEAD: +#+HTML_HEAD: +#+HTML_HEAD: +#+HTML_HEAD: + +#+LATEX_HEADER: \usepackage{a4} +#+LATEX_HEADER: \usepackage[french]{babel} + +# #+PROPERTY: header-args :session :exports both + +Le 27 Janvier 1986, veille du décollage de la navette /Challenger/, 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. + +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. + + +* Première analyse + +** Chargement des données +Nous commençons donc par charger ces données: +#+begin_src python :results value :session first :exports both +import numpy as np +import pandas as pd +data = pd.read_csv("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 +#+end_example + +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 +température (en Fahrenheit) et la pression (en psi), et enfin le +nombre de dysfonctionnements relevés. + +** 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 +dysfonctionnements, nous nous concentrons sur les expériences où au +moins un joint a été défectueux. + +#+begin_src python :results value :session first :exports both +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 + +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 first +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) +#+end_src + +#+RESULTS: +[[file:freq_temp_python.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 +dysfonctionnements d'un joint. + +** Estimation de l'influence de la température + +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. + +#+begin_src python :results value :session first :exports both +import statsmodels.api as sm + +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() +#+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 + +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. + +** Estimation de la probabilité de dysfonctionnant des joints toriques +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 first +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.plot(x="Temperature",y="Frequency",kind="line",ylim=[0,1]) +plt.scatter(x=data["Temperature"],y=data["Frequency"]) +plt.grid(True) + +plt.savefig(matplot_lib_filename) +print(matplot_lib_filename) +#+end_src + +#+RESULTS: +[[file:proba_estimate_python.png]] + +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 +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 first :exports both +data = pd.read_csv("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 +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. + + + +* Nouvelle analyse + +Nous commençons donc par charger ces données: +#+begin_src python :results value :session second :exports both +import numpy as np +import pandas as pd +import matplotlib.pyplot as plt +import statsmodels.api as sm + +data = pd.read_csv("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 +#+end_example + + +Tracons, pour ces 3 jeux de données, le nombre de Malfunction en fonction de la température sur le même graphique, une couleur par pression. On le nomme "P_temp_python.png". + +#+begin_src python :results output file :var matplot_lib_filename="P_temp_python.png" :exports both :session second + +dataP50 = data[data.Pressure==50] +dataP100 = data[data.Pressure==100] +dataP200 = data[data.Pressure==200] + +plt.clf() +dataP50.plot(x="Temperature",y="Malfunction",kind="scatter",ylim=[-1,6],color='blue') +dataP100.plot(x="Temperature",y="Malfunction",kind="scatter",ylim=[-1,6],color='red',ax=plt.gca()) +dataP200.plot(x="Temperature",y="Malfunction",kind="scatter",ylim=[-1,6],color='green',ax=plt.gca()) +plt.grid(True) +plt.savefig(matplot_lib_filename) +print(matplot_lib_filename) +#+end_src + +#+RESULTS: +[[file:P_temp_python.png]] + + +** Regression logistique + +Nous allons maintenant effectuer 4 régressions logistiques + +*** Régression logistique sur les données P50 + +#+begin_src python :results value :session second :exports both + +# Filtrer les données et créer une copie +dataP50 = data[data.Pressure == 50].copy() + +# Calcul des colonnes nécessaires +dataP50["Success"] = dataP50["Count"] - dataP50["Malfunction"] +dataP50["Intercept"] = 1 + +# Création du modèle avec une instance de sm.families.links.Logit() +logmodelP50 = sm.GLM( + dataP50['Malfunction'], + dataP50[['Intercept', 'Temperature']], + family=sm.families.Binomial(link=sm.families.links.Logit()) +).fit() + +# Résumé du modèle +logmodelP50.summary() + +# Générer des données pour la prédiction (de 25 à 90) +data_pred50 = pd.DataFrame({ + 'Temperature': np.linspace(25, 90, 100), + 'Intercept': 1 +}) + +# Ajouter les prédictions au DataFrame +data_pred50['Malfunction_Prob'] = logmodelP50.predict(data_pred50[['Intercept', 'Temperature']]) + +# Tracer les points observés +plt.figure() # Nouvelle figure +plt.scatter(dataP50['Temperature'], dataP50['Malfunction'], color='blue', label='Observations') + +# Tracer la courbe de régression logistique +plt.plot(data_pred50['Temperature'], data_pred50['Malfunction_Prob'], color='red', label='Régression logistique') + +# Ajouter des labels et une légende +plt.xlabel('Température') +plt.ylabel('Probabilité de dysfonctionnement') +plt.title('Régression logistique sur les données P50') +plt.legend() +plt.grid(True) + +# Créer une image regression_logistique_P50.png +plt.savefig("regression_logistique_P50.png") +print("regression_logistique_P50.png") +#+end_src + +#+RESULTS: +[[file:regression_logistique_P50.png]] + + +*** Régression logistique sur les données P200 + +#+begin_src python :results value :session second :exports both +#+begin_src python :results value :session second :exports both +# Filtrer les données et créer une copie +dataP200 = data[data.Pressure == 200].copy() + +# Calcul des colonnes nécessaires +dataP200["Success"] = dataP200["Count"] - dataP200["Malfunction"] +dataP200["Intercept"] = 1 +dataP200["Proportion"] = dataP200["Malfunction"] / dataP200["Count"] + +# Création du modèle avec une instance de sm.families.links.Logit() +logmodelP200 = sm.GLM( + dataP200['Proportion'], + dataP200[['Intercept', 'Temperature']], + family=sm.families.Binomial(link=sm.families.links.Logit()) +).fit() + +# Résumé du modèle +logmodelP200.summary() + +# Générer des données pour la prédiction (de 25 à 90) +data_pred200 = pd.DataFrame({ + 'Temperature': np.linspace(25, 90, 100), + 'Intercept': 1 +}) +data_pred200['Proportion_Pred'] = logmodelP200.predict(data_pred200[['Intercept', 'Temperature']]) + +# Tracer les points observés +plt.figure() # Nouvelle figure +plt.scatter(dataP200['Temperature'], dataP200['Proportion'], color='blue', label='Observations') + +# Tracer la courbe de régression logistique +plt.plot(data_pred200['Temperature'], data_pred200['Proportion_Pred'], color='red', label='Régression logistique') + +# Ajouter des labels et une légende +plt.xlabel('Température') +plt.ylabel('Proportion de dysfonctionnements') +plt.title('Régression logistique sur les données P200') +plt.legend() +plt.grid(True) + +# Créer une image regression_logistique_P200.png +plt.savefig("regression_logistique_P200.png") +print("regression_logistique_P200.png") +#+end_src + + + +*** Régression logistique sur toutes données + +#+begin_src python :results value :session second :exports both +# Filtrer les données et créer une copie +dataAll = data.copy() +# Calcul des colonnes nécessaires +dataAll["Success"] = dataAll["Count"] - dataAll["Malfunction"] +dataAll["Intercept"] = 1 +# Création du modèle avec une instance de sm.families.links.Logit() +logmodelAll = sm.GLM( + dataAll['Malfunction'], + dataAll[['Intercept', 'Temperature']], + family=sm.families.Binomial(link=sm.families.links.Logit()) +).fit() +# Résumé du modèle +logmodelAll.summary() +# Générer des données pour la prédiction (de 25 à 90) +data_predAll = pd.DataFrame({ + 'Temperature': np.linspace(25, 90, 100), + 'Intercept': 1 +}) +# Ajouter les prédictions au DataFrame +data_predAll['Malfunction_Prob'] = logmodelAll.predict(data_predAll[['Intercept', 'Temperature']]) +# Tracer les points observés +plt.figure() # Nouvelle figure +plt.scatter(dataAll['Temperature'], dataAll['Malfunction'], color='blue', label='Observations') +# Tracer la courbe de régression logistique +plt.plot(data_predAll['Temperature'], data_predAll['Malfunction_Prob'], color='red', label='Régression logistique') +# Ajouter des labels et une légende +plt.xlabel('Température') +plt.ylabel('Probabilité de dysfonctionnement') +plt.title('Régression logistique sur toutes les données') +# Ajouter une légende +plt.legend() +# Créer une image regression_logistique_All.png +plt.savefig("regression_logistique_All.png") +print("regression_logistique_All.png") +#+end_src + +#+RESULTS: +[[file:regression_logistique_All.png]] + + -- 2.18.1