{
"cells": [
{
"cell_type": "markdown",
"metadata": {
"hideCode": true,
"hidePrompt": true
},
"source": [
"\n",
"# Risk Analysis of the Space Shuttle: Pre-Challenger Prediction of Failure\n",
"\n",
"*Marie-Gabrielle Dondon*\n",
"\n",
"*11/10/2018*\n",
"\n",
"*This is a SAS version based on the Rmd document of Arnaud Legrand available at [https://app-learninglab.inria.fr](https://applearninglab.inria.fr/gitlab/moocrr-session1/moocrr-reproducibility-study/blob/5e298cd5b653bcb31b6f9fedcb7f75e800c2dc98/src/R/challenger.Rmd) and the SAS program of Michael Friendly available at [https://dokumen.tips](https://dokumen.tips/documents/categorical-data-analysis-with-graphics.html) (cf. page 52).*\n",
"\n",
"---"
]
},
{
"cell_type": "markdown",
"metadata": {
"hideCode": true,
"hidePrompt": true
},
"source": [
"In this document we reperform some of the analysis provided in \n",
"*Risk Analysis of the Space Shuttle: Pre-Challenger Prediction of Failure* by *Siddhartha R. Dalal, Edward B. Fowlkes, Bruce Hoadley* published in *Journal of the American Statistical Association*, Vol. 84, No. 408 (Dec., 1989), pp. 945-957 and available at http://www.jstor.org/stable/2290069. \n",
"\n",
"On the fourth page of this article, they indicate that the maximum likelihood estimates of the logistic regression using only temperature are: $\\hat{\\alpha}$ = **5.085** and $\\hat{\\beta}$ = **-0.1156** and their asymptotic standard errors are $s_{\\hat{\\alpha}}$ = **3.052** and $s_{\\hat{\\beta}}$ = **0.047**. The Goodness of fit indicated for this model was $G^2$ = **18.086** with **21** degrees of freedom. Our goal is to reproduce the computation behind these values."
]
},
{
"cell_type": "markdown",
"metadata": {
"hideCode": true,
"hidePrompt": true
},
"source": [
"The DATA step below reads the data on the number of O-ring failures and temperature for the 23 flights forwhich information was available before the Challenger launch. Our interest here is in predicting the likelihood of failures at low temperatures."
]
},
{
"cell_type": "code",
"execution_count": 1,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"SAS Connection established. Subprocess id is 3488\n",
"\n"
]
},
{
"data": {
"text/html": [
"\r\n",
"\r\n",
"
\r\n",
"\r\n",
"\r\n",
"\r\n",
"\r\n"
],
"text/plain": [
""
]
},
"execution_count": 2,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"data nasa;\n",
" set nasa;\n",
" Freq_Malfunction = Malfunction/Count;\n",
"run;\n",
"\n",
"proc gplot data=nasa;\n",
" plot Freq_Malfunction * Temperature / haxis=axis1 hminor=0\n",
" vaxis=axis2 vminor=0;\n",
" symbol1 v=dot i=none c=blue h=2;\n",
" axis1 order=(50 to 85 by 5);\n",
" axis2 order=(0 to 1 by 0.2) label=(angle=90 'Estimated Failure probability');\n",
"run;\n",
"quit;"
]
},
{
"cell_type": "markdown",
"metadata": {
"hideCode": true,
"hidePrompt": true
},
"source": [
"To obtain predicted probabilities for observations not in the original sample, create an additional data set which contains values for the independent variables in the extrapolation sample, and join these observations to the actual data set. The response variable (Malfunction) will be missing for the extrapolation sample.\n",
"\n",
"Obtain predicted values for 30-80 degrees"
]
},
{
"cell_type": "code",
"execution_count": 3,
"metadata": {},
"outputs": [
{
"data": {
"text/html": [
"\r\n",
"\r\n",
"\r\n",
"\r\n",
"\r\n",
"Sortie SAS\r\n",
"\r\n",
"\r\n",
"\r\n",
"
Sortie SAS
\r\n",
"\r\n",
"
\r\n",
"
NASA Space Shuttle O-Ring Failures
\r\n",
"
\r\n",
"
La PRINT Procédure
\r\n",
"\r\n",
"
Table WORK.NASA2
\r\n",
"
\r\n",
"
\r\n",
"\r\n",
"
\r\n",
"
Flight
\r\n",
"
Date
\r\n",
"
Count
\r\n",
"
Temperature
\r\n",
"
Pressure
\r\n",
"
Malfunction
\r\n",
"
Freq_Malfunction
\r\n",
"
\r\n",
"\r\n",
"\r\n",
"
\r\n",
"
1
\r\n",
"
04/12/81
\r\n",
"
6
\r\n",
"
66
\r\n",
"
50
\r\n",
"
0
\r\n",
"
0.00000
\r\n",
"
\r\n",
"
\r\n",
"
2
\r\n",
"
11/12/81
\r\n",
"
6
\r\n",
"
70
\r\n",
"
50
\r\n",
"
1
\r\n",
"
0.16667
\r\n",
"
\r\n",
"
\r\n",
"
3
\r\n",
"
03/22/82
\r\n",
"
6
\r\n",
"
69
\r\n",
"
50
\r\n",
"
0
\r\n",
"
0.00000
\r\n",
"
\r\n",
"
\r\n",
"
5
\r\n",
"
11/11/82
\r\n",
"
6
\r\n",
"
68
\r\n",
"
50
\r\n",
"
0
\r\n",
"
0.00000
\r\n",
"
\r\n",
"
\r\n",
"
6
\r\n",
"
04/04/83
\r\n",
"
6
\r\n",
"
67
\r\n",
"
50
\r\n",
"
0
\r\n",
"
0.00000
\r\n",
"
\r\n",
"
\r\n",
"
7
\r\n",
"
06/18/82
\r\n",
"
6
\r\n",
"
72
\r\n",
"
50
\r\n",
"
0
\r\n",
"
0.00000
\r\n",
"
\r\n",
"
\r\n",
"
8
\r\n",
"
08/30/83
\r\n",
"
6
\r\n",
"
73
\r\n",
"
100
\r\n",
"
0
\r\n",
"
0.00000
\r\n",
"
\r\n",
"
\r\n",
"
9
\r\n",
"
11/28/83
\r\n",
"
6
\r\n",
"
70
\r\n",
"
100
\r\n",
"
0
\r\n",
"
0.00000
\r\n",
"
\r\n",
"
\r\n",
"
41B
\r\n",
"
02/03/84
\r\n",
"
6
\r\n",
"
57
\r\n",
"
200
\r\n",
"
1
\r\n",
"
0.16667
\r\n",
"
\r\n",
"
\r\n",
"
41C
\r\n",
"
04/06/84
\r\n",
"
6
\r\n",
"
63
\r\n",
"
200
\r\n",
"
1
\r\n",
"
0.16667
\r\n",
"
\r\n",
"
\r\n",
"
41D
\r\n",
"
08/30/84
\r\n",
"
6
\r\n",
"
70
\r\n",
"
200
\r\n",
"
1
\r\n",
"
0.16667
\r\n",
"
\r\n",
"
\r\n",
"
41G
\r\n",
"
10/05/84
\r\n",
"
6
\r\n",
"
78
\r\n",
"
200
\r\n",
"
0
\r\n",
"
0.00000
\r\n",
"
\r\n",
"
\r\n",
"
51A
\r\n",
"
11/08/84
\r\n",
"
6
\r\n",
"
67
\r\n",
"
200
\r\n",
"
0
\r\n",
"
0.00000
\r\n",
"
\r\n",
"
\r\n",
"
51C
\r\n",
"
01/24/85
\r\n",
"
6
\r\n",
"
53
\r\n",
"
200
\r\n",
"
2
\r\n",
"
0.33333
\r\n",
"
\r\n",
"
\r\n",
"
51D
\r\n",
"
04/12/85
\r\n",
"
6
\r\n",
"
67
\r\n",
"
200
\r\n",
"
0
\r\n",
"
0.00000
\r\n",
"
\r\n",
"
\r\n",
"
51B
\r\n",
"
04/29/85
\r\n",
"
6
\r\n",
"
75
\r\n",
"
200
\r\n",
"
0
\r\n",
"
0.00000
\r\n",
"
\r\n",
"
\r\n",
"
51G
\r\n",
"
06/17/85
\r\n",
"
6
\r\n",
"
70
\r\n",
"
200
\r\n",
"
0
\r\n",
"
0.00000
\r\n",
"
\r\n",
"
\r\n",
"
51F
\r\n",
"
07/29/85
\r\n",
"
6
\r\n",
"
81
\r\n",
"
200
\r\n",
"
0
\r\n",
"
0.00000
\r\n",
"
\r\n",
"
\r\n",
"
51I
\r\n",
"
08/27/85
\r\n",
"
6
\r\n",
"
76
\r\n",
"
200
\r\n",
"
0
\r\n",
"
0.00000
\r\n",
"
\r\n",
"
\r\n",
"
51J
\r\n",
"
10/03/85
\r\n",
"
6
\r\n",
"
79
\r\n",
"
200
\r\n",
"
0
\r\n",
"
0.00000
\r\n",
"
\r\n",
"
\r\n",
"
61A
\r\n",
"
10/30/85
\r\n",
"
6
\r\n",
"
75
\r\n",
"
200
\r\n",
"
2
\r\n",
"
0.33333
\r\n",
"
\r\n",
"
\r\n",
"
61B
\r\n",
"
11/26/85
\r\n",
"
6
\r\n",
"
76
\r\n",
"
200
\r\n",
"
0
\r\n",
"
0.00000
\r\n",
"
\r\n",
"
\r\n",
"
61C
\r\n",
"
01/12/86
\r\n",
"
6
\r\n",
"
58
\r\n",
"
200
\r\n",
"
1
\r\n",
"
0.16667
\r\n",
"
\r\n",
"
\r\n",
"
\r\n",
"
.
\r\n",
"
.
\r\n",
"
31
\r\n",
"
.
\r\n",
"
.
\r\n",
"
.
\r\n",
"
\r\n",
"
\r\n",
"
\r\n",
"
.
\r\n",
"
.
\r\n",
"
30
\r\n",
"
.
\r\n",
"
.
\r\n",
"
.
\r\n",
"
\r\n",
"
\r\n",
"
\r\n",
"
.
\r\n",
"
.
\r\n",
"
35
\r\n",
"
.
\r\n",
"
.
\r\n",
"
.
\r\n",
"
\r\n",
"
\r\n",
"
\r\n",
"
.
\r\n",
"
.
\r\n",
"
40
\r\n",
"
.
\r\n",
"
.
\r\n",
"
.
\r\n",
"
\r\n",
"
\r\n",
"
\r\n",
"
.
\r\n",
"
.
\r\n",
"
45
\r\n",
"
.
\r\n",
"
.
\r\n",
"
.
\r\n",
"
\r\n",
"
\r\n",
"
\r\n",
"
.
\r\n",
"
.
\r\n",
"
50
\r\n",
"
.
\r\n",
"
.
\r\n",
"
.
\r\n",
"
\r\n",
"
\r\n",
"
\r\n",
"
.
\r\n",
"
.
\r\n",
"
55
\r\n",
"
.
\r\n",
"
.
\r\n",
"
.
\r\n",
"
\r\n",
"
\r\n",
"
\r\n",
"
.
\r\n",
"
.
\r\n",
"
60
\r\n",
"
.
\r\n",
"
.
\r\n",
"
.
\r\n",
"
\r\n",
"
\r\n",
"
\r\n",
"
.
\r\n",
"
.
\r\n",
"
65
\r\n",
"
.
\r\n",
"
.
\r\n",
"
.
\r\n",
"
\r\n",
"
\r\n",
"
\r\n",
"
.
\r\n",
"
.
\r\n",
"
70
\r\n",
"
.
\r\n",
"
.
\r\n",
"
.
\r\n",
"
\r\n",
"
\r\n",
"
\r\n",
"
.
\r\n",
"
.
\r\n",
"
75
\r\n",
"
.
\r\n",
"
.
\r\n",
"
.
\r\n",
"
\r\n",
"
\r\n",
"
\r\n",
"
.
\r\n",
"
.
\r\n",
"
80
\r\n",
"
.
\r\n",
"
.
\r\n",
"
.
\r\n",
"
\r\n",
"\r\n",
"
\r\n",
"\r\n",
"\r\n",
"\r\n",
"\r\n"
],
"text/plain": [
""
]
},
"execution_count": 3,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"data temp;\n",
" input Temperature @@;\n",
" cards;\n",
"31 30 35 40 45 50 55 60 65 70 75 80\n",
";\n",
"run;\n",
"\n",
"data nasa2;\n",
" set nasa temp;\n",
"run;\n",
"\n",
"proc print data=nasa2;\n",
" id Flight;\n",
"run;"
]
},
{
"cell_type": "markdown",
"metadata": {
"hideCode": true,
"hidePrompt": true
},
"source": [
"In the PROC LOGISTIC step, we use the events/trials syntax to indicate the number of failures and number of trials. The observations in the extrapolation sample are not used in fitting the model, yet the procedure produces predicted probabilities and logits (as long as the independent variable(s) are non-missing)."
]
},
{
"cell_type": "code",
"execution_count": 4,
"metadata": {},
"outputs": [
{
"data": {
"text/html": [
"\r\n",
"\r\n",
"\r\n",
"\r\n",
"\r\n",
"Sortie SAS\r\n",
"\r\n",
"\r\n",
"\r\n",
"
Sortie SAS
\r\n",
"\r\n",
"
\r\n",
"
NASA Space Shuttle O-Ring Failures
\r\n",
"
\r\n",
"
La LOGISTIC Procédure
\r\n",
"\r\n",
"
Informations sur le modèle
\r\n",
"
\r\n",
"
\r\n",
"\r\n",
"
\r\n",
"
Informations sur le modèle
\r\n",
"
\r\n",
"\r\n",
"\r\n",
"
\r\n",
"
Table
\r\n",
"
WORK.NASA2
\r\n",
"
\r\n",
"
\r\n",
"
\r\n",
"
Variable de réponse (Evénements)
\r\n",
"
Malfunction
\r\n",
"
Number of O-ring failures
\r\n",
"
\r\n",
"
\r\n",
"
Variable de réponse (Expériences)
\r\n",
"
Count
\r\n",
"
\r\n",
"
\r\n",
"
\r\n",
"
Modèle
\r\n",
"
logit binaire
\r\n",
"
\r\n",
"
\r\n",
"
\r\n",
"
Technique d'optimisation
\r\n",
"
Score de Fisher
\r\n",
"
\r\n",
"
\r\n",
"\r\n",
"
\r\n",
"\r\n",
"\r\n",
"
Synthèse des observations
\r\n",
"
\r\n",
"
\r\n",
"\r\n",
"
\r\n",
"
Nb d'observations lues
\r\n",
"
35
\r\n",
"
\r\n",
"
\r\n",
"
Nb d'observations utilisées
\r\n",
"
23
\r\n",
"
\r\n",
"
\r\n",
"
Somme des fréquences lues
\r\n",
"
138
\r\n",
"
\r\n",
"
\r\n",
"
Somme des fréquences utilis
\r\n",
"
138
\r\n",
"
\r\n",
"\r\n",
"
\r\n",
"\r\n",
"\r\n",
"
Profil de réponse
\r\n",
"
\r\n",
"
\r\n",
"\r\n",
"
\r\n",
"
Profil de réponse
\r\n",
"
\r\n",
"
\r\n",
"
Valeur ordonnée
\r\n",
"
Résultat binaire
\r\n",
"
Fréquence totale
\r\n",
"
\r\n",
"\r\n",
"\r\n",
"
\r\n",
"
1
\r\n",
"
Evénemen
\r\n",
"
9
\r\n",
"
\r\n",
"
\r\n",
"
2
\r\n",
"
Non-évén
\r\n",
"
129
\r\n",
"
\r\n",
"\r\n",
"
\r\n",
"\r\n",
"
\r\n",
"
Note:12 observations were deleted due to missing values for the response or explanatory variables.
\r\n",
"
\r\n",
"\r\n",
"
Etat de convergence
\r\n",
"
\r\n",
"
\r\n",
"\r\n",
"
\r\n",
"
Etat de convergence du modèle
\r\n",
"
\r\n",
"\r\n",
"\r\n",
"
\r\n",
"
Critère de convergence (GCONV=1E-8) respecté.
\r\n",
"
\r\n",
"\r\n",
"
\r\n",
"\r\n",
"\r\n",
"
Qualité de l'ajustement
\r\n",
"
\r\n",
"
\r\n",
"\r\n",
"
\r\n",
"
Statistique d'adéquation de la déviance et de Pearson
\r\n",
"
\r\n",
"
\r\n",
"
Critère
\r\n",
"
Valeur
\r\n",
"
DDL
\r\n",
"
Valeur/DDL
\r\n",
"
Pr > Khi-2
\r\n",
"
\r\n",
"\r\n",
"\r\n",
"
\r\n",
"
Ecart
\r\n",
"
18.0863
\r\n",
"
21
\r\n",
"
0.8613
\r\n",
"
0.6435
\r\n",
"
\r\n",
"
\r\n",
"
Pearson
\r\n",
"
29.9803
\r\n",
"
21
\r\n",
"
1.4276
\r\n",
"
0.0924
\r\n",
"
\r\n",
"\r\n",
"
\r\n",
"
\r\n",
"
Nombre d'observations d'événements/expériences : 23
\r\n",
"
\r\n",
"\r\n",
"\r\n",
"
Statistiques d'ajustement
\r\n",
"
\r\n",
"
\r\n",
"\r\n",
"
\r\n",
"
Statistiques d'ajustement du modèle
\r\n",
"
\r\n",
"
\r\n",
"
Critère
\r\n",
"
Constante uniquement
\r\n",
"
Constante et covariables
\r\n",
"
\r\n",
"
\r\n",
"
Log-vraisemblance
\r\n",
"
Log-vraisemblance complète
\r\n",
"
\r\n",
"\r\n",
"\r\n",
"
\r\n",
"
AIC
\r\n",
"
68.540
\r\n",
"
64.396
\r\n",
"
35.647
\r\n",
"
\r\n",
"
\r\n",
"
SC
\r\n",
"
71.468
\r\n",
"
70.251
\r\n",
"
41.501
\r\n",
"
\r\n",
"
\r\n",
"
-2 Log
\r\n",
"
66.540
\r\n",
"
60.396
\r\n",
"
31.647
\r\n",
"
\r\n",
"\r\n",
"
\r\n",
"\r\n",
"\r\n",
"
Tests globaux
\r\n",
"
\r\n",
"
\r\n",
"\r\n",
"
\r\n",
"
Test de l'hypothèse nulle globale : BETA=0
\r\n",
"
\r\n",
"
\r\n",
"
Test
\r\n",
"
Khi-2
\r\n",
"
DDL
\r\n",
"
Pr > Khi-2
\r\n",
"
\r\n",
"\r\n",
"\r\n",
"
\r\n",
"
Rapport de vrais
\r\n",
"
6.1440
\r\n",
"
1
\r\n",
"
0.0132
\r\n",
"
\r\n",
"
\r\n",
"
Score
\r\n",
"
6.7696
\r\n",
"
1
\r\n",
"
0.0093
\r\n",
"
\r\n",
"
\r\n",
"
Wald
\r\n",
"
6.0435
\r\n",
"
1
\r\n",
"
0.0140
\r\n",
"
\r\n",
"\r\n",
"
\r\n",
"\r\n",
"\r\n",
"
Valeurs estimées du paramètre
\r\n",
"
\r\n",
"
\r\n",
"\r\n",
"
\r\n",
"
Estimations par l'analyse du maximum de vraisemblance
\r\n",
"
\r\n",
"
\r\n",
"
Paramètre
\r\n",
"
DDL
\r\n",
"
Estimation
\r\n",
"
Erreur type
\r\n",
"
Khi-2 de Wald
\r\n",
"
Pr > Khi-2
\r\n",
"
\r\n",
"\r\n",
"\r\n",
"
\r\n",
"
Intercept
\r\n",
"
1
\r\n",
"
5.0850
\r\n",
"
3.0525
\r\n",
"
2.7751
\r\n",
"
0.0957
\r\n",
"
\r\n",
"
\r\n",
"
Temperature
\r\n",
"
1
\r\n",
"
-0.1156
\r\n",
"
0.0470
\r\n",
"
6.0435
\r\n",
"
0.0140
\r\n",
"
\r\n",
"\r\n",
"
\r\n",
"\r\n",
"\r\n",
"
Rapports de cotes
\r\n",
"
\r\n",
"
\r\n",
"\r\n",
"
\r\n",
"
Estimations des rapports de cotes
\r\n",
"
\r\n",
"
\r\n",
"
Effet
\r\n",
"
Valeur estimée du point
\r\n",
"
95% Intervalle de confiance de Wald
\r\n",
"
\r\n",
"\r\n",
"\r\n",
"
\r\n",
"
Temperature
\r\n",
"
0.891
\r\n",
"
0.812
\r\n",
"
0.977
\r\n",
"
\r\n",
"\r\n",
"
\r\n",
"\r\n",
"\r\n",
"
Statistiques d'association
\r\n",
"
\r\n",
"
\r\n",
"\r\n",
"
\r\n",
"
Association des probabilités prédites et des réponses observées
\r\n",
"
\r\n",
"\r\n",
"\r\n",
"
\r\n",
"
Pourcentage concordant
\r\n",
"
65.4
\r\n",
"
D de Somers
\r\n",
"
0.382
\r\n",
"
\r\n",
"
\r\n",
"
Pourcentage discordant
\r\n",
"
27.1
\r\n",
"
Gamma
\r\n",
"
0.413
\r\n",
"
\r\n",
"
\r\n",
"
Pourcentage lié
\r\n",
"
7.5
\r\n",
"
Tau-a
\r\n",
"
0.047
\r\n",
"
\r\n",
"
\r\n",
"
Paires
\r\n",
"
1161
\r\n",
"
c
\r\n",
"
0.691
\r\n",
"
\r\n",
"\r\n",
"
\r\n",
"\r\n",
"\r\n",
"\r\n",
"\r\n"
],
"text/plain": [
""
]
},
"execution_count": 4,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"ods output OddsRatios=OddsRatios ParameterEstimates=ParameterEstimates GoodnessOfFit=GoodnessOfFit;\n",
"proc logistic data=nasa2 nosimple;\n",
" model Malfunction/Count = Temperature / scale=none;\n",
" output out=results p=predict l=lower u=upper;\n",
"run;"
]
},
{
"cell_type": "markdown",
"metadata": {
"hideCode": true,
"hidePrompt": true
},
"source": [
"The printed output indicates that the 12 new observations were not used in the analysis."
]
},
{
"cell_type": "code",
"execution_count": 5,
"metadata": {},
"outputs": [
{
"data": {
"text/html": [
"\r\n",
"\r\n",
"\r\n",
"\r\n",
"\r\n",
"Sortie SAS\r\n",
"\r\n",
"\r\n",
"\r\n",
"
Sortie SAS
\r\n",
"\r\n",
"
\r\n",
"
Odds Ratios
\r\n",
"
\r\n",
"
La PRINT Procédure
\r\n",
"\r\n",
"
Table WORK.ODDSRATIOS
\r\n",
"
\r\n",
"
\r\n",
"\r\n",
"
\r\n",
"
Effet
\r\n",
"
Estimation du rapport de cotes
\r\n",
"
Borne inférieure de l'IC à 95% pour le rapport de cotes
\r\n",
"
Borne supérieure de l'IC à 95% pour le rapport de cotes
\r\n",
"
\r\n",
"\r\n",
"\r\n",
"
\r\n",
"
Temperature
\r\n",
"
0.891
\r\n",
"
0.812
\r\n",
"
0.977
\r\n",
"
\r\n",
"\r\n",
"
\r\n",
"\r\n",
"\r\n",
"\r\n",
"\r\n"
],
"text/plain": [
""
]
},
"execution_count": 5,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"proc print data=OddsRatios label;\n",
" id Effect;\n",
" title \"Odds Ratios\";\n",
"run;"
]
},
{
"cell_type": "markdown",
"metadata": {
"hideCode": true,
"hidePrompt": true
},
"source": [
"The odds ratio, 0.891, is interpreted to mean that each increase of 1° in temperature decreases the odds of a failure by 11%!"
]
},
{
"cell_type": "code",
"execution_count": 6,
"metadata": {},
"outputs": [
{
"data": {
"text/html": [
"\r\n",
"\r\n",
"\r\n",
"\r\n",
"\r\n",
"Sortie SAS\r\n",
"\r\n",
"\r\n",
"\r\n",
"
Sortie SAS
\r\n",
"\r\n",
"
\r\n",
"
Parameter Estimates
\r\n",
"
\r\n",
"
La PRINT Procédure
\r\n",
"\r\n",
"
Table WORK.PARAMETERESTIMATES
\r\n",
"
\r\n",
"
\r\n",
"\r\n",
"
\r\n",
"
Variable
\r\n",
"
DDL
\r\n",
"
Estimation
\r\n",
"
Erreur type
\r\n",
"
Khi-2 de Wald
\r\n",
"
Pr > Khi-2
\r\n",
"
\r\n",
"\r\n",
"\r\n",
"
\r\n",
"
Intercept
\r\n",
"
1
\r\n",
"
5.08498
\r\n",
"
3.05249
\r\n",
"
2.7751
\r\n",
"
0.0957
\r\n",
"
\r\n",
"
\r\n",
"
Temperature
\r\n",
"
1
\r\n",
"
-0.11560
\r\n",
"
0.04702
\r\n",
"
6.0435
\r\n",
"
0.0140
\r\n",
"
\r\n",
"\r\n",
"
\r\n",
"\r\n",
"\r\n",
"\r\n",
"\r\n"
],
"text/plain": [
""
]
},
"execution_count": 6,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"proc print data=ParameterEstimates label;\n",
" id variable;\n",
" format estimate stderr 8.5;\n",
" var _numeric_;\n",
" title \"Parameter Estimates\";\n",
"run;"
]
},
{
"cell_type": "markdown",
"metadata": {
"hideCode": true,
"hidePrompt": true
},
"source": [
"The maximum likelyhood estimator of the intercept and of Temperature are thus $\\hat{\\alpha}$ = **5.0850** and $\\hat{\\beta}$ = **-0.1156** and their standard errors are $s_{\\hat{\\alpha}}$ = **3.052** and $s_{\\hat{\\beta}}$ = **0.04702**."
]
},
{
"cell_type": "code",
"execution_count": 7,
"metadata": {},
"outputs": [
{
"data": {
"text/html": [
"\r\n",
"\r\n",
"\r\n",
"\r\n",
"\r\n",
"Sortie SAS\r\n",
"\r\n",
"\r\n",
"\r\n",
"
Sortie SAS
\r\n",
"\r\n",
"
\r\n",
"
Goodness Of Fit
\r\n",
"
\r\n",
"
La PRINT Procédure
\r\n",
"\r\n",
"
Table WORK.GOODNESSOFFIT
\r\n",
"
\r\n",
"
\r\n",
"\r\n",
"
\r\n",
"
Critère
\r\n",
"
DDL
\r\n",
"
Khi-2
\r\n",
"
Khi-2/DDL
\r\n",
"
Pr > Khi-2
\r\n",
"
\r\n",
"\r\n",
"\r\n",
"
\r\n",
"
Ecart
\r\n",
"
21
\r\n",
"
18.0863
\r\n",
"
0.8613
\r\n",
"
0.6435
\r\n",
"
\r\n",
"
\r\n",
"
Pearson
\r\n",
"
21
\r\n",
"
29.9803
\r\n",
"
1.4276
\r\n",
"
0.0924
\r\n",
"
\r\n",
"\r\n",
"
\r\n",
"\r\n",
"\r\n",
"\r\n",
"\r\n"
],
"text/plain": [
""
]
},
"execution_count": 7,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"proc print data=GoodnessOfFit label;\n",
" id criterion;\n",
" title \"Goodness Of Fit\";\n",
"run;"
]
},
{
"cell_type": "markdown",
"metadata": {
"hideCode": true,
"hidePrompt": true
},
"source": [
"The Residual deviance corresponds to the Goodness of fit $G^2$ = **18.086** with **21** degrees of freedom."
]
},
{
"cell_type": "markdown",
"metadata": {
"hideCode": true,
"hidePrompt": true
},
"source": [
"The output data set results contains the predicted probability of a failure at each temperature and upper and\n",
"lower confidence 95% limits for this probability."
]
},
{
"cell_type": "code",
"execution_count": 8,
"metadata": {},
"outputs": [
{
"data": {
"text/html": [
"\r\n",
"\r\n",
"\r\n",
"\r\n",
"\r\n",
"Sortie SAS\r\n",
"\r\n",
"\r\n",
"\r\n",
"
Sortie SAS
\r\n",
"\r\n",
"
\r\n",
"
Predicted probabilities of a failure
\r\n",
"
\r\n",
"
La PRINT Procédure
\r\n",
"\r\n",
"
Table WORK.RESULTS
\r\n",
"
\r\n",
"
\r\n",
"\r\n",
"
\r\n",
"
Flight
\r\n",
"
Date
\r\n",
"
Count
\r\n",
"
Temperature
\r\n",
"
Pressure
\r\n",
"
Malfunction
\r\n",
"
Freq_Malfunction
\r\n",
"
predict
\r\n",
"
lower
\r\n",
"
upper
\r\n",
"
\r\n",
"\r\n",
"\r\n",
"
\r\n",
"
1
\r\n",
"
04/12/81
\r\n",
"
6
\r\n",
"
66
\r\n",
"
50
\r\n",
"
0
\r\n",
"
0.00000
\r\n",
"
0.07278
\r\n",
"
0.03716
\r\n",
"
0.13767
\r\n",
"
\r\n",
"
\r\n",
"
2
\r\n",
"
11/12/81
\r\n",
"
6
\r\n",
"
70
\r\n",
"
50
\r\n",
"
1
\r\n",
"
0.16667
\r\n",
"
0.04711
\r\n",
"
0.02044
\r\n",
"
0.10482
\r\n",
"
\r\n",
"
\r\n",
"
3
\r\n",
"
03/22/82
\r\n",
"
6
\r\n",
"
69
\r\n",
"
50
\r\n",
"
0
\r\n",
"
0.00000
\r\n",
"
0.05258
\r\n",
"
0.02406
\r\n",
"
0.11104
\r\n",
"
\r\n",
"
\r\n",
"
5
\r\n",
"
11/11/82
\r\n",
"
6
\r\n",
"
68
\r\n",
"
50
\r\n",
"
0
\r\n",
"
0.00000
\r\n",
"
0.05864
\r\n",
"
0.02809
\r\n",
"
0.11837
\r\n",
"
\r\n",
"
\r\n",
"
6
\r\n",
"
04/04/83
\r\n",
"
6
\r\n",
"
67
\r\n",
"
50
\r\n",
"
0
\r\n",
"
0.00000
\r\n",
"
0.06536
\r\n",
"
0.03248
\r\n",
"
0.12713
\r\n",
"
\r\n",
"
\r\n",
"
7
\r\n",
"
06/18/82
\r\n",
"
6
\r\n",
"
72
\r\n",
"
50
\r\n",
"
0
\r\n",
"
0.00000
\r\n",
"
0.03775
\r\n",
"
0.01447
\r\n",
"
0.09485
\r\n",
"
\r\n",
"
\r\n",
"
8
\r\n",
"
08/30/83
\r\n",
"
6
\r\n",
"
73
\r\n",
"
100
\r\n",
"
0
\r\n",
"
0.00000
\r\n",
"
0.03377
\r\n",
"
0.01209
\r\n",
"
0.09077
\r\n",
"
\r\n",
"
\r\n",
"
9
\r\n",
"
11/28/83
\r\n",
"
6
\r\n",
"
70
\r\n",
"
100
\r\n",
"
0
\r\n",
"
0.00000
\r\n",
"
0.04711
\r\n",
"
0.02044
\r\n",
"
0.10482
\r\n",
"
\r\n",
"
\r\n",
"
41B
\r\n",
"
02/03/84
\r\n",
"
6
\r\n",
"
57
\r\n",
"
200
\r\n",
"
1
\r\n",
"
0.16667
\r\n",
"
0.18179
\r\n",
"
0.07703
\r\n",
"
0.37163
\r\n",
"
\r\n",
"
\r\n",
"
41C
\r\n",
"
04/06/84
\r\n",
"
6
\r\n",
"
63
\r\n",
"
200
\r\n",
"
1
\r\n",
"
0.16667
\r\n",
"
0.09994
\r\n",
"
0.05182
\r\n",
"
0.18408
\r\n",
"
\r\n",
"
\r\n",
"
41D
\r\n",
"
08/30/84
\r\n",
"
6
\r\n",
"
70
\r\n",
"
200
\r\n",
"
1
\r\n",
"
0.16667
\r\n",
"
0.04711
\r\n",
"
0.02044
\r\n",
"
0.10482
\r\n",
"
\r\n",
"
\r\n",
"
41G
\r\n",
"
10/05/84
\r\n",
"
6
\r\n",
"
78
\r\n",
"
200
\r\n",
"
0
\r\n",
"
0.00000
\r\n",
"
0.01923
\r\n",
"
0.00468
\r\n",
"
0.07557
\r\n",
"
\r\n",
"
\r\n",
"
51A
\r\n",
"
11/08/84
\r\n",
"
6
\r\n",
"
67
\r\n",
"
200
\r\n",
"
0
\r\n",
"
0.00000
\r\n",
"
0.06536
\r\n",
"
0.03248
\r\n",
"
0.12713
\r\n",
"
\r\n",
"
\r\n",
"
51C
\r\n",
"
01/24/85
\r\n",
"
6
\r\n",
"
53
\r\n",
"
200
\r\n",
"
2
\r\n",
"
0.33333
\r\n",
"
0.26079
\r\n",
"
0.09049
\r\n",
"
0.55575
\r\n",
"
\r\n",
"
\r\n",
"
51D
\r\n",
"
04/12/85
\r\n",
"
6
\r\n",
"
67
\r\n",
"
200
\r\n",
"
0
\r\n",
"
0.00000
\r\n",
"
0.06536
\r\n",
"
0.03248
\r\n",
"
0.12713
\r\n",
"
\r\n",
"
\r\n",
"
51B
\r\n",
"
04/29/85
\r\n",
"
6
\r\n",
"
75
\r\n",
"
200
\r\n",
"
0
\r\n",
"
0.00000
\r\n",
"
0.02699
\r\n",
"
0.00833
\r\n",
"
0.08385
\r\n",
"
\r\n",
"
\r\n",
"
51G
\r\n",
"
06/17/85
\r\n",
"
6
\r\n",
"
70
\r\n",
"
200
\r\n",
"
0
\r\n",
"
0.00000
\r\n",
"
0.04711
\r\n",
"
0.02044
\r\n",
"
0.10482
\r\n",
"
\r\n",
"
\r\n",
"
51F
\r\n",
"
07/29/85
\r\n",
"
6
\r\n",
"
81
\r\n",
"
200
\r\n",
"
0
\r\n",
"
0.00000
\r\n",
"
0.01367
\r\n",
"
0.00259
\r\n",
"
0.06887
\r\n",
"
\r\n",
"
\r\n",
"
51I
\r\n",
"
08/27/85
\r\n",
"
6
\r\n",
"
76
\r\n",
"
200
\r\n",
"
0
\r\n",
"
0.00000
\r\n",
"
0.02411
\r\n",
"
0.00689
\r\n",
"
0.08086
\r\n",
"
\r\n",
"
\r\n",
"
51J
\r\n",
"
10/03/85
\r\n",
"
6
\r\n",
"
79
\r\n",
"
200
\r\n",
"
0
\r\n",
"
0.00000
\r\n",
"
0.01717
\r\n",
"
0.00385
\r\n",
"
0.07319
\r\n",
"
\r\n",
"
\r\n",
"
61A
\r\n",
"
10/30/85
\r\n",
"
6
\r\n",
"
75
\r\n",
"
200
\r\n",
"
2
\r\n",
"
0.33333
\r\n",
"
0.02699
\r\n",
"
0.00833
\r\n",
"
0.08385
\r\n",
"
\r\n",
"
\r\n",
"
61B
\r\n",
"
11/26/85
\r\n",
"
6
\r\n",
"
76
\r\n",
"
200
\r\n",
"
0
\r\n",
"
0.00000
\r\n",
"
0.02411
\r\n",
"
0.00689
\r\n",
"
0.08086
\r\n",
"
\r\n",
"
\r\n",
"
61C
\r\n",
"
01/12/86
\r\n",
"
6
\r\n",
"
58
\r\n",
"
200
\r\n",
"
1
\r\n",
"
0.16667
\r\n",
"
0.16522
\r\n",
"
0.07334
\r\n",
"
0.33107
\r\n",
"
\r\n",
"
\r\n",
"
\r\n",
"
.
\r\n",
"
.
\r\n",
"
31
\r\n",
"
.
\r\n",
"
.
\r\n",
"
.
\r\n",
"
0.81777
\r\n",
"
0.15960
\r\n",
"
0.99066
\r\n",
"
\r\n",
"
\r\n",
"
\r\n",
"
.
\r\n",
"
.
\r\n",
"
30
\r\n",
"
.
\r\n",
"
.
\r\n",
"
.
\r\n",
"
0.83437
\r\n",
"
0.16307
\r\n",
"
0.99238
\r\n",
"
\r\n",
"
\r\n",
"
\r\n",
"
.
\r\n",
"
.
\r\n",
"
35
\r\n",
"
.
\r\n",
"
.
\r\n",
"
.
\r\n",
"
0.73864
\r\n",
"
0.14615
\r\n",
"
0.97902
\r\n",
"
\r\n",
"
\r\n",
"
\r\n",
"
.
\r\n",
"
.
\r\n",
"
40
\r\n",
"
.
\r\n",
"
.
\r\n",
"
.
\r\n",
"
0.61323
\r\n",
"
0.13016
\r\n",
"
0.94382
\r\n",
"
\r\n",
"
\r\n",
"
\r\n",
"
.
\r\n",
"
.
\r\n",
"
45
\r\n",
"
.
\r\n",
"
.
\r\n",
"
.
\r\n",
"
0.47076
\r\n",
"
0.11487
\r\n",
"
0.85910
\r\n",
"
\r\n",
"
\r\n",
"
\r\n",
"
.
\r\n",
"
.
\r\n",
"
50
\r\n",
"
.
\r\n",
"
.
\r\n",
"
.
\r\n",
"
0.33290
\r\n",
"
0.09979
\r\n",
"
0.69198
\r\n",
"
\r\n",
"
\r\n",
"
\r\n",
"
.
\r\n",
"
.
\r\n",
"
55
\r\n",
"
.
\r\n",
"
.
\r\n",
"
.
\r\n",
"
0.21873
\r\n",
"
0.08398
\r\n",
"
0.46092
\r\n",
"
\r\n",
"
\r\n",
"
\r\n",
"
.
\r\n",
"
.
\r\n",
"
60
\r\n",
"
.
\r\n",
"
.
\r\n",
"
.
\r\n",
"
0.13574
\r\n",
"
0.06538
\r\n",
"
0.26070
\r\n",
"
\r\n",
"
\r\n",
"
\r\n",
"
.
\r\n",
"
.
\r\n",
"
65
\r\n",
"
.
\r\n",
"
.
\r\n",
"
.
\r\n",
"
0.08098
\r\n",
"
0.04202
\r\n",
"
0.15038
\r\n",
"
\r\n",
"
\r\n",
"
\r\n",
"
.
\r\n",
"
.
\r\n",
"
70
\r\n",
"
.
\r\n",
"
.
\r\n",
"
.
\r\n",
"
0.04711
\r\n",
"
0.02044
\r\n",
"
0.10482
\r\n",
"
\r\n",
"
\r\n",
"
\r\n",
"
.
\r\n",
"
.
\r\n",
"
75
\r\n",
"
.
\r\n",
"
.
\r\n",
"
.
\r\n",
"
0.02699
\r\n",
"
0.00833
\r\n",
"
0.08385
\r\n",
"
\r\n",
"
\r\n",
"
\r\n",
"
.
\r\n",
"
.
\r\n",
"
80
\r\n",
"
.
\r\n",
"
.
\r\n",
"
.
\r\n",
"
0.01532
\r\n",
"
0.00316
\r\n",
"
0.07097
\r\n",
"
\r\n",
"\r\n",
"
\r\n",
"\r\n",
"\r\n",
"\r\n",
"\r\n"
],
"text/plain": [
""
]
},
"execution_count": 8,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"proc print data=results;\n",
" id Flight;\n",
" title \"Predicted probabilities of a failure\";\n",
"run;"
]
},
{
"cell_type": "markdown",
"metadata": {
"hideCode": true,
"hidePrompt": true
},
"source": [
"We can plot the predicted and observed values as shown below. A vertical reference line at 31° is used to highlight the conditions at the Challenger launch.\n",
"\n",
"The graph is shown in Figure below. There’s not much data at low temperatures (the confidence band is quite\n",
"wide), but the predicted probability of failure is uncomfortably high. Would you take a ride on Challenger when\n",
"the weather is cold?"
]
},
{
"cell_type": "code",
"execution_count": 9,
"metadata": {},
"outputs": [
{
"data": {
"text/html": [
"\r\n",
"\r\n",
"\r\n",
"\r\n",
"\r\n",
"Sortie SAS\r\n",
"\r\n",
"\r\n",
"\r\n",
"