diff --git a/module2/exo5/exo5_fr_COR.Rmd b/module2/exo5/exo5_fr_COR.Rmd new file mode 100644 index 0000000000000000000000000000000000000000..0cb1894023f67ae4d4a552e704fbec7379e58e70 --- /dev/null +++ b/module2/exo5/exo5_fr_COR.Rmd @@ -0,0 +1,146 @@ +--- +title: "Analyse du risque de d?faillance des joints toriques de la navette Challenger - CORRECTION" +author: "Arnaud Legrand & cor.FMG" +date: "14/12/2021" +output: html_document +--- + +```{r setup, include=FALSE} +knitr::opts_chunk$set(echo = TRUE) +``` + +# CORRECTION PERSONNELLE +On va essayer de refaire l'analyse avec les éléments qui nous semblent "légers" : + +\- la disparition des données de malfonctions jugées inutiles, qui vont réduire +mécaniquement l’incertitude et biaisé la prise de décision + +\- la sortie arbitraire du facteur pression dans la regression logistique + +\- ou alternativement (1), tenter de réduire la difficulté d’analyse des +résultats de la régression par une approche plus conservatrice consistant à +détecter toute malfonction comme la survenue d’un évènement à éviter, à partir +du recodage de la variable en une variable binaire/booléenne 0/1 pure (=\>1 pour +la moindre malfonction de joint) =\> permettrait d’avoir **exp(coefficient +regression log) = odds-ratio qui associe la variable explicative à la variable +expliquée**) + + +## Chargement et visualisation du dataset +On commence de la même façon par le chargement des données /création du dataset de travail, mais nous ne retirons pas les données jugées valables arbitrairement. Nous conservons toutes les données de disfonctionnement des joints. Cela augmente la taille de l'échantillon et devrait réduire l'incertitude générale. + +On les affiches ensuite +```{r} +data_f = read.csv("shuttle.csv",header=T) +data_f +``` +On poursuit pour l'instant la même procédure d'analyse. +On va donc représenter les variations de dysfonctionnement des joints avec la température, sur les données non-tromquées +```{r} +plot(data=data_f, Malfunction/Count ~ Temperature, ylim=c(0,1)) +``` + +C'est pas super évident, mais il semble malgré tout que la température "réduit" les dysfonctionnement. + +Ce qui est visible en revanche, si on essaie d'agrandir les échelles, c'est que nous n'avons qu'une plage très ressérées d'observations, représentant des évènements de températures assez hohmogènes + + ```{r} + plot(data=data_f, Malfunction/Count~Temperature, xlim=c(30, 90), ylim=c(0,1)) + ``` + + Il était aussi possible de réaliser une visualisation plus précise, qui intégrait directement la regression logistique + + ```{r} + library(ggplot2) + ggplot(data_f,aes(x=Temperature, y=Malfunction/Count)) + geom_point(alpha=.3, size=3) + + theme_bw() + geom_smooth(method = "glm", method.args = list(family = "binomial")) + xlim(30,90) + ``` + +On voit ici que la tendancce que nous donne la regression n'est déjà plus aussi sécurisante que les simple plots... + + + +### Estimation de l'impact de la température + Avec le modèle complet + + ```{r} +logistic_reg1 = glm(data=data_f, Malfunction/Count ~ Temperature, weights=Count, + family=binomial(link='logit')) +summary(logistic_reg1) +``` +**La température est devenue significative !** +L'analyse des coefficients restent cependant toujours compliqués (pas une véritable variable binaire qui permettrait d'avoir des odds-ratios par exponentielle des coefficients), bien que le signe de trompe pas : relation négative. Donc, plus la température monte, plus la proportion de défaillance des joints toriques baisse. + + +##### Essai twoby et Epi +Juste pour le fun +```{r} + +library(Epi) +malf_c <- (data_f$Malfunction/data_f$Count) +twoby2(1-malf_c, 1-data_f$Temperature) +``` +Ben ça marche pas, pas des variables binaires... :-) + +```{r} +exp(-0.11560) +``` +Envrion 0.9... Cela caractérise bien la perte de proba à mesure que la température augmente (et inversement, 1,15 pour si la température chute) + +### Prédiction/Estimation de la probabilité de dysfonctionnant des joints toriques +On reprend le code proposé en modifiant les variables de calcul. +Ici, on crée une séquence par "0,5" de 30 à 90 (abscisses axe), et on crée un jeu de données prédites par le modèle de la regression logistique (ordonnées axe). +Ensuite, on le représente graphiqement +```{r} +tempv = seq(from=30, to=90, by = .5) +rmv <- predict(logistic_reg1,list(Temperature=tempv),type="response") +plot(tempv,rmv,type="l",ylim=c(0,1)) +points(data=data_f, Malfunction/Count ~ Temperature) +``` + +Là, cela devient évident! Avec des température de 30 degré farenheit, le risque de défaillance est très grand ! Nous pourrions nous en tenir là, et estimer les defauts par l'absence d'un dataset complet... + + +## Intégration de la variable pression + +```{r} +logistic_reg2 = glm(data=data_f, Malfunction/Count ~ Temperature + Pressure + Temperature:Pressure, weights=Count, + family=binomial(link='logit')) +summary(logistic_reg2) +``` +Pas d'interaction, ni d'effet simple : le retrait de la pression n'était donc pas dommageable selon les données disponibles + +## Modèle avec variables binaires : approche plus conservatrice +L'idée est de recoder la variable "incident" (ici la défaillance d'un joint torique) en une variablel binaire/bouléenne. Ainsi, on comptabilise le risque de défaillance d'au moins 1 joint torique (avec l'idée qu'un seul suffit à mettre en danger la navette) + +Création de la variable binaire + +```{r} +data_f$Malfunction_b <- ifelse(data_f$Malfunction>=1, 1,0) +``` + +On refait la regression + ```{r} +logistic_reg3 = glm(data=data_f, Malfunction_b ~ Temperature, weights=Count, + family=binomial(link='logit')) +summary(logistic_reg3) +``` + +### Analyse des odds-ratios +La librairie twoby ne fonctonnera pas. +Mais ici, l'exp(coef) = OR! + +```{r} +exp(0.23216) +``` +On trouve qu'il y a 1,3 fois plus de chance de défaillance avec la baisse de la temperature. +Cette approche, comme on s'y attendait, montre plus facilemnt encore le risque qu'encourt la navette avec une plus forte probabilité signifiative qu'au moins 1 joint torique ne lâche avec la baisse de la température lors du lancement. + +```{r} +tempv = seq(from=30, to=90, by = .5) +rmv <- predict(logistic_reg3,list(Temperature=tempv),type="response") +plot(tempv,rmv,type="l",ylim=c(0,1)) +points(data=data_f, Malfunction_b ~ Temperature) +``` + +**En fait, si on en croit les prédictions du modèle, l'accident (défaillance d'au moins 1 joint torique) était inévitable en dessous de 50°F...** diff --git a/module2/exo5/exo5_fr_COR.html b/module2/exo5/exo5_fr_COR.html new file mode 100644 index 0000000000000000000000000000000000000000..d359836e17164faf9fed9840dc0eddc569ecc95e --- /dev/null +++ b/module2/exo5/exo5_fr_COR.html @@ -0,0 +1,438 @@ + + + + + + + + + + + + + + +Analyse du risque de d?faillance des joints toriques de la navette Challenger - CORRECTION + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + +
+ + + + + + + +
+

CORRECTION PERSONNELLE

+

On va essayer de refaire l’analyse avec les éléments qui nous semblent “légers” :

+

- la disparition des données de malfonctions jugées inutiles, qui vont réduire mécaniquement l’incertitude et biaisé la prise de décision

+

- la sortie arbitraire du facteur pression dans la regression logistique

+

- ou alternativement (1), tenter de réduire la difficulté d’analyse des résultats de la régression par une approche plus conservatrice consistant à détecter toute malfonction comme la survenue d’un évènement à éviter, à partir du recodage de la variable en une variable binaire/booléenne 0/1 pure (=>1 pour la moindre malfonction de joint) => permettrait d’avoir exp(coefficient regression log) = odds-ratio qui associe la variable explicative à la variable expliquée)

+
+

Chargement et visualisation du dataset

+

On commence de la même façon par le chargement des données /création du dataset de travail, mais nous ne retirons pas les données jugées valables arbitrairement. Nous conservons toutes les données de disfonctionnement des joints. Cela augmente la taille de l’échantillon et devrait réduire l’incertitude générale.

+

On les affiches ensuite

+
data_f = read.csv("shuttle.csv",header=T)
+data_f
+
##        Date Count Temperature Pressure Malfunction
+## 1   4/12/81     6          66       50           0
+## 2  11/12/81     6          70       50           1
+## 3   3/22/82     6          69       50           0
+## 4  11/11/82     6          68       50           0
+## 5   4/04/83     6          67       50           0
+## 6   6/18/82     6          72       50           0
+## 7   8/30/83     6          73      100           0
+## 8  11/28/83     6          70      100           0
+## 9   2/03/84     6          57      200           1
+## 10  4/06/84     6          63      200           1
+## 11  8/30/84     6          70      200           1
+## 12 10/05/84     6          78      200           0
+## 13 11/08/84     6          67      200           0
+## 14  1/24/85     6          53      200           2
+## 15  4/12/85     6          67      200           0
+## 16  4/29/85     6          75      200           0
+## 17  6/17/85     6          70      200           0
+## 18  7/29/85     6          81      200           0
+## 19  8/27/85     6          76      200           0
+## 20 10/03/85     6          79      200           0
+## 21 10/30/85     6          75      200           2
+## 22 11/26/85     6          76      200           0
+## 23  1/12/86     6          58      200           1
+

On poursuit pour l’instant la même procédure d’analyse. On va donc représenter les variations de dysfonctionnement des joints avec la température, sur les données non-tromquées

+
plot(data=data_f, Malfunction/Count ~ Temperature, ylim=c(0,1))
+

+

C’est pas super évident, mais il semble malgré tout que la température “réduit” les dysfonctionnement.

+

Ce qui est visible en revanche, si on essaie d’agrandir les échelles, c’est que nous n’avons qu’une plage très ressérées d’observations, représentant des évènements de températures assez hohmogènes

+
plot(data=data_f, Malfunction/Count~Temperature, xlim=c(30, 90), ylim=c(0,1))
+

+

Il était aussi possible de réaliser une visualisation plus précise, qui intégrait directement la regression logistique

+
library(ggplot2)
+
## Registered S3 methods overwritten by 'tibble':
+##   method     from  
+##   format.tbl pillar
+##   print.tbl  pillar
+
ggplot(data_f,aes(x=Temperature, y=Malfunction/Count)) + geom_point(alpha=.3, size=3) +
+theme_bw() + geom_smooth(method = "glm", method.args = list(family = "binomial")) + xlim(30,90)
+
## `geom_smooth()` using formula 'y ~ x'
+
## Warning in eval(family$initialize): non-integer #successes in a binomial glm!
+

+

On voit ici que la tendancce que nous donne la regression n’est déjà plus aussi sécurisante que les simple plots…

+
+

Estimation de l’impact de la température

+

Avec le modèle complet

+
logistic_reg1 = glm(data=data_f, Malfunction/Count ~ Temperature, weights=Count, 
+                  family=binomial(link='logit'))
+summary(logistic_reg1)
+
## 
+## Call:
+## glm(formula = Malfunction/Count ~ Temperature, family = binomial(link = "logit"), 
+##     data = data_f, weights = Count)
+## 
+## Deviance Residuals: 
+##      Min        1Q    Median        3Q       Max  
+## -0.95227  -0.78299  -0.54117  -0.04379   2.65152  
+## 
+## Coefficients:
+##             Estimate Std. Error z value Pr(>|z|)  
+## (Intercept)  5.08498    3.05247   1.666   0.0957 .
+## Temperature -0.11560    0.04702  -2.458   0.0140 *
+## ---
+## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
+## 
+## (Dispersion parameter for binomial family taken to be 1)
+## 
+##     Null deviance: 24.230  on 22  degrees of freedom
+## Residual deviance: 18.086  on 21  degrees of freedom
+## AIC: 35.647
+## 
+## Number of Fisher Scoring iterations: 5
+

La température est devenue significative ! L’analyse des coefficients restent cependant toujours compliqués (pas une véritable variable binaire qui permettrait d’avoir des odds-ratios par exponentielle des coefficients), bien que le signe de trompe pas : relation négative. Donc, plus la température monte, plus la proportion de défaillance des joints toriques baisse.

+
+
Essai twoby et Epi
+

Juste pour le fun

+
library(Epi)
+malf_c <- (data_f$Malfunction/data_f$Count)
+twoby2(1-malf_c, 1-data_f$Temperature)
+
## 2 by 2 table analysis: 
+## ------------------------------------------------------ 
+## Outcome   : -80 
+## Comparing : 0.666666666666667 vs. 0.833333333333333 
+## 
+##                   -80 -78    P(-80) 95% conf. interval
+## 0.666666666666667   0   0       NaN       NaN      NaN
+## 0.833333333333333   0   0       NaN       NaN      NaN
+## 
+##                                   95% conf. interval
+##              Relative Risk:   NaN       NaN      NaN
+##          Sample Odds Ratio:   NaN       NaN      NaN
+## Conditional MLE Odds Ratio:     0         0      Inf
+##     Probability difference:   NaN       NaN      NaN
+## 
+##              Exact P-value: 1.0000 
+##         Asymptotic P-value:   NaN 
+## ------------------------------------------------------
+

Ben ça marche pas, pas des variables binaires… :-)

+
exp(-0.11560)
+
## [1] 0.8908315
+

Envrion 0.9… Cela caractérise bien la perte de proba à mesure que la température augmente (et inversement, 1,15 pour si la température chute)

+
+
+
+

Prédiction/Estimation de la probabilité de dysfonctionnant des joints toriques

+

On reprend le code proposé en modifiant les variables de calcul. Ici, on crée une séquence par “0,5” de 30 à 90 (abscisses axe), et on crée un jeu de données prédites par le modèle de la regression logistique (ordonnées axe). Ensuite, on le représente graphiqement

+
tempv = seq(from=30, to=90, by = .5)
+rmv <- predict(logistic_reg1,list(Temperature=tempv),type="response")
+plot(tempv,rmv,type="l",ylim=c(0,1))
+points(data=data_f, Malfunction/Count ~ Temperature)
+

+

Là, cela devient évident! Avec des température de 30 degré farenheit, le risque de défaillance est très grand ! Nous pourrions nous en tenir là, et estimer les defauts par l’absence d’un dataset complet…

+
+
+
+

Intégration de la variable pression

+
logistic_reg2 = glm(data=data_f, Malfunction/Count ~ Temperature + Pressure + Temperature:Pressure, weights=Count, 
+                   family=binomial(link='logit'))
+summary(logistic_reg2)
+
## 
+## Call:
+## glm(formula = Malfunction/Count ~ Temperature + Pressure + Temperature:Pressure, 
+##     family = binomial(link = "logit"), data = data_f, weights = Count)
+## 
+## Deviance Residuals: 
+##     Min       1Q   Median       3Q      Max  
+## -1.0279  -0.6600  -0.5141  -0.1744   2.3724  
+## 
+## Coefficients:
+##                        Estimate Std. Error z value Pr(>|z|)
+## (Intercept)          -21.765547  46.285607  -0.470    0.638
+## Temperature            0.252463   0.662989   0.381    0.703
+## Pressure               0.130857   0.232643   0.562    0.574
+## Temperature:Pressure  -0.001769   0.003336  -0.530    0.596
+## 
+## (Dispersion parameter for binomial family taken to be 1)
+## 
+##     Null deviance: 24.230  on 22  degrees of freedom
+## Residual deviance: 16.257  on 19  degrees of freedom
+## AIC: 37.817
+## 
+## Number of Fisher Scoring iterations: 6
+

Pas d’interaction, ni d’effet simple : le retrait de la pression n’était donc pas dommageable selon les données disponibles

+
+
+

Modèle avec variables binaires : approche plus conservatrice

+

L’idée est de recoder la variable “incident” (ici la défaillance d’un joint torique) en une variablel binaire/bouléenne. Ainsi, on comptabilise le risque de défaillance d’au moins 1 joint torique (avec l’idée qu’un seul suffit à mettre en danger la navette)

+

Création de la variable binaire

+
data_f$Malfunction_b <- ifelse(data_f$Malfunction>=1, 1,0)
+

On refait la regression

+
logistic_reg3 = glm(data=data_f, Malfunction_b ~ Temperature, weights=Count, 
+                  family=binomial(link='logit'))
+summary(logistic_reg3)
+
## 
+## Call:
+## glm(formula = Malfunction_b ~ Temperature, family = binomial(link = "logit"), 
+##     data = data_f, weights = Count)
+## 
+## Deviance Residuals: 
+##     Min       1Q   Median       3Q      Max  
+## -2.5992  -1.8647  -0.9266   1.1081   5.4318  
+## 
+## Coefficients:
+##             Estimate Std. Error z value Pr(>|z|)    
+## (Intercept) 15.04290    3.01231   4.994 5.92e-07 ***
+## Temperature -0.23216    0.04419  -5.254 1.49e-07 ***
+## ---
+## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
+## 
+## (Dispersion parameter for binomial family taken to be 1)
+## 
+##     Null deviance: 169.60  on 22  degrees of freedom
+## Residual deviance: 121.89  on 21  degrees of freedom
+## AIC: 125.89
+## 
+## Number of Fisher Scoring iterations: 5
+
+

Analyse des odds-ratios

+

La librairie twoby ne fonctonnera pas. Mais ici, l’exp(coef) = OR!

+
exp(0.23216)
+
## [1] 1.261322
+

On trouve qu’il y a 1,3 fois plus de chance de défaillance avec la baisse de la temperature. Cette approche, comme on s’y attendait, montre plus facilemnt encore le risque qu’encourt la navette avec une plus forte probabilité signifiative qu’au moins 1 joint torique ne lâche avec la baisse de la température lors du lancement.

+
tempv = seq(from=30, to=90, by = .5)
+rmv <- predict(logistic_reg3,list(Temperature=tempv),type="response")
+plot(tempv,rmv,type="l",ylim=c(0,1))
+points(data=data_f, Malfunction_b ~ Temperature)
+

+

En fait, si on en croit les prédictions du modèle, l’accident (défaillance d’au moins 1 joint torique) était inévitable en dessous de 50°F…

+
+
+
+ + + + +
+ + + + + + + + + + + + + + +