The Pima dataset - predicting diabetes
library (MASS)head (Pima.te)
npreg glu bp skin bmi ped age type
1 6 148 72 35 33.6 0.627 50 Yes
2 1 85 66 29 26.6 0.351 31 No
3 1 89 66 23 28.1 0.167 21 No
4 3 78 50 32 31.0 0.248 26 Yes
5 2 197 70 45 30.5 0.158 53 Yes
6 5 166 72 19 25.8 0.587 51 Yes
$ glu_bin <- (round (Pima.te$ glu/ 10 ))* 10 $ typeTF <- as.numeric (Pima.te$ type== "Yes" )<- glm (typeTF ~ glu, data= Pima.te, family= "binomial" )<- aggregate (typeTF ~ glu_bin, data= Pima.te, FUN= mean)<- seq (60 ,200 ,10 )<- predict (Pima.glm, newdata= data.frame (glu= xs))<- exp (zs)/ (1 + exp (zs))plot (x= Pima.te$ glu, y= Pima.te$ typeTF, pch= 16 , col= rgb (0 ,0 ,1 ,.2 ), xlim= c (50 ,200 ))#points(aggr.props, ylim=c(0,1)) lines (x= xs,y= ps, col= "red" )text (x= xs-2 , y= as.vector (aggr.props$ typeTF* .8 +.5 * .2 ), label= round (aggr.props$ typeTF,2 ))segments (x0= xs-5 , aggr.props$ typeTF, xs+ 5 ,aggr.props$ typeTF)abline (v= seq (55 ,195 ,10 ), col= "gray" )
Model summary
Call:
glm(formula = typeTF ~ glu, family = "binomial", data = Pima.te)
Deviance Residuals:
Min 1Q Median 3Q Max
-2.2343 -0.7270 -0.4985 0.6663 2.3268
Coefficients:
Estimate Std. Error z value Pr(>|z|)
(Intercept) -5.946808 0.659839 -9.013 <2e-16 ***
glu 0.042421 0.005165 8.213 <2e-16 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
(Dispersion parameter for binomial family taken to be 1)
Null deviance: 420.30 on 331 degrees of freedom
Residual deviance: 325.99 on 330 degrees of freedom
AIC: 329.99
Number of Fisher Scoring iterations: 4
Odds Ratios
Odds Ratio = Odds(Event | Case A) / Odds(Event | Case B) We look at odds ratio when the cases differe by an increase of 1 unit in the predictor \(x\) . \(\beta_1\) is the log of the odds ratio associated with a 1 unit increase in \(x\) .
The log(odds ratio)=0.042421
The odds ratio = exp(0.042421)
Interpret the slope
exp (Pima.glm$ coefficients[2 ])
Odd ratio for a 1 unit increase in blood glucose is 1.0433 so if we hold all other variables constant and increase blood glucose by 1 unit, we expect an increase in the Odds by 4.33%
odds ratio vs risk ratio
Risk is Prob(Event | case A) / Pr(Event | Case B) which is not the same as odds ratio. Consider this example:
Predict probability for glu 200 vs 201
predict (Pima.glm, type= "response" , newdata= data.frame (glu= c (200 ,201 )))#The risk ratio is 0.9295508 / 0.9267217 #The Odds ratio is 0.9295508 / (1 -0.9295508 )) / (0.9267217 / (1 -0.9267217 ))
Compare glu 100 to 101
predict (Pima.glm, type= "response" , newdata= data.frame (glu= c (100 ,101 )))#The risk ratio is 0.1594550 / 0.1538511 #The Odds ratio is 0.1594550 / (1 -0.1594550 )) / (0.1538511 / (1 -0.1538511 ))
The risk ratios are different but the odds ratio is a constant 1.04333
Log Likelihood
log likelihood calculation Let’s compare the true Y values and the predicted Y values
<- Pima.te$ typeTF<- predict (Pima.glm, type= "response" )<- y * y.hat + (1 - y)* (1 - y.hat)data.frame (y,y.hat,pY)
y y.hat pY
1 1 0.58212363 0.58212363
2 0 0.08778183 0.91221817
3 0 0.10235379 0.89764621
4 1 0.06673426 0.06673426
5 1 0.91759616 0.91759616
6 1 0.74933615 0.74933615
7 1 0.28067169 0.28067169
8 0 0.17115736 0.82884264
9 0 0.35394014 0.64605986
10 1 0.28931541 0.28931541
11 0 0.13800342 0.86199658
12 0 0.21033272 0.78966728
13 0 0.09852148 0.90147852
14 0 0.31616589 0.68383411
15 0 0.17115736 0.82884264
16 1 0.16522309 0.16522309
17 1 0.10631759 0.10631759
18 1 0.22477052 0.22477052
19 0 0.84408847 0.15591153
20 0 0.18997325 0.81002675
21 1 0.78704085 0.78704085
22 1 0.84408847 0.84408847
23 0 0.05045417 0.94954583
24 0 0.17115736 0.82884264
25 0 0.15945498 0.84054502
26 0 0.09852148 0.90147852
27 0 0.60260691 0.39739309
28 0 0.05467737 0.94532263
29 0 0.18353106 0.81646894
30 0 0.14840944 0.85159056
31 1 0.21033272 0.21033272
32 0 0.12821718 0.87178282
33 0 0.56135307 0.43864693
34 1 0.15385113 0.15385113
35 0 0.38355137 0.61644863
36 0 0.12821718 0.87178282
37 0 0.23224851 0.76775149
38 0 0.23989830 0.76010170
39 0 0.08122139 0.91877861
40 0 0.15945498 0.84054502
41 0 0.18997325 0.81002675
42 0 0.15385113 0.84614887
43 1 0.45572761 0.45572761
44 0 0.32540821 0.67459179
45 0 0.07511087 0.92488913
46 0 0.11465216 0.88534784
47 0 0.11902908 0.88097092
48 1 0.31616589 0.31616589
49 0 0.07511087 0.92488913
50 0 0.35394014 0.64605986
51 0 0.54036528 0.45963472
52 0 0.10235379 0.89764621
53 1 0.69858083 0.69858083
54 0 0.13800342 0.86199658
55 1 0.71613927 0.71613927
56 0 0.19658713 0.80341287
57 1 0.09852148 0.09852148
58 1 0.27218732 0.27218732
59 1 0.80091481 0.80091481
60 1 0.77984422 0.77984422
61 0 0.08444380 0.91555620
62 0 0.15385113 0.84614887
63 0 0.11902908 0.88097092
64 0 0.18997325 0.81002675
65 0 0.20337345 0.79662655
66 0 0.18997325 0.81002675
67 0 0.10631759 0.89368241
68 1 0.66169683 0.66169683
69 0 0.63264996 0.36735004
70 1 0.62273686 0.62273686
71 0 0.09852148 0.90147852
72 1 0.72468315 0.72468315
73 0 0.61272001 0.38727999
74 0 0.16522309 0.83477691
75 0 0.24771876 0.75228124
76 1 0.17725952 0.17725952
77 0 0.05923201 0.94076799
78 1 0.83842391 0.83842391
79 1 0.38355137 0.38355137
80 1 0.37357216 0.37357216
81 1 0.21033272 0.21033272
82 1 0.21033272 0.21033272
83 0 0.23989830 0.76010170
84 0 0.20337345 0.79662655
85 0 0.22477052 0.77522948
86 1 0.91433121 0.91433121
87 0 0.07511087 0.92488913
88 0 0.57176996 0.42823004
89 1 0.21033272 0.21033272
90 0 0.34430105 0.65569895
91 1 0.08778183 0.08778183
92 0 0.51923334 0.48076666
93 0 0.15385113 0.84614887
94 0 0.09481750 0.90518250
95 1 0.71613927 0.71613927
96 0 0.91759616 0.08240384
97 0 0.27218732 0.72781268
98 1 0.43477468 0.43477468
99 0 0.05691189 0.94308811
100 1 0.84959008 0.84959008
101 1 0.83842391 0.83842391
102 0 0.11041602 0.88958398
103 0 0.11041602 0.88958398
104 1 0.28931541 0.28931541
105 0 0.56135307 0.43864693
106 1 0.86513981 0.86513981
107 0 0.74128419 0.25871581
108 0 0.33478843 0.66521157
109 0 0.22477052 0.77522948
110 0 0.10631759 0.89368241
111 1 0.23989830 0.23989830
112 0 0.22477052 0.77522948
113 1 0.65213658 0.65213658
114 0 0.12821718 0.87178282
115 0 0.13303414 0.86696586
116 0 0.37357216 0.62642784
117 1 0.15945498 0.15945498
118 0 0.20337345 0.79662655
119 0 0.15385113 0.84614887
120 1 0.18997325 0.18997325
121 0 0.17725952 0.82274048
122 0 0.20337345 0.79662655
123 0 0.42438023 0.57561977
124 0 0.45572761 0.54427239
125 1 0.28931541 0.28931541
126 0 0.13303414 0.86696586
127 0 0.06673426 0.93326574
128 1 0.19658713 0.19658713
129 0 0.61272001 0.38727999
130 1 0.56135307 0.56135307
131 0 0.35394014 0.64605986
132 1 0.54036528 0.54036528
133 0 0.29811500 0.70188500
134 1 0.70743729 0.70743729
135 1 0.37357216 0.37357216
136 1 0.33478843 0.33478843
137 1 0.65213658 0.65213658
138 0 0.23989830 0.76010170
139 1 0.21033272 0.21033272
140 0 0.25570840 0.74429160
141 1 0.62273686 0.62273686
142 0 0.23224851 0.76775149
143 0 0.67112685 0.32887315
144 1 0.31616589 0.31616589
145 1 0.16522309 0.16522309
146 0 0.18353106 0.81646894
147 0 0.28067169 0.71932831
148 0 0.09481750 0.90518250
149 0 0.12821718 0.87178282
150 0 0.39362965 0.60637035
151 0 0.12821718 0.87178282
152 0 0.35394014 0.64605986
153 0 0.48742974 0.51257026
154 0 0.14840944 0.85159056
155 1 0.34430105 0.34430105
156 1 0.91433121 0.91433121
157 1 0.88802830 0.88802830
158 0 0.17115736 0.82884264
159 0 0.57176996 0.42823004
160 0 0.14840944 0.85159056
161 0 0.07511087 0.92488913
162 1 0.80091481 0.80091481
163 0 0.08444380 0.91555620
164 0 0.14312768 0.85687232
165 0 0.09481750 0.90518250
166 0 0.11902908 0.88097092
167 0 0.18353106 0.81646894
168 0 0.10631759 0.89368241
169 0 0.34430105 0.65569895
170 0 0.28931541 0.71068459
171 1 0.26386534 0.26386534
172 1 0.18353106 0.18353106
173 0 0.15385113 0.84614887
174 0 0.40379928 0.59620072
175 0 0.36369940 0.63630060
176 0 0.13303414 0.86696586
177 1 0.90382283 0.90382283
178 1 0.45572761 0.45572761
179 0 0.05252569 0.94747431
180 1 0.79406435 0.79406435
181 0 0.16522309 0.83477691
182 0 0.23224851 0.76775149
183 0 0.52981267 0.47018733
184 1 0.80091481 0.80091481
185 1 0.54036528 0.54036528
186 1 0.38355137 0.38355137
187 0 0.28931541 0.71068459
188 0 0.12354978 0.87645022
189 0 0.16522309 0.83477691
190 1 0.61272001 0.61272001
191 1 0.84959008 0.84959008
192 1 0.12821718 0.12821718
193 0 0.10235379 0.89764621
194 0 0.07221652 0.92778348
195 0 0.10631759 0.89368241
196 1 0.88802830 0.88802830
197 0 0.27218732 0.72781268
198 1 0.84408847 0.84408847
199 1 0.17725952 0.17725952
200 0 0.29811500 0.70188500
201 0 0.07811146 0.92188854
202 0 0.11041602 0.88958398
203 0 0.43477468 0.56522532
204 0 0.29811500 0.70188500
205 0 0.05691189 0.94308811
206 0 0.09852148 0.90147852
207 0 0.33478843 0.66521157
208 0 0.13800342 0.86199658
209 0 0.54036528 0.45963472
210 0 0.46626792 0.53373208
211 0 0.41405224 0.58594776
212 1 0.68042096 0.68042096
213 0 0.32540821 0.67459179
214 0 0.08444380 0.91555620
215 1 0.44522680 0.44522680
216 0 0.48742974 0.51257026
217 0 0.80091481 0.19908519
218 0 0.08122139 0.91877861
219 0 0.10235379 0.89764621
220 0 0.14840944 0.85159056
221 1 0.34430105 0.34430105
222 0 0.07511087 0.92488913
223 0 0.64245214 0.35754786
224 0 0.27218732 0.72781268
225 0 0.08444380 0.91555620
226 0 0.12354978 0.87645022
227 0 0.13303414 0.86696586
228 1 0.08444380 0.08444380
229 0 0.14840944 0.85159056
230 1 0.72468315 0.72468315
231 1 0.55088182 0.55088182
232 0 0.38355137 0.61644863
233 0 0.04469448 0.95530552
234 0 0.09481750 0.90518250
235 0 0.31616589 0.68383411
236 0 0.06414031 0.93585969
237 0 0.36369940 0.63630060
238 1 0.37357216 0.37357216
239 1 0.10631759 0.10631759
240 0 0.08444380 0.91555620
241 0 0.09852148 0.90147852
242 1 0.87473555 0.87473555
243 1 0.87931034 0.87931034
244 0 0.40379928 0.59620072
245 0 0.26386534 0.73613466
246 0 0.08444380 0.91555620
247 0 0.09852148 0.90147852
248 0 0.08444380 0.91555620
249 0 0.17115736 0.82884264
250 0 0.14840944 0.85159056
251 0 0.14840944 0.85159056
252 0 0.22477052 0.77522948
253 0 0.14312768 0.85687232
254 0 0.52981267 0.47018733
255 0 0.28931541 0.71068459
256 0 0.20337345 0.79662655
257 1 0.82043312 0.82043312
258 1 0.22477052 0.22477052
259 0 0.23224851 0.76775149
260 0 0.07811146 0.92188854
261 0 0.32540821 0.67459179
262 0 0.10235379 0.89764621
263 0 0.20337345 0.79662655
264 0 0.33478843 0.66521157
265 1 0.84959008 0.84959008
266 0 0.11465216 0.88534784
267 0 0.62273686 0.37726314
268 1 0.80759260 0.80759260
269 0 0.18353106 0.81646894
270 1 0.47683844 0.47683844
271 0 0.04469448 0.95530552
272 1 0.23224851 0.23224851
273 0 0.12354978 0.87645022
274 0 0.10631759 0.89368241
275 0 0.12354978 0.87645022
276 0 0.16522309 0.83477691
277 0 0.37357216 0.62642784
278 1 0.13800342 0.13800342
279 0 0.15385113 0.84614887
280 0 0.17115736 0.82884264
281 1 0.83842391 0.83842391
282 1 0.45572761 0.45572761
283 0 0.27218732 0.72781268
284 1 0.65213658 0.65213658
285 0 0.15945498 0.84054502
286 0 0.23224851 0.76775149
287 1 0.55088182 0.55088182
288 1 0.22477052 0.22477052
289 0 0.14312768 0.85687232
290 0 0.74128419 0.25871581
291 0 0.04469448 0.95530552
292 0 0.32540821 0.67459179
293 1 0.71613927 0.71613927
294 0 0.12821718 0.87178282
295 0 0.38355137 0.61644863
296 0 0.19658713 0.80341287
297 1 0.51923334 0.51923334
298 1 0.77247470 0.77247470
299 0 0.07221652 0.92778348
300 0 0.36369940 0.63630060
301 1 0.11902908 0.11902908
302 0 0.35394014 0.64605986
303 1 0.38355137 0.38355137
304 0 0.43477468 0.56522532
305 1 0.87931034 0.87931034
306 1 0.80091481 0.80091481
307 0 0.12354978 0.87645022
308 0 0.20337345 0.79662655
309 0 0.27218732 0.72781268
310 0 0.26386534 0.73613466
311 0 0.50863673 0.49136327
312 1 0.80759260 0.80759260
313 0 0.18997325 0.81002675
314 0 0.35394014 0.64605986
315 0 0.03956489 0.96043511
316 0 0.14840944 0.85159056
317 1 0.29811500 0.29811500
318 0 0.16522309 0.83477691
319 0 0.21033272 0.78966728
320 0 0.63264996 0.36735004
321 0 0.15385113 0.84614887
322 1 0.57176996 0.57176996
323 1 0.87931034 0.87931034
324 0 0.30706659 0.69293341
325 0 0.20337345 0.79662655
326 1 0.84959008 0.84959008
327 1 0.37357216 0.37357216
328 0 0.09852148 0.90147852
329 1 0.77984422 0.77984422
330 0 0.15945498 0.84054502
331 0 0.30706659 0.69293341
332 0 0.11902908 0.88097092
The likelihood is the product of predicted probabilities (of the actual label) For actual 1s, we use y.hat For actual 0s we use 1-y.hat
<- prod (y * y.hat + (1 - y)* (1 - y.hat))
Maximizing likelihood is equivalent to maximising the log of likelihood, and because logarithms turn products into sums, the log likelihood is easier to work with mathematically. Because likelihood < 1, log likelihood is going to be negative.
<- sum (log (y * y.hat + (1 - y)* (1 - y.hat)))
The closer to zero the log likelihood is, the better; The closer to -infinity, the “worse”. The residual deviance is defined as -2 times the log likelihood. It turns a big negative log likelihood into a large positive number (easier to interpret)
For each observation, we can translate \(P(Y=y | model)\) into a log-likelihood by just taking the log of that decimal. Furthermore, we could multiply each of these by -2 to see which data points contribute more (or less) to the total “residual deviance”
<- log (pY)data.frame (y,y.hat,pY, ll, - 2 * ll)
y y.hat pY ll X.2...ll
1 1 0.58212363 0.58212363 -0.54107244 1.08214488
2 0 0.08778183 0.91221817 -0.09187610 0.18375220
3 0 0.10235379 0.89764621 -0.10797926 0.21595853
4 1 0.06673426 0.06673426 -2.70703677 5.41407354
5 1 0.91759616 0.91759616 -0.08599790 0.17199581
6 1 0.74933615 0.74933615 -0.28856760 0.57713520
7 1 0.28067169 0.28067169 -1.27056964 2.54113928
8 0 0.17115736 0.82884264 -0.18772496 0.37544992
9 0 0.35394014 0.64605986 -0.43686311 0.87372622
10 1 0.28931541 0.28931541 -1.24023781 2.48047561
11 0 0.13800342 0.86199658 -0.14850398 0.29700795
12 0 0.21033272 0.78966728 -0.23614358 0.47228716
13 0 0.09852148 0.90147852 -0.10371906 0.20743812
14 0 0.31616589 0.68383411 -0.38003992 0.76007985
15 0 0.17115736 0.82884264 -0.18772496 0.37544992
16 1 0.16522309 0.16522309 -1.80045868 3.60091737
17 1 0.10631759 0.10631759 -2.24132451 4.48264903
18 1 0.22477052 0.22477052 -1.49267529 2.98535058
19 0 0.84408847 0.15591153 -1.85846657 3.71693313
20 0 0.18997325 0.81002675 -0.21068801 0.42137601
21 1 0.78704085 0.78704085 -0.23947513 0.47895026
22 1 0.84408847 0.84408847 -0.16949796 0.33899593
23 0 0.05045417 0.94954583 -0.05177149 0.10354297
24 0 0.17115736 0.82884264 -0.18772496 0.37544992
25 0 0.15945498 0.84054502 -0.17370476 0.34740952
26 0 0.09852148 0.90147852 -0.10371906 0.20743812
27 0 0.60260691 0.39739309 -0.92282935 1.84565869
28 0 0.05467737 0.94532263 -0.05622900 0.11245800
29 0 0.18353106 0.81646894 -0.20276640 0.40553280
30 0 0.14840944 0.85159056 -0.16064944 0.32129887
31 1 0.21033272 0.21033272 -1.55906464 3.11812928
32 0 0.12821718 0.87178282 -0.13721495 0.27442989
33 0 0.56135307 0.43864693 -0.82406045 1.64812091
34 1 0.15385113 0.15385113 -1.87176985 3.74353970
35 0 0.38355137 0.61644863 -0.48378028 0.96756057
36 0 0.12821718 0.87178282 -0.13721495 0.27442989
37 0 0.23224851 0.76775149 -0.26428918 0.52857837
38 0 0.23989830 0.76010170 -0.27430304 0.54860609
39 0 0.08122139 0.91877861 -0.08471009 0.16942018
40 0 0.15945498 0.84054502 -0.17370476 0.34740952
41 0 0.18997325 0.81002675 -0.21068801 0.42137601
42 0 0.15385113 0.84614887 -0.16705996 0.33411992
43 1 0.45572761 0.45572761 -0.78586000 1.57172000
44 0 0.32540821 0.67459179 -0.39364753 0.78729506
45 0 0.07511087 0.92488913 -0.07808141 0.15616282
46 0 0.11465216 0.88534784 -0.12177467 0.24354935
47 0 0.11902908 0.88097092 -0.12673066 0.25346133
48 1 0.31616589 0.31616589 -1.15148823 2.30297645
49 0 0.07511087 0.92488913 -0.07808141 0.15616282
50 0 0.35394014 0.64605986 -0.43686311 0.87372622
51 0 0.54036528 0.45963472 -0.77732320 1.55464640
52 0 0.10235379 0.89764621 -0.10797926 0.21595853
53 1 0.69858083 0.69858083 -0.35870439 0.71740879
54 0 0.13800342 0.86199658 -0.14850398 0.29700795
55 1 0.71613927 0.71613927 -0.33388062 0.66776123
56 0 0.19658713 0.80341287 -0.21888653 0.43777307
57 1 0.09852148 0.09852148 -2.31748072 4.63496144
58 1 0.27218732 0.27218732 -1.30126478 2.60252956
59 1 0.80091481 0.80091481 -0.22200070 0.44400140
60 1 0.77984422 0.77984422 -0.24866110 0.49732220
61 0 0.08444380 0.91555620 -0.08822352 0.17644705
62 0 0.15385113 0.84614887 -0.16705996 0.33411992
63 0 0.11902908 0.88097092 -0.12673066 0.25346133
64 0 0.18997325 0.81002675 -0.21068801 0.42137601
65 0 0.20337345 0.79662655 -0.22736928 0.45473857
66 0 0.18997325 0.81002675 -0.21068801 0.42137601
67 0 0.10631759 0.89368241 -0.11240482 0.22480963
68 1 0.66169683 0.66169683 -0.41294778 0.82589557
69 0 0.63264996 0.36735004 -1.00144010 2.00288021
70 1 0.62273686 0.62273686 -0.47363122 0.94726245
71 0 0.09852148 0.90147852 -0.10371906 0.20743812
72 1 0.72468315 0.72468315 -0.32202075 0.64404150
73 0 0.61272001 0.38727999 -0.94860735 1.89721471
74 0 0.16522309 0.83477691 -0.18059076 0.36118152
75 0 0.24771876 0.75228124 -0.28464504 0.56929008
76 1 0.17725952 0.17725952 -1.73014042 3.46028084
77 0 0.05923201 0.94076799 -0.06105873 0.12211746
78 1 0.83842391 0.83842391 -0.17623145 0.35246290
79 1 0.38355137 0.38355137 -0.95828172 1.91656344
80 1 0.37357216 0.37357216 -0.98464410 1.96928820
81 1 0.21033272 0.21033272 -1.55906464 3.11812928
82 1 0.21033272 0.21033272 -1.55906464 3.11812928
83 0 0.23989830 0.76010170 -0.27430304 0.54860609
84 0 0.20337345 0.79662655 -0.22736928 0.45473857
85 0 0.22477052 0.77522948 -0.25459620 0.50919239
86 1 0.91433121 0.91433121 -0.08956240 0.17912480
87 0 0.07511087 0.92488913 -0.07808141 0.15616282
88 0 0.57176996 0.42823004 -0.84809476 1.69618951
89 1 0.21033272 0.21033272 -1.55906464 3.11812928
90 0 0.34430105 0.65569895 -0.42205351 0.84410702
91 1 0.08778183 0.08778183 -2.43290071 4.86580141
92 0 0.51923334 0.48076666 -0.73237323 1.46474646
93 0 0.15385113 0.84614887 -0.16705996 0.33411992
94 0 0.09481750 0.90518250 -0.09961869 0.19923739
95 1 0.71613927 0.71613927 -0.33388062 0.66776123
96 0 0.91759616 0.08240384 -2.49612319 4.99224637
97 0 0.27218732 0.72781268 -0.31771157 0.63542314
98 1 0.43477468 0.43477468 -0.83292736 1.66585472
99 0 0.05691189 0.94308811 -0.05859557 0.11719113
100 1 0.84959008 0.84959008 -0.16300131 0.32600262
101 1 0.83842391 0.83842391 -0.17623145 0.35246290
102 0 0.11041602 0.88958398 -0.11700136 0.23400272
103 0 0.11041602 0.88958398 -0.11700136 0.23400272
104 1 0.28931541 0.28931541 -1.24023781 2.48047561
105 0 0.56135307 0.43864693 -0.82406045 1.64812091
106 1 0.86513981 0.86513981 -0.14486415 0.28972831
107 0 0.74128419 0.25871581 -1.35202509 2.70405018
108 0 0.33478843 0.66521157 -0.40765015 0.81530029
109 0 0.22477052 0.77522948 -0.25459620 0.50919239
110 0 0.10631759 0.89368241 -0.11240482 0.22480963
111 1 0.23989830 0.23989830 -1.42754018 2.85508036
112 0 0.22477052 0.77522948 -0.25459620 0.50919239
113 1 0.65213658 0.65213658 -0.42750126 0.85500253
114 0 0.12821718 0.87178282 -0.13721495 0.27442989
115 0 0.13303414 0.86696586 -0.14275568 0.28551136
116 0 0.37357216 0.62642784 -0.46772169 0.93544337
117 1 0.15945498 0.15945498 -1.83599367 3.67198734
118 0 0.20337345 0.79662655 -0.22736928 0.45473857
119 0 0.15385113 0.84614887 -0.16705996 0.33411992
120 1 0.18997325 0.18997325 -1.66087201 3.32174402
121 0 0.17725952 0.82274048 -0.19511446 0.39022892
122 0 0.20337345 0.79662655 -0.22736928 0.45473857
123 0 0.42438023 0.57561977 -0.55230795 1.10461591
124 0 0.45572761 0.54427239 -0.60830543 1.21661087
125 1 0.28931541 0.28931541 -1.24023781 2.48047561
126 0 0.13303414 0.86696586 -0.14275568 0.28551136
127 0 0.06673426 0.93326574 -0.06906530 0.13813060
128 1 0.19658713 0.19658713 -1.62664955 3.25329911
129 0 0.61272001 0.38727999 -0.94860735 1.89721471
130 1 0.56135307 0.56135307 -0.57740521 1.15481042
131 0 0.35394014 0.64605986 -0.43686311 0.87372622
132 1 0.54036528 0.54036528 -0.61550992 1.23101983
133 0 0.29811500 0.70188500 -0.35398570 0.70797141
134 1 0.70743729 0.70743729 -0.34610629 0.69221257
135 1 0.37357216 0.37357216 -0.98464410 1.96928820
136 1 0.33478843 0.33478843 -1.09425649 2.18851297
137 1 0.65213658 0.65213658 -0.42750126 0.85500253
138 0 0.23989830 0.76010170 -0.27430304 0.54860609
139 1 0.21033272 0.21033272 -1.55906464 3.11812928
140 0 0.25570840 0.74429160 -0.29532238 0.59064477
141 1 0.62273686 0.62273686 -0.47363122 0.94726245
142 0 0.23224851 0.76775149 -0.26428918 0.52857837
143 0 0.67112685 0.32887315 -1.11208316 2.22416631
144 1 0.31616589 0.31616589 -1.15148823 2.30297645
145 1 0.16522309 0.16522309 -1.80045868 3.60091737
146 0 0.18353106 0.81646894 -0.20276640 0.40553280
147 0 0.28067169 0.71932831 -0.32943741 0.65887482
148 0 0.09481750 0.90518250 -0.09961869 0.19923739
149 0 0.12821718 0.87178282 -0.13721495 0.27442989
150 0 0.39362965 0.60637035 -0.50026434 1.00052868
151 0 0.12821718 0.87178282 -0.13721495 0.27442989
152 0 0.35394014 0.64605986 -0.43686311 0.87372622
153 0 0.48742974 0.51257026 -0.66831749 1.33663498
154 0 0.14840944 0.85159056 -0.16064944 0.32129887
155 1 0.34430105 0.34430105 -1.06623887 2.13247774
156 1 0.91433121 0.91433121 -0.08956240 0.17912480
157 1 0.88802830 0.88802830 -0.11875167 0.23750334
158 0 0.17115736 0.82884264 -0.18772496 0.37544992
159 0 0.57176996 0.42823004 -0.84809476 1.69618951
160 0 0.14840944 0.85159056 -0.16064944 0.32129887
161 0 0.07511087 0.92488913 -0.07808141 0.15616282
162 1 0.80091481 0.80091481 -0.22200070 0.44400140
163 0 0.08444380 0.91555620 -0.08822352 0.17644705
164 0 0.14312768 0.85687232 -0.15446635 0.30893271
165 0 0.09481750 0.90518250 -0.09961869 0.19923739
166 0 0.11902908 0.88097092 -0.12673066 0.25346133
167 0 0.18353106 0.81646894 -0.20276640 0.40553280
168 0 0.10631759 0.89368241 -0.11240482 0.22480963
169 0 0.34430105 0.65569895 -0.42205351 0.84410702
170 0 0.28931541 0.71068459 -0.34152656 0.68305312
171 1 0.26386534 0.26386534 -1.33231640 2.66463280
172 1 0.18353106 0.18353106 -1.69537138 3.39074277
173 0 0.15385113 0.84614887 -0.16705996 0.33411992
174 0 0.40379928 0.59620072 -0.51717789 1.03435577
175 0 0.36369940 0.63630060 -0.45208418 0.90416837
176 0 0.13303414 0.86696586 -0.14275568 0.28551136
177 1 0.90382283 0.90382283 -0.10112192 0.20224384
178 1 0.45572761 0.45572761 -0.78586000 1.57172000
179 0 0.05252569 0.94747431 -0.05395546 0.10791092
180 1 0.79406435 0.79406435 -0.23059078 0.46118155
181 0 0.16522309 0.83477691 -0.18059076 0.36118152
182 0 0.23224851 0.76775149 -0.26428918 0.52857837
183 0 0.52981267 0.47018733 -0.75462409 1.50924818
184 1 0.80091481 0.80091481 -0.22200070 0.44400140
185 1 0.54036528 0.54036528 -0.61550992 1.23101983
186 1 0.38355137 0.38355137 -0.95828172 1.91656344
187 0 0.28931541 0.71068459 -0.34152656 0.68305312
188 0 0.12354978 0.87645022 -0.13187537 0.26375074
189 0 0.16522309 0.83477691 -0.18059076 0.36118152
190 1 0.61272001 0.61272001 -0.48984720 0.97969441
191 1 0.84959008 0.84959008 -0.16300131 0.32600262
192 1 0.12821718 0.12821718 -2.05402974 4.10805948
193 0 0.10235379 0.89764621 -0.10797926 0.21595853
194 0 0.07221652 0.92778348 -0.07495690 0.14991379
195 0 0.10631759 0.89368241 -0.11240482 0.22480963
196 1 0.88802830 0.88802830 -0.11875167 0.23750334
197 0 0.27218732 0.72781268 -0.31771157 0.63542314
198 1 0.84408847 0.84408847 -0.16949796 0.33899593
199 1 0.17725952 0.17725952 -1.73014042 3.46028084
200 0 0.29811500 0.70188500 -0.35398570 0.70797141
201 0 0.07811146 0.92188854 -0.08133095 0.16266190
202 0 0.11041602 0.88958398 -0.11700136 0.23400272
203 0 0.43477468 0.56522532 -0.57053083 1.14106166
204 0 0.29811500 0.70188500 -0.35398570 0.70797141
205 0 0.05691189 0.94308811 -0.05859557 0.11719113
206 0 0.09852148 0.90147852 -0.10371906 0.20743812
207 0 0.33478843 0.66521157 -0.40765015 0.81530029
208 0 0.13800342 0.86199658 -0.14850398 0.29700795
209 0 0.54036528 0.45963472 -0.77732320 1.55464640
210 0 0.46626792 0.53373208 -0.62786129 1.25572258
211 0 0.41405224 0.58594776 -0.53452464 1.06904928
212 1 0.68042096 0.68042096 -0.38504362 0.77008724
213 0 0.32540821 0.67459179 -0.39364753 0.78729506
214 0 0.08444380 0.91555620 -0.08822352 0.17644705
215 1 0.44522680 0.44522680 -0.80917145 1.61834291
216 0 0.48742974 0.51257026 -0.66831749 1.33663498
217 0 0.80091481 0.19908519 -1.61402243 3.22804487
218 0 0.08122139 0.91877861 -0.08471009 0.16942018
219 0 0.10235379 0.89764621 -0.10797926 0.21595853
220 0 0.14840944 0.85159056 -0.16064944 0.32129887
221 1 0.34430105 0.34430105 -1.06623887 2.13247774
222 0 0.07511087 0.92488913 -0.07808141 0.15616282
223 0 0.64245214 0.35754786 -1.02848605 2.05697210
224 0 0.27218732 0.72781268 -0.31771157 0.63542314
225 0 0.08444380 0.91555620 -0.08822352 0.17644705
226 0 0.12354978 0.87645022 -0.13187537 0.26375074
227 0 0.13303414 0.86696586 -0.14275568 0.28551136
228 1 0.08444380 0.08444380 -2.47166911 4.94333822
229 0 0.14840944 0.85159056 -0.16064944 0.32129887
230 1 0.72468315 0.72468315 -0.32202075 0.64404150
231 1 0.55088182 0.55088182 -0.59623497 1.19246993
232 0 0.38355137 0.61644863 -0.48378028 0.96756057
233 0 0.04469448 0.95530552 -0.04572407 0.09144814
234 0 0.09481750 0.90518250 -0.09961869 0.19923739
235 0 0.31616589 0.68383411 -0.38003992 0.76007985
236 0 0.06414031 0.93585969 -0.06628972 0.13257945
237 0 0.36369940 0.63630060 -0.45208418 0.90416837
238 1 0.37357216 0.37357216 -0.98464410 1.96928820
239 1 0.10631759 0.10631759 -2.24132451 4.48264903
240 0 0.08444380 0.91555620 -0.08822352 0.17644705
241 0 0.09852148 0.90147852 -0.10371906 0.20743812
242 1 0.87473555 0.87473555 -0.13383367 0.26766734
243 1 0.87931034 0.87931034 -0.12861738 0.25723476
244 0 0.40379928 0.59620072 -0.51717789 1.03435577
245 0 0.26386534 0.73613466 -0.30634221 0.61268442
246 0 0.08444380 0.91555620 -0.08822352 0.17644705
247 0 0.09852148 0.90147852 -0.10371906 0.20743812
248 0 0.08444380 0.91555620 -0.08822352 0.17644705
249 0 0.17115736 0.82884264 -0.18772496 0.37544992
250 0 0.14840944 0.85159056 -0.16064944 0.32129887
251 0 0.14840944 0.85159056 -0.16064944 0.32129887
252 0 0.22477052 0.77522948 -0.25459620 0.50919239
253 0 0.14312768 0.85687232 -0.15446635 0.30893271
254 0 0.52981267 0.47018733 -0.75462409 1.50924818
255 0 0.28931541 0.71068459 -0.34152656 0.68305312
256 0 0.20337345 0.79662655 -0.22736928 0.45473857
257 1 0.82043312 0.82043312 -0.19792288 0.39584576
258 1 0.22477052 0.22477052 -1.49267529 2.98535058
259 0 0.23224851 0.76775149 -0.26428918 0.52857837
260 0 0.07811146 0.92188854 -0.08133095 0.16266190
261 0 0.32540821 0.67459179 -0.39364753 0.78729506
262 0 0.10235379 0.89764621 -0.10797926 0.21595853
263 0 0.20337345 0.79662655 -0.22736928 0.45473857
264 0 0.33478843 0.66521157 -0.40765015 0.81530029
265 1 0.84959008 0.84959008 -0.16300131 0.32600262
266 0 0.11465216 0.88534784 -0.12177467 0.24354935
267 0 0.62273686 0.37726314 -0.97481235 1.94962471
268 1 0.80759260 0.80759260 -0.21369756 0.42739511
269 0 0.18353106 0.81646894 -0.20276640 0.40553280
270 1 0.47683844 0.47683844 -0.74057755 1.48115510
271 0 0.04469448 0.95530552 -0.04572407 0.09144814
272 1 0.23224851 0.23224851 -1.45994730 2.91989460
273 0 0.12354978 0.87645022 -0.13187537 0.26375074
274 0 0.10631759 0.89368241 -0.11240482 0.22480963
275 0 0.12354978 0.87645022 -0.13187537 0.26375074
276 0 0.16522309 0.83477691 -0.18059076 0.36118152
277 0 0.37357216 0.62642784 -0.46772169 0.93544337
278 1 0.13800342 0.13800342 -1.98047681 3.96095362
279 0 0.15385113 0.84614887 -0.16705996 0.33411992
280 0 0.17115736 0.82884264 -0.18772496 0.37544992
281 1 0.83842391 0.83842391 -0.17623145 0.35246290
282 1 0.45572761 0.45572761 -0.78586000 1.57172000
283 0 0.27218732 0.72781268 -0.31771157 0.63542314
284 1 0.65213658 0.65213658 -0.42750126 0.85500253
285 0 0.15945498 0.84054502 -0.17370476 0.34740952
286 0 0.23224851 0.76775149 -0.26428918 0.52857837
287 1 0.55088182 0.55088182 -0.59623497 1.19246993
288 1 0.22477052 0.22477052 -1.49267529 2.98535058
289 0 0.14312768 0.85687232 -0.15446635 0.30893271
290 0 0.74128419 0.25871581 -1.35202509 2.70405018
291 0 0.04469448 0.95530552 -0.04572407 0.09144814
292 0 0.32540821 0.67459179 -0.39364753 0.78729506
293 1 0.71613927 0.71613927 -0.33388062 0.66776123
294 0 0.12821718 0.87178282 -0.13721495 0.27442989
295 0 0.38355137 0.61644863 -0.48378028 0.96756057
296 0 0.19658713 0.80341287 -0.21888653 0.43777307
297 1 0.51923334 0.51923334 -0.65540191 1.31080382
298 1 0.77247470 0.77247470 -0.25815602 0.51631205
299 0 0.07221652 0.92778348 -0.07495690 0.14991379
300 0 0.36369940 0.63630060 -0.45208418 0.90416837
301 1 0.11902908 0.11902908 -2.12838742 4.25677484
302 0 0.35394014 0.64605986 -0.43686311 0.87372622
303 1 0.38355137 0.38355137 -0.95828172 1.91656344
304 0 0.43477468 0.56522532 -0.57053083 1.14106166
305 1 0.87931034 0.87931034 -0.12861738 0.25723476
306 1 0.80091481 0.80091481 -0.22200070 0.44400140
307 0 0.12354978 0.87645022 -0.13187537 0.26375074
308 0 0.20337345 0.79662655 -0.22736928 0.45473857
309 0 0.27218732 0.72781268 -0.31771157 0.63542314
310 0 0.26386534 0.73613466 -0.30634221 0.61268442
311 0 0.50863673 0.49136327 -0.71057156 1.42114312
312 1 0.80759260 0.80759260 -0.21369756 0.42739511
313 0 0.18997325 0.81002675 -0.21068801 0.42137601
314 0 0.35394014 0.64605986 -0.43686311 0.87372622
315 0 0.03956489 0.96043511 -0.04036886 0.08073772
316 0 0.14840944 0.85159056 -0.16064944 0.32129887
317 1 0.29811500 0.29811500 -1.21027597 2.42055194
318 0 0.16522309 0.83477691 -0.18059076 0.36118152
319 0 0.21033272 0.78966728 -0.23614358 0.47228716
320 0 0.63264996 0.36735004 -1.00144010 2.00288021
321 0 0.15385113 0.84614887 -0.16705996 0.33411992
322 1 0.57176996 0.57176996 -0.55901853 1.11803706
323 1 0.87931034 0.87931034 -0.12861738 0.25723476
324 0 0.30706659 0.69293341 -0.36682137 0.73364274
325 0 0.20337345 0.79662655 -0.22736928 0.45473857
326 1 0.84959008 0.84959008 -0.16300131 0.32600262
327 1 0.37357216 0.37357216 -0.98464410 1.96928820
328 0 0.09852148 0.90147852 -0.10371906 0.20743812
329 1 0.77984422 0.77984422 -0.24866110 0.49732220
330 0 0.15945498 0.84054502 -0.17370476 0.34740952
331 0 0.30706659 0.69293341 -0.36682137 0.73364274
332 0 0.11902908 0.88097092 -0.12673066 0.25346133
Prediction
Suppose someone has glucose level 150; what is the estimated probability of diabetes? \[\hat{logodds}(Y+i=1) = \hat{\beta_0}+\hat{\beta_1}X_i\]
#Predict log odds predict (Pima.glm, newdata= data.frame (glu= 150 ))
#Predicts probability predict (Pima.glm, newdata= data.frame (glu= 150 ), type= "response" )
#calculation .150 <- sum (Pima.glm$ coefficients * c (1 , 150 ))#Estimated probability = e^log-odds / (1+e^log-odds) #or # = 1/ (1 + e^ -logodds) exp (log.odds.150 ) / (1 + exp (log.odds.150 ))#or 1 / (1 + exp (- log.odds.150 ))
The null model
The null model would be what we compare our “better” model to - it is a logistic model that does not have any predictors. Just like the null model for linear regression predicts \(\hat{y}=\bar{y}\) the null model predicts \(\hat{y}=\) the proportion of 1s in the data set.
<- glm (typeTF ~ 1 , data= Pima.te, family= "binomial" )summary (pima.glm.null)
Call:
glm(formula = typeTF ~ 1, family = "binomial", data = Pima.te)
Deviance Residuals:
Min 1Q Median 3Q Max
-0.8921 -0.8921 -0.8921 1.4925 1.4925
Coefficients:
Estimate Std. Error z value Pr(>|z|)
(Intercept) -0.7158 0.1169 -6.125 9.07e-10 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
(Dispersion parameter for binomial family taken to be 1)
Null deviance: 420.3 on 331 degrees of freedom
Residual deviance: 420.3 on 331 degrees of freedom
AIC: 422.3
Number of Fisher Scoring iterations: 4
<- pima.glm.null$ coefficients[1 ]exp (log.odds) / (1 + exp (log.odds))
Let’s just double check that. We can calculate the coefficient manually
#What is the proportion of 1s in the data? <- mean (Pima.te$ typeTF))#convert this to a log-odds log (prop1/ (1 - prop1))
This is the same as the constant in the null model.
#null deviance calculation #log likelihood = sum [y * y.hat + (1-y)*(1-yhat)] #deviance = -2*log likelihood <- nrow (Pima.te)<- sum (Pima.te$ typeTF)- 2 * (m* log (m/ n) + (n- m)* log ((n- m)/ n))
Results of logistic regression: confusion matrix
See https://en.wikipedia.org/wiki/Confusion_matrix
<- glm (typeTF ~ glu, data= Pima.te, family= "binomial" )<- table (Pima.te$ typeTF, factor (predict (pima_logit, type= "response" ) >= .5 , levels= c (FALSE , TRUE )))colnames (confusion.table) <- c ("Predict 0" ,"Predict 1" )rownames (confusion.table) <- c ("No diabetes" ,"Yes diabetes" )#confusion.table <- confusion.table[2 ,2 ]<- confusion.table[1 ,1 ]<- confusion.table[1 ,2 ]<- confusion.table[2 ,1 ]t (addmargins (confusion.table))
No diabetes Yes diabetes Sum
Predict 0 200 53 253
Predict 1 23 56 79
Sum 223 109 332
Statistics for a predictor:
Accuracy of a classifier - the total # correct predictions / total predictions
<- (TP+ TN)/ (TP+ TN+ FP+ FN))
Sensitivity (aka recall, hit rate, true positive rate, aka POWER) - what proportion of actual positive cases are classified as positive
<- TP / (TP+ FN))
Specificity (selectivity, true negative rate, \(1-\alpha\) ) - what proportion of actual negative cases are classified as negative
<- TN / (TN+ FP))
Precision (positive predictive value) - what proportion of positive predictions are correct
Negative Predictive Value (NPV) - what proportion of negative predictions are correct
What was the Type 1 error rate? The False positive rate
Let’s just change the theshold from \(p=.5\) to some other number to see what happens to this matrix.
<- table (Pima.te$ typeTF, factor (predict (pima_logit, type= "response" ) >= .6 , levels= c (FALSE , TRUE )))colnames (confusion.table) <- c ("Predict 0" ,"Predict 1" )rownames (confusion.table) <- c ("No diabetes" ,"Yes diabetes" )#confusion.table t (addmargins (confusion.table))
No diabetes Yes diabetes Sum
Predict 0 210 61 271
Predict 1 13 48 61
Sum 223 109 332
#column proportions prop.table (t (confusion.table), 2 )
No diabetes Yes diabetes
Predict 0 0.94170404 0.55963303
Predict 1 0.05829596 0.44036697
If I look at column proportions, the first column (bottom row) is the alpha, type 1 error rate i.e. FPR. Column 2, in row 2 gives power.
The accuracy is (210 + 48 ) / 332
By adjusting the positive prediction threshold we can improve/worsen statistics such as accuracy, specificity, sensitivity, etc.
<- glm (typeTF ~ glu, data= Pima.te, family= "binomial" )<- 0 <- seq (0.01 , 0.99 , .01 )for (i in 1 : length (thresh)){<- table (Pima.te$ typeTF, factor (predict (pima_logit, type= "response" )>= thresh[i], levels= c (FALSE , TRUE )))colnames (confusion.table) <- c ("Predict 0" ,"Predict 1" )rownames (confusion.table) <- c ("No diabetes" ,"Yes diabetes" )#confusion.table <- confusion.table[2 ,2 ]<- confusion.table[1 ,1 ]<- confusion.table[1 ,2 ]<- confusion.table[2 ,1 ]= (TP+ TN)/ (TP+ FP+ TN+ FN)plot (thresh, acc, type= "l" , main= "Accuracy by threshold" , ylim= c (0 ,1 ))
In particular, we can look at the relationship between True positive rate and False Positive Rate (what proportion of actual negatives are incorrectly predicted to be positive) (FP / (FP+TN))
<- 0 ; FPR <- 0 ; <- seq (0 , 1 , .001 )for (i in 1 : length (thresh)){<- table (Pima.te$ typeTF, factor (predict (pima_logit, type= "response" )>= thresh[i], levels= c (FALSE , TRUE )))colnames (confusion.table) <- c ("Predict 0" ,"Predict 1" )rownames (confusion.table) <- c ("No diabetes" ,"Yes diabetes" )#confusion.table <- confusion.table[2 ,2 ]<- confusion.table[1 ,1 ]<- confusion.table[1 ,2 ]<- confusion.table[2 ,1 ]<- TP/ (TP+ FN) #power <- FP/ (FP+ TN) #alpha / significance / type 1 error plot (thresh, TPR, type= "n" , main= "TPR and FPR by threshold" , xlim= c (0 ,1 ), ylim= c (0 ,1 ))lines (thresh, TPR, lty= 1 )lines (thresh, FPR, lty= 2 )abline (h= 0.05 , col= "red" , lty= 3 )legend (x= .8 , y= 1 , legend = c ("TPR" , "FPR" ), lty= c (1 ,2 ))
To put things into the language of hypothesis testing we could ask this: What threshold would we use to achieve a 5% significance level (false positive rate) and what would the power of the classifier be?
<- min (which (FPR<= 0.05 ))as.numeric (FPR[optimalThresh])as.numeric (TPR[optimalThresh])
TPR FPR
1 1.000000000 1.000000000
2 1.000000000 1.000000000
3 1.000000000 1.000000000
4 1.000000000 1.000000000
5 1.000000000 1.000000000
6 1.000000000 1.000000000
7 1.000000000 1.000000000
8 1.000000000 1.000000000
9 1.000000000 1.000000000
10 1.000000000 1.000000000
11 1.000000000 1.000000000
12 1.000000000 1.000000000
13 1.000000000 1.000000000
14 1.000000000 1.000000000
15 1.000000000 1.000000000
16 1.000000000 1.000000000
17 1.000000000 1.000000000
18 1.000000000 1.000000000
19 1.000000000 1.000000000
20 1.000000000 1.000000000
21 1.000000000 1.000000000
22 1.000000000 1.000000000
23 1.000000000 1.000000000
24 1.000000000 1.000000000
25 1.000000000 1.000000000
26 1.000000000 1.000000000
27 1.000000000 1.000000000
28 1.000000000 1.000000000
29 1.000000000 1.000000000
30 1.000000000 1.000000000
31 1.000000000 1.000000000
32 1.000000000 1.000000000
33 1.000000000 1.000000000
34 1.000000000 1.000000000
35 1.000000000 1.000000000
36 1.000000000 1.000000000
37 1.000000000 1.000000000
38 1.000000000 1.000000000
39 1.000000000 1.000000000
40 1.000000000 1.000000000
41 1.000000000 0.995515695
42 1.000000000 0.995515695
43 1.000000000 0.995515695
44 1.000000000 0.995515695
45 1.000000000 0.995515695
46 1.000000000 0.982062780
47 1.000000000 0.982062780
48 1.000000000 0.982062780
49 1.000000000 0.982062780
50 1.000000000 0.982062780
51 1.000000000 0.982062780
52 1.000000000 0.977578475
53 1.000000000 0.977578475
54 1.000000000 0.973094170
55 1.000000000 0.973094170
56 1.000000000 0.968609865
57 1.000000000 0.968609865
58 1.000000000 0.959641256
59 1.000000000 0.959641256
60 1.000000000 0.959641256
61 1.000000000 0.955156951
62 1.000000000 0.955156951
63 1.000000000 0.955156951
64 1.000000000 0.955156951
65 1.000000000 0.955156951
66 1.000000000 0.950672646
67 1.000000000 0.950672646
68 0.990825688 0.946188341
69 0.990825688 0.946188341
70 0.990825688 0.946188341
71 0.990825688 0.946188341
72 0.990825688 0.946188341
73 0.990825688 0.946188341
74 0.990825688 0.937219731
75 0.990825688 0.937219731
76 0.990825688 0.937219731
77 0.990825688 0.914798206
78 0.990825688 0.914798206
79 0.990825688 0.914798206
80 0.990825688 0.905829596
81 0.990825688 0.905829596
82 0.990825688 0.905829596
83 0.990825688 0.896860987
84 0.990825688 0.896860987
85 0.990825688 0.896860987
86 0.981651376 0.865470852
87 0.981651376 0.865470852
88 0.981651376 0.865470852
89 0.972477064 0.860986547
90 0.972477064 0.860986547
91 0.972477064 0.860986547
92 0.972477064 0.860986547
93 0.972477064 0.860986547
94 0.972477064 0.860986547
95 0.972477064 0.860986547
96 0.972477064 0.843049327
97 0.972477064 0.843049327
98 0.972477064 0.843049327
99 0.972477064 0.843049327
100 0.963302752 0.811659193
101 0.963302752 0.811659193
102 0.963302752 0.811659193
103 0.963302752 0.811659193
104 0.963302752 0.789237668
105 0.963302752 0.789237668
106 0.963302752 0.789237668
107 0.963302752 0.789237668
108 0.944954128 0.766816143
109 0.944954128 0.766816143
110 0.944954128 0.766816143
111 0.944954128 0.766816143
112 0.944954128 0.753363229
113 0.944954128 0.753363229
114 0.944954128 0.753363229
115 0.944954128 0.753363229
116 0.944954128 0.744394619
117 0.944954128 0.744394619
118 0.944954128 0.744394619
119 0.944954128 0.744394619
120 0.944954128 0.744394619
121 0.935779817 0.726457399
122 0.935779817 0.726457399
123 0.935779817 0.726457399
124 0.935779817 0.726457399
125 0.935779817 0.704035874
126 0.935779817 0.704035874
127 0.935779817 0.704035874
128 0.935779817 0.704035874
129 0.935779817 0.704035874
130 0.926605505 0.677130045
131 0.926605505 0.677130045
132 0.926605505 0.677130045
133 0.926605505 0.677130045
134 0.926605505 0.677130045
135 0.926605505 0.659192825
136 0.926605505 0.659192825
137 0.926605505 0.659192825
138 0.926605505 0.659192825
139 0.926605505 0.659192825
140 0.917431193 0.645739910
141 0.917431193 0.645739910
142 0.917431193 0.645739910
143 0.917431193 0.645739910
144 0.917431193 0.645739910
145 0.917431193 0.632286996
146 0.917431193 0.632286996
147 0.917431193 0.632286996
148 0.917431193 0.632286996
149 0.917431193 0.632286996
150 0.917431193 0.596412556
151 0.917431193 0.596412556
152 0.917431193 0.596412556
153 0.917431193 0.596412556
154 0.917431193 0.596412556
155 0.908256881 0.565022422
156 0.908256881 0.565022422
157 0.908256881 0.565022422
158 0.908256881 0.565022422
159 0.908256881 0.565022422
160 0.908256881 0.565022422
161 0.899082569 0.547085202
162 0.899082569 0.547085202
163 0.899082569 0.547085202
164 0.899082569 0.547085202
165 0.899082569 0.547085202
166 0.899082569 0.547085202
167 0.880733945 0.524663677
168 0.880733945 0.524663677
169 0.880733945 0.524663677
170 0.880733945 0.524663677
171 0.880733945 0.524663677
172 0.880733945 0.524663677
173 0.880733945 0.497757848
174 0.880733945 0.497757848
175 0.880733945 0.497757848
176 0.880733945 0.497757848
177 0.880733945 0.497757848
178 0.880733945 0.497757848
179 0.862385321 0.493273543
180 0.862385321 0.493273543
181 0.862385321 0.493273543
182 0.862385321 0.493273543
183 0.862385321 0.493273543
184 0.862385321 0.493273543
185 0.853211009 0.475336323
186 0.853211009 0.475336323
187 0.853211009 0.475336323
188 0.853211009 0.475336323
189 0.853211009 0.475336323
190 0.853211009 0.475336323
191 0.844036697 0.452914798
192 0.844036697 0.452914798
193 0.844036697 0.452914798
194 0.844036697 0.452914798
195 0.844036697 0.452914798
196 0.844036697 0.452914798
197 0.844036697 0.452914798
198 0.834862385 0.443946188
199 0.834862385 0.443946188
200 0.834862385 0.443946188
201 0.834862385 0.443946188
202 0.834862385 0.443946188
203 0.834862385 0.443946188
204 0.834862385 0.443946188
205 0.834862385 0.408071749
206 0.834862385 0.408071749
207 0.834862385 0.408071749
208 0.834862385 0.408071749
209 0.834862385 0.408071749
210 0.834862385 0.408071749
211 0.834862385 0.408071749
212 0.788990826 0.399103139
213 0.788990826 0.399103139
214 0.788990826 0.399103139
215 0.788990826 0.399103139
216 0.788990826 0.399103139
217 0.788990826 0.399103139
218 0.788990826 0.399103139
219 0.788990826 0.399103139
220 0.788990826 0.399103139
221 0.788990826 0.399103139
222 0.788990826 0.399103139
223 0.788990826 0.399103139
224 0.788990826 0.399103139
225 0.788990826 0.399103139
226 0.761467890 0.381165919
227 0.761467890 0.381165919
228 0.761467890 0.381165919
229 0.761467890 0.381165919
230 0.761467890 0.381165919
231 0.761467890 0.381165919
232 0.761467890 0.381165919
233 0.761467890 0.381165919
234 0.752293578 0.358744395
235 0.752293578 0.358744395
236 0.752293578 0.358744395
237 0.752293578 0.358744395
238 0.752293578 0.358744395
239 0.752293578 0.358744395
240 0.752293578 0.358744395
241 0.743119266 0.345291480
242 0.743119266 0.345291480
243 0.743119266 0.345291480
244 0.743119266 0.345291480
245 0.743119266 0.345291480
246 0.743119266 0.345291480
247 0.743119266 0.345291480
248 0.743119266 0.345291480
249 0.743119266 0.340807175
250 0.743119266 0.340807175
251 0.743119266 0.340807175
252 0.743119266 0.340807175
253 0.743119266 0.340807175
254 0.743119266 0.340807175
255 0.743119266 0.340807175
256 0.743119266 0.340807175
257 0.743119266 0.336322870
258 0.743119266 0.336322870
259 0.743119266 0.336322870
260 0.743119266 0.336322870
261 0.743119266 0.336322870
262 0.743119266 0.336322870
263 0.743119266 0.336322870
264 0.743119266 0.336322870
265 0.733944954 0.327354260
266 0.733944954 0.327354260
267 0.733944954 0.327354260
268 0.733944954 0.327354260
269 0.733944954 0.327354260
270 0.733944954 0.327354260
271 0.733944954 0.327354260
272 0.733944954 0.327354260
273 0.733944954 0.327354260
274 0.724770642 0.304932735
275 0.724770642 0.304932735
276 0.724770642 0.304932735
277 0.724770642 0.304932735
278 0.724770642 0.304932735
279 0.724770642 0.304932735
280 0.724770642 0.304932735
281 0.724770642 0.304932735
282 0.715596330 0.300448430
283 0.715596330 0.300448430
284 0.715596330 0.300448430
285 0.715596330 0.300448430
286 0.715596330 0.300448430
287 0.715596330 0.300448430
288 0.715596330 0.300448430
289 0.715596330 0.300448430
290 0.715596330 0.300448430
291 0.688073394 0.286995516
292 0.688073394 0.286995516
293 0.688073394 0.286995516
294 0.688073394 0.286995516
295 0.688073394 0.286995516
296 0.688073394 0.286995516
297 0.688073394 0.286995516
298 0.688073394 0.286995516
299 0.688073394 0.286995516
300 0.678899083 0.273542601
301 0.678899083 0.273542601
302 0.678899083 0.273542601
303 0.678899083 0.273542601
304 0.678899083 0.273542601
305 0.678899083 0.273542601
306 0.678899083 0.273542601
307 0.678899083 0.273542601
308 0.678899083 0.273542601
309 0.678899083 0.264573991
310 0.678899083 0.264573991
311 0.678899083 0.264573991
312 0.678899083 0.264573991
313 0.678899083 0.264573991
314 0.678899083 0.264573991
315 0.678899083 0.264573991
316 0.678899083 0.264573991
317 0.678899083 0.264573991
318 0.660550459 0.255605381
319 0.660550459 0.255605381
320 0.660550459 0.255605381
321 0.660550459 0.255605381
322 0.660550459 0.255605381
323 0.660550459 0.255605381
324 0.660550459 0.255605381
325 0.660550459 0.255605381
326 0.660550459 0.255605381
327 0.660550459 0.237668161
328 0.660550459 0.237668161
329 0.660550459 0.237668161
330 0.660550459 0.237668161
331 0.660550459 0.237668161
332 0.660550459 0.237668161
333 0.660550459 0.237668161
334 0.660550459 0.237668161
335 0.660550459 0.237668161
336 0.651376147 0.224215247
337 0.651376147 0.224215247
338 0.651376147 0.224215247
339 0.651376147 0.224215247
340 0.651376147 0.224215247
341 0.651376147 0.224215247
342 0.651376147 0.224215247
343 0.651376147 0.224215247
344 0.651376147 0.224215247
345 0.651376147 0.224215247
346 0.633027523 0.215246637
347 0.633027523 0.215246637
348 0.633027523 0.215246637
349 0.633027523 0.215246637
350 0.633027523 0.215246637
351 0.633027523 0.215246637
352 0.633027523 0.215246637
353 0.633027523 0.215246637
354 0.633027523 0.215246637
355 0.633027523 0.188340807
356 0.633027523 0.188340807
357 0.633027523 0.188340807
358 0.633027523 0.188340807
359 0.633027523 0.188340807
360 0.633027523 0.188340807
361 0.633027523 0.188340807
362 0.633027523 0.188340807
363 0.633027523 0.188340807
364 0.633027523 0.188340807
365 0.633027523 0.174887892
366 0.633027523 0.174887892
367 0.633027523 0.174887892
368 0.633027523 0.174887892
369 0.633027523 0.174887892
370 0.633027523 0.174887892
371 0.633027523 0.174887892
372 0.633027523 0.174887892
373 0.633027523 0.174887892
374 0.633027523 0.174887892
375 0.596330275 0.165919283
376 0.596330275 0.165919283
377 0.596330275 0.165919283
378 0.596330275 0.165919283
379 0.596330275 0.165919283
380 0.596330275 0.165919283
381 0.596330275 0.165919283
382 0.596330275 0.165919283
383 0.596330275 0.165919283
384 0.596330275 0.165919283
385 0.568807339 0.152466368
386 0.568807339 0.152466368
387 0.568807339 0.152466368
388 0.568807339 0.152466368
389 0.568807339 0.152466368
390 0.568807339 0.152466368
391 0.568807339 0.152466368
392 0.568807339 0.152466368
393 0.568807339 0.152466368
394 0.568807339 0.152466368
395 0.568807339 0.147982063
396 0.568807339 0.147982063
397 0.568807339 0.147982063
398 0.568807339 0.147982063
399 0.568807339 0.147982063
400 0.568807339 0.147982063
401 0.568807339 0.147982063
402 0.568807339 0.147982063
403 0.568807339 0.147982063
404 0.568807339 0.147982063
405 0.568807339 0.139013453
406 0.568807339 0.139013453
407 0.568807339 0.139013453
408 0.568807339 0.139013453
409 0.568807339 0.139013453
410 0.568807339 0.139013453
411 0.568807339 0.139013453
412 0.568807339 0.139013453
413 0.568807339 0.139013453
414 0.568807339 0.139013453
415 0.568807339 0.139013453
416 0.568807339 0.134529148
417 0.568807339 0.134529148
418 0.568807339 0.134529148
419 0.568807339 0.134529148
420 0.568807339 0.134529148
421 0.568807339 0.134529148
422 0.568807339 0.134529148
423 0.568807339 0.134529148
424 0.568807339 0.134529148
425 0.568807339 0.134529148
426 0.568807339 0.130044843
427 0.568807339 0.130044843
428 0.568807339 0.130044843
429 0.568807339 0.130044843
430 0.568807339 0.130044843
431 0.568807339 0.130044843
432 0.568807339 0.130044843
433 0.568807339 0.130044843
434 0.568807339 0.130044843
435 0.568807339 0.130044843
436 0.559633028 0.121076233
437 0.559633028 0.121076233
438 0.559633028 0.121076233
439 0.559633028 0.121076233
440 0.559633028 0.121076233
441 0.559633028 0.121076233
442 0.559633028 0.121076233
443 0.559633028 0.121076233
444 0.559633028 0.121076233
445 0.559633028 0.121076233
446 0.559633028 0.121076233
447 0.550458716 0.121076233
448 0.550458716 0.121076233
449 0.550458716 0.121076233
450 0.550458716 0.121076233
451 0.550458716 0.121076233
452 0.550458716 0.121076233
453 0.550458716 0.121076233
454 0.550458716 0.121076233
455 0.550458716 0.121076233
456 0.550458716 0.121076233
457 0.522935780 0.116591928
458 0.522935780 0.116591928
459 0.522935780 0.116591928
460 0.522935780 0.116591928
461 0.522935780 0.116591928
462 0.522935780 0.116591928
463 0.522935780 0.116591928
464 0.522935780 0.116591928
465 0.522935780 0.116591928
466 0.522935780 0.116591928
467 0.522935780 0.116591928
468 0.522935780 0.112107623
469 0.522935780 0.112107623
470 0.522935780 0.112107623
471 0.522935780 0.112107623
472 0.522935780 0.112107623
473 0.522935780 0.112107623
474 0.522935780 0.112107623
475 0.522935780 0.112107623
476 0.522935780 0.112107623
477 0.522935780 0.112107623
478 0.513761468 0.112107623
479 0.513761468 0.112107623
480 0.513761468 0.112107623
481 0.513761468 0.112107623
482 0.513761468 0.112107623
483 0.513761468 0.112107623
484 0.513761468 0.112107623
485 0.513761468 0.112107623
486 0.513761468 0.112107623
487 0.513761468 0.112107623
488 0.513761468 0.112107623
489 0.513761468 0.103139013
490 0.513761468 0.103139013
491 0.513761468 0.103139013
492 0.513761468 0.103139013
493 0.513761468 0.103139013
494 0.513761468 0.103139013
495 0.513761468 0.103139013
496 0.513761468 0.103139013
497 0.513761468 0.103139013
498 0.513761468 0.103139013
499 0.513761468 0.103139013
500 0.513761468 0.103139013
501 0.513761468 0.103139013
502 0.513761468 0.103139013
503 0.513761468 0.103139013
504 0.513761468 0.103139013
505 0.513761468 0.103139013
506 0.513761468 0.103139013
507 0.513761468 0.103139013
508 0.513761468 0.103139013
509 0.513761468 0.103139013
510 0.513761468 0.098654709
511 0.513761468 0.098654709
512 0.513761468 0.098654709
513 0.513761468 0.098654709
514 0.513761468 0.098654709
515 0.513761468 0.098654709
516 0.513761468 0.098654709
517 0.513761468 0.098654709
518 0.513761468 0.098654709
519 0.513761468 0.098654709
520 0.513761468 0.098654709
521 0.504587156 0.094170404
522 0.504587156 0.094170404
523 0.504587156 0.094170404
524 0.504587156 0.094170404
525 0.504587156 0.094170404
526 0.504587156 0.094170404
527 0.504587156 0.094170404
528 0.504587156 0.094170404
529 0.504587156 0.094170404
530 0.504587156 0.094170404
531 0.504587156 0.085201794
532 0.504587156 0.085201794
533 0.504587156 0.085201794
534 0.504587156 0.085201794
535 0.504587156 0.085201794
536 0.504587156 0.085201794
537 0.504587156 0.085201794
538 0.504587156 0.085201794
539 0.504587156 0.085201794
540 0.504587156 0.085201794
541 0.504587156 0.085201794
542 0.486238532 0.076233184
543 0.486238532 0.076233184
544 0.486238532 0.076233184
545 0.486238532 0.076233184
546 0.486238532 0.076233184
547 0.486238532 0.076233184
548 0.486238532 0.076233184
549 0.486238532 0.076233184
550 0.486238532 0.076233184
551 0.486238532 0.076233184
552 0.467889908 0.076233184
553 0.467889908 0.076233184
554 0.467889908 0.076233184
555 0.467889908 0.076233184
556 0.467889908 0.076233184
557 0.467889908 0.076233184
558 0.467889908 0.076233184
559 0.467889908 0.076233184
560 0.467889908 0.076233184
561 0.467889908 0.076233184
562 0.467889908 0.076233184
563 0.458715596 0.067264574
564 0.458715596 0.067264574
565 0.458715596 0.067264574
566 0.458715596 0.067264574
567 0.458715596 0.067264574
568 0.458715596 0.067264574
569 0.458715596 0.067264574
570 0.458715596 0.067264574
571 0.458715596 0.067264574
572 0.458715596 0.067264574
573 0.449541284 0.058295964
574 0.449541284 0.058295964
575 0.449541284 0.058295964
576 0.449541284 0.058295964
577 0.449541284 0.058295964
578 0.449541284 0.058295964
579 0.449541284 0.058295964
580 0.449541284 0.058295964
581 0.449541284 0.058295964
582 0.449541284 0.058295964
583 0.449541284 0.058295964
584 0.440366972 0.058295964
585 0.440366972 0.058295964
586 0.440366972 0.058295964
587 0.440366972 0.058295964
588 0.440366972 0.058295964
589 0.440366972 0.058295964
590 0.440366972 0.058295964
591 0.440366972 0.058295964
592 0.440366972 0.058295964
593 0.440366972 0.058295964
594 0.440366972 0.058295964
595 0.440366972 0.058295964
596 0.440366972 0.058295964
597 0.440366972 0.058295964
598 0.440366972 0.058295964
599 0.440366972 0.058295964
600 0.440366972 0.058295964
601 0.440366972 0.058295964
602 0.440366972 0.058295964
603 0.440366972 0.058295964
604 0.440366972 0.053811659
605 0.440366972 0.053811659
606 0.440366972 0.053811659
607 0.440366972 0.053811659
608 0.440366972 0.053811659
609 0.440366972 0.053811659
610 0.440366972 0.053811659
611 0.440366972 0.053811659
612 0.440366972 0.053811659
613 0.440366972 0.053811659
614 0.431192661 0.044843049
615 0.431192661 0.044843049
616 0.431192661 0.044843049
617 0.431192661 0.044843049
618 0.431192661 0.044843049
619 0.431192661 0.044843049
620 0.431192661 0.044843049
621 0.431192661 0.044843049
622 0.431192661 0.044843049
623 0.431192661 0.044843049
624 0.412844037 0.040358744
625 0.412844037 0.040358744
626 0.412844037 0.040358744
627 0.412844037 0.040358744
628 0.412844037 0.040358744
629 0.412844037 0.040358744
630 0.412844037 0.040358744
631 0.412844037 0.040358744
632 0.412844037 0.040358744
633 0.412844037 0.040358744
634 0.412844037 0.031390135
635 0.412844037 0.031390135
636 0.412844037 0.031390135
637 0.412844037 0.031390135
638 0.412844037 0.031390135
639 0.412844037 0.031390135
640 0.412844037 0.031390135
641 0.412844037 0.031390135
642 0.412844037 0.031390135
643 0.412844037 0.031390135
644 0.412844037 0.026905830
645 0.412844037 0.026905830
646 0.412844037 0.026905830
647 0.412844037 0.026905830
648 0.412844037 0.026905830
649 0.412844037 0.026905830
650 0.412844037 0.026905830
651 0.412844037 0.026905830
652 0.412844037 0.026905830
653 0.412844037 0.026905830
654 0.385321101 0.026905830
655 0.385321101 0.026905830
656 0.385321101 0.026905830
657 0.385321101 0.026905830
658 0.385321101 0.026905830
659 0.385321101 0.026905830
660 0.385321101 0.026905830
661 0.385321101 0.026905830
662 0.385321101 0.026905830
663 0.376146789 0.026905830
664 0.376146789 0.026905830
665 0.376146789 0.026905830
666 0.376146789 0.026905830
667 0.376146789 0.026905830
668 0.376146789 0.026905830
669 0.376146789 0.026905830
670 0.376146789 0.026905830
671 0.376146789 0.026905830
672 0.376146789 0.026905830
673 0.376146789 0.022421525
674 0.376146789 0.022421525
675 0.376146789 0.022421525
676 0.376146789 0.022421525
677 0.376146789 0.022421525
678 0.376146789 0.022421525
679 0.376146789 0.022421525
680 0.376146789 0.022421525
681 0.376146789 0.022421525
682 0.366972477 0.022421525
683 0.366972477 0.022421525
684 0.366972477 0.022421525
685 0.366972477 0.022421525
686 0.366972477 0.022421525
687 0.366972477 0.022421525
688 0.366972477 0.022421525
689 0.366972477 0.022421525
690 0.366972477 0.022421525
691 0.366972477 0.022421525
692 0.366972477 0.022421525
693 0.366972477 0.022421525
694 0.366972477 0.022421525
695 0.366972477 0.022421525
696 0.366972477 0.022421525
697 0.366972477 0.022421525
698 0.366972477 0.022421525
699 0.366972477 0.022421525
700 0.357798165 0.022421525
701 0.357798165 0.022421525
702 0.357798165 0.022421525
703 0.357798165 0.022421525
704 0.357798165 0.022421525
705 0.357798165 0.022421525
706 0.357798165 0.022421525
707 0.357798165 0.022421525
708 0.357798165 0.022421525
709 0.348623853 0.022421525
710 0.348623853 0.022421525
711 0.348623853 0.022421525
712 0.348623853 0.022421525
713 0.348623853 0.022421525
714 0.348623853 0.022421525
715 0.348623853 0.022421525
716 0.348623853 0.022421525
717 0.348623853 0.022421525
718 0.321100917 0.022421525
719 0.321100917 0.022421525
720 0.321100917 0.022421525
721 0.321100917 0.022421525
722 0.321100917 0.022421525
723 0.321100917 0.022421525
724 0.321100917 0.022421525
725 0.321100917 0.022421525
726 0.302752294 0.022421525
727 0.302752294 0.022421525
728 0.302752294 0.022421525
729 0.302752294 0.022421525
730 0.302752294 0.022421525
731 0.302752294 0.022421525
732 0.302752294 0.022421525
733 0.302752294 0.022421525
734 0.302752294 0.022421525
735 0.302752294 0.022421525
736 0.302752294 0.022421525
737 0.302752294 0.022421525
738 0.302752294 0.022421525
739 0.302752294 0.022421525
740 0.302752294 0.022421525
741 0.302752294 0.022421525
742 0.302752294 0.022421525
743 0.302752294 0.013452915
744 0.302752294 0.013452915
745 0.302752294 0.013452915
746 0.302752294 0.013452915
747 0.302752294 0.013452915
748 0.302752294 0.013452915
749 0.302752294 0.013452915
750 0.302752294 0.013452915
751 0.293577982 0.013452915
752 0.293577982 0.013452915
753 0.293577982 0.013452915
754 0.293577982 0.013452915
755 0.293577982 0.013452915
756 0.293577982 0.013452915
757 0.293577982 0.013452915
758 0.293577982 0.013452915
759 0.293577982 0.013452915
760 0.293577982 0.013452915
761 0.293577982 0.013452915
762 0.293577982 0.013452915
763 0.293577982 0.013452915
764 0.293577982 0.013452915
765 0.293577982 0.013452915
766 0.293577982 0.013452915
767 0.293577982 0.013452915
768 0.293577982 0.013452915
769 0.293577982 0.013452915
770 0.293577982 0.013452915
771 0.293577982 0.013452915
772 0.293577982 0.013452915
773 0.293577982 0.013452915
774 0.284403670 0.013452915
775 0.284403670 0.013452915
776 0.284403670 0.013452915
777 0.284403670 0.013452915
778 0.284403670 0.013452915
779 0.284403670 0.013452915
780 0.284403670 0.013452915
781 0.266055046 0.013452915
782 0.266055046 0.013452915
783 0.266055046 0.013452915
784 0.266055046 0.013452915
785 0.266055046 0.013452915
786 0.266055046 0.013452915
787 0.266055046 0.013452915
788 0.266055046 0.013452915
789 0.256880734 0.013452915
790 0.256880734 0.013452915
791 0.256880734 0.013452915
792 0.256880734 0.013452915
793 0.256880734 0.013452915
794 0.256880734 0.013452915
795 0.256880734 0.013452915
796 0.247706422 0.013452915
797 0.247706422 0.013452915
798 0.247706422 0.013452915
799 0.247706422 0.013452915
800 0.247706422 0.013452915
801 0.247706422 0.013452915
802 0.211009174 0.008968610
803 0.211009174 0.008968610
804 0.211009174 0.008968610
805 0.211009174 0.008968610
806 0.211009174 0.008968610
807 0.211009174 0.008968610
808 0.211009174 0.008968610
809 0.192660550 0.008968610
810 0.192660550 0.008968610
811 0.192660550 0.008968610
812 0.192660550 0.008968610
813 0.192660550 0.008968610
814 0.192660550 0.008968610
815 0.192660550 0.008968610
816 0.192660550 0.008968610
817 0.192660550 0.008968610
818 0.192660550 0.008968610
819 0.192660550 0.008968610
820 0.192660550 0.008968610
821 0.192660550 0.008968610
822 0.183486239 0.008968610
823 0.183486239 0.008968610
824 0.183486239 0.008968610
825 0.183486239 0.008968610
826 0.183486239 0.008968610
827 0.183486239 0.008968610
828 0.183486239 0.008968610
829 0.183486239 0.008968610
830 0.183486239 0.008968610
831 0.183486239 0.008968610
832 0.183486239 0.008968610
833 0.183486239 0.008968610
834 0.183486239 0.008968610
835 0.183486239 0.008968610
836 0.183486239 0.008968610
837 0.183486239 0.008968610
838 0.183486239 0.008968610
839 0.183486239 0.008968610
840 0.155963303 0.008968610
841 0.155963303 0.008968610
842 0.155963303 0.008968610
843 0.155963303 0.008968610
844 0.155963303 0.008968610
845 0.155963303 0.008968610
846 0.137614679 0.004484305
847 0.137614679 0.004484305
848 0.137614679 0.004484305
849 0.137614679 0.004484305
850 0.137614679 0.004484305
851 0.100917431 0.004484305
852 0.100917431 0.004484305
853 0.100917431 0.004484305
854 0.100917431 0.004484305
855 0.100917431 0.004484305
856 0.100917431 0.004484305
857 0.100917431 0.004484305
858 0.100917431 0.004484305
859 0.100917431 0.004484305
860 0.100917431 0.004484305
861 0.100917431 0.004484305
862 0.100917431 0.004484305
863 0.100917431 0.004484305
864 0.100917431 0.004484305
865 0.100917431 0.004484305
866 0.100917431 0.004484305
867 0.091743119 0.004484305
868 0.091743119 0.004484305
869 0.091743119 0.004484305
870 0.091743119 0.004484305
871 0.091743119 0.004484305
872 0.091743119 0.004484305
873 0.091743119 0.004484305
874 0.091743119 0.004484305
875 0.091743119 0.004484305
876 0.082568807 0.004484305
877 0.082568807 0.004484305
878 0.082568807 0.004484305
879 0.082568807 0.004484305
880 0.082568807 0.004484305
881 0.055045872 0.004484305
882 0.055045872 0.004484305
883 0.055045872 0.004484305
884 0.055045872 0.004484305
885 0.055045872 0.004484305
886 0.055045872 0.004484305
887 0.055045872 0.004484305
888 0.055045872 0.004484305
889 0.055045872 0.004484305
890 0.036697248 0.004484305
891 0.036697248 0.004484305
892 0.036697248 0.004484305
893 0.036697248 0.004484305
894 0.036697248 0.004484305
895 0.036697248 0.004484305
896 0.036697248 0.004484305
897 0.036697248 0.004484305
898 0.036697248 0.004484305
899 0.036697248 0.004484305
900 0.036697248 0.004484305
901 0.036697248 0.004484305
902 0.036697248 0.004484305
903 0.036697248 0.004484305
904 0.036697248 0.004484305
905 0.027522936 0.004484305
906 0.027522936 0.004484305
907 0.027522936 0.004484305
908 0.027522936 0.004484305
909 0.027522936 0.004484305
910 0.027522936 0.004484305
911 0.027522936 0.004484305
912 0.027522936 0.004484305
913 0.027522936 0.004484305
914 0.027522936 0.004484305
915 0.027522936 0.004484305
916 0.009174312 0.004484305
917 0.009174312 0.004484305
918 0.009174312 0.004484305
919 0.000000000 0.000000000
920 0.000000000 0.000000000
921 0.000000000 0.000000000
922 0.000000000 0.000000000
923 0.000000000 0.000000000
924 0.000000000 0.000000000
925 0.000000000 0.000000000
926 0.000000000 0.000000000
927 0.000000000 0.000000000
928 0.000000000 0.000000000
929 0.000000000 0.000000000
930 0.000000000 0.000000000
931 0.000000000 0.000000000
932 0.000000000 0.000000000
933 0.000000000 0.000000000
934 0.000000000 0.000000000
935 0.000000000 0.000000000
936 0.000000000 0.000000000
937 0.000000000 0.000000000
938 0.000000000 0.000000000
939 0.000000000 0.000000000
940 0.000000000 0.000000000
941 0.000000000 0.000000000
942 0.000000000 0.000000000
943 0.000000000 0.000000000
944 0.000000000 0.000000000
945 0.000000000 0.000000000
946 0.000000000 0.000000000
947 0.000000000 0.000000000
948 0.000000000 0.000000000
949 0.000000000 0.000000000
950 0.000000000 0.000000000
951 0.000000000 0.000000000
952 0.000000000 0.000000000
953 0.000000000 0.000000000
954 0.000000000 0.000000000
955 0.000000000 0.000000000
956 0.000000000 0.000000000
957 0.000000000 0.000000000
958 0.000000000 0.000000000
959 0.000000000 0.000000000
960 0.000000000 0.000000000
961 0.000000000 0.000000000
962 0.000000000 0.000000000
963 0.000000000 0.000000000
964 0.000000000 0.000000000
965 0.000000000 0.000000000
966 0.000000000 0.000000000
967 0.000000000 0.000000000
968 0.000000000 0.000000000
969 0.000000000 0.000000000
970 0.000000000 0.000000000
971 0.000000000 0.000000000
972 0.000000000 0.000000000
973 0.000000000 0.000000000
974 0.000000000 0.000000000
975 0.000000000 0.000000000
976 0.000000000 0.000000000
977 0.000000000 0.000000000
978 0.000000000 0.000000000
979 0.000000000 0.000000000
980 0.000000000 0.000000000
981 0.000000000 0.000000000
982 0.000000000 0.000000000
983 0.000000000 0.000000000
984 0.000000000 0.000000000
985 0.000000000 0.000000000
986 0.000000000 0.000000000
987 0.000000000 0.000000000
988 0.000000000 0.000000000
989 0.000000000 0.000000000
990 0.000000000 0.000000000
991 0.000000000 0.000000000
992 0.000000000 0.000000000
993 0.000000000 0.000000000
994 0.000000000 0.000000000
995 0.000000000 0.000000000
996 0.000000000 0.000000000
997 0.000000000 0.000000000
998 0.000000000 0.000000000
999 0.000000000 0.000000000
1000 0.000000000 0.000000000
1001 0.000000000 0.000000000
A scatterplot of these rates is called the ROC curve. We start at (0,0) and go to (1,1)
plot (c (1 ,FPR,0 ),c (1 , TPR,0 ), type= "s" , xlab= "False Positive Rate" , ylab= "True Positive Rate" )
The PRROC package produces pretty ROC curves.
#install.packages("PRROC") library (PRROC)
Loading required package: rlang
<- roc.curve (scores.class0 = predict (pima_logit, type= "response" ), weights.class0= Pima.te$ typeTF,curve= TRUE )plot (PRROC_obj)
The area under the curve (AUC) is a common metric measuring how good your classifier is. An AUC = 0.5 (where the ROC curve is a diagonal line) is a “dumb” predictor. AUC closer to 1 represents a discerning “good” predictor.
Comparing using AIC and BIC
As we’ll see next week, the AIC is a useful statistic to compare logistic models. \[AIC = \text{residual deviance} + 2k\text{ (k is the number of coefficients fit)}\] \[BIC = \text{residual deviance} + ln(n)*k\] When sample size is large then BIC will give a larger penalty per predictor. \(ln(8) > 2\) so when \(n \geq 8\) BIC gives a heftier penalty per coefficient in the model than AIC does and thus will favor simpler models.
When using either of these comparative statistics, the lower the value the better.
<- glm (typeTF ~ 1 + glu,data= Pima.te, family= "binomial" )<- glm (typeTF ~ 1 + glu + bmi, data= Pima.te, family= "binomial" )summary (pima.logit1)
Call:
glm(formula = typeTF ~ 1 + glu, family = "binomial", data = Pima.te)
Deviance Residuals:
Min 1Q Median 3Q Max
-2.2343 -0.7270 -0.4985 0.6663 2.3268
Coefficients:
Estimate Std. Error z value Pr(>|z|)
(Intercept) -5.946808 0.659839 -9.013 <2e-16 ***
glu 0.042421 0.005165 8.213 <2e-16 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
(Dispersion parameter for binomial family taken to be 1)
Null deviance: 420.30 on 331 degrees of freedom
Residual deviance: 325.99 on 330 degrees of freedom
AIC: 329.99
Number of Fisher Scoring iterations: 4
Call:
glm(formula = typeTF ~ 1 + glu + bmi, family = "binomial", data = Pima.te)
Deviance Residuals:
Min 1Q Median 3Q Max
-2.2240 -0.7216 -0.4421 0.6158 2.3544
Coefficients:
Estimate Std. Error z value Pr(>|z|)
(Intercept) -8.025302 0.942896 -8.511 < 2e-16 ***
glu 0.039311 0.005257 7.478 7.57e-14 ***
bmi 0.072645 0.020621 3.523 0.000427 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
(Dispersion parameter for binomial family taken to be 1)
Null deviance: 420.30 on 331 degrees of freedom
Residual deviance: 312.49 on 329 degrees of freedom
AIC: 318.49
Number of Fisher Scoring iterations: 4
Example 2: A Logistic model with a categorical predictor
#Let's use age as a predictor #Only for age > 50 $ age50 <- as.numeric (Pima.te$ age > 50 )<- aggregate (typeTF ~ age50, data= Pima.te, FUN= mean)
age50 typeTF
1 0 0.3129032
2 1 0.5454545
Among those 50 or younger, proportion with diabetes is .3129032 among those 51+ it is .545454
#Log-odds for the two groups <- diabetes.aggr$ typeTF / (1 - diabetes.aggr$ typeTF)<- log (odds)
Notice the log-odds are -.7865 and .1823 respectively. If we were to build a model we would say that \[ \hat{LO} = -.7865812 + .9689 \times age_{>50}\]
<- glm (typeTF ~ age50, data= Pima.te, family= "binomial" )summary (Pima50)
Call:
glm(formula = typeTF ~ age50, family = "binomial", data = Pima.te)
Deviance Residuals:
Min 1Q Median 3Q Max
-1.2557 -0.8663 -0.8663 1.5244 1.5244
Coefficients:
Estimate Std. Error z value Pr(>|z|)
(Intercept) -0.7866 0.1225 -6.422 1.35e-10 ***
age50 0.9689 0.4454 2.176 0.0296 *
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
(Dispersion parameter for binomial family taken to be 1)
Null deviance: 420.30 on 331 degrees of freedom
Residual deviance: 415.59 on 330 degrees of freedom
AIC: 419.59
Number of Fisher Scoring iterations: 4
Example 3: Iris dataset
Another classic example is the iris dataset. This famous dataset gives measurements in centimeters of the sepal and petal width and length for 3 species of iris: setosa, versicolor and virginica. The dataset has 150 rows, with 50 samples of each.
Because the logistic model handles binary response variables, what we will do is create a logistic model to predict whether an iris is virginica Let’s build a full model with all 4 predictors.
<- glm (Species== "virginica" ~ 1 + Sepal.Length + Sepal.Width+ Petal.Length+ Petal.Width, data= iris, family= "binomial" )
Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
Call:
glm(formula = Species == "virginica" ~ 1 + Sepal.Length + Sepal.Width +
Petal.Length + Petal.Width, family = "binomial", data = iris)
Deviance Residuals:
Min 1Q Median 3Q Max
-2.01105 -0.00065 0.00000 0.00048 1.78065
Coefficients:
Estimate Std. Error z value Pr(>|z|)
(Intercept) -42.638 25.708 -1.659 0.0972 .
Sepal.Length -2.465 2.394 -1.030 0.3032
Sepal.Width -6.681 4.480 -1.491 0.1359
Petal.Length 9.429 4.737 1.990 0.0465 *
Petal.Width 18.286 9.743 1.877 0.0605 .
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
(Dispersion parameter for binomial family taken to be 1)
Null deviance: 190.954 on 149 degrees of freedom
Residual deviance: 11.899 on 145 degrees of freedom
AIC: 21.899
Number of Fisher Scoring iterations: 12
plot (predict (virginica.model), iris$ Species== "virginica" , col= iris$ Species)<- range (predict (virginica.model))<- seq (rng[1 ],rng[2 ], length.out= 100 )<- 1 / (1 + exp (- xs))lines (xs,ys)
Let’s just contrast this plot with a model that only uses sepal length.
<- glm (Species== "virginica" ~ 1 + Sepal.Length, data= iris, family= "binomial" )summary (virginica.model.simple)
Call:
glm(formula = Species == "virginica" ~ 1 + Sepal.Length, family = "binomial",
data = iris)
Deviance Residuals:
Min 1Q Median 3Q Max
-1.9870 -0.5520 -0.2614 0.5832 2.7001
Coefficients:
Estimate Std. Error z value Pr(>|z|)
(Intercept) -16.3198 2.6581 -6.140 8.27e-10 ***
Sepal.Length 2.5921 0.4316 6.006 1.90e-09 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
(Dispersion parameter for binomial family taken to be 1)
Null deviance: 190.95 on 149 degrees of freedom
Residual deviance: 117.35 on 148 degrees of freedom
AIC: 121.35
Number of Fisher Scoring iterations: 5
plot (jitter (predict (virginica.model.simple)), iris$ Species== "virginica" , col= iris$ Species)<- range (predict (virginica.model.simple))<- seq (rng[1 ],rng[2 ], length.out= 100 )<- 1 / (1 + exp (- xs))lines (xs,ys)
How well we did? TRUE represents virginica
table (predict = predict (virginica.model, type= "response" )> .5 , actual = iris$ Species== "virginica" )
actual
predict FALSE TRUE
FALSE 99 1
TRUE 1 49
column prop table will show us type 1&2 error rates and power:
prop.table (table (predict = predict (virginica.model, type= "response" )> .5 , actual = iris$ Species== "virginica" ), margin= 2 )
actual
predict FALSE TRUE
FALSE 0.99 0.02
TRUE 0.01 0.98
Can we achieve a 5% Type 1 error rate? Adjust the threshold. We can find the 5th percentile for the predictions for the “other”
quantile (predict (virginica.model, type= "response" )[iris$ Species!= "virginica" ], .95 )
prop.table (table (predict = predict (virginica.model, type= "response" )> 0.00504067 , actual = iris$ Species== "virginica" ), margin= 2 )
actual
predict FALSE TRUE
FALSE 0.95 0.00
TRUE 0.05 1.00
We require very little evidence to claim the flower is virginica, a probability of 0.00504 is enough. In this case we have a 5% type 1 error rate. But note that the power has increased to 100%.
When a logistic model cannot converge
If your model is so good that it correctly labels everyhing as 1 or zero this is what the output looks like
<- glm (Species== "setosa" ~ 1 + Sepal.Length + Sepal.Width+ Petal.Length+ Petal.Width, data= iris, family= "binomial" )
Warning: glm.fit: algorithm did not converge
Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
Call:
glm(formula = Species == "setosa" ~ 1 + Sepal.Length + Sepal.Width +
Petal.Length + Petal.Width, family = "binomial", data = iris)
Deviance Residuals:
Min 1Q Median 3Q Max
-3.185e-05 -2.100e-08 -2.100e-08 2.100e-08 3.173e-05
Coefficients:
Estimate Std. Error z value Pr(>|z|)
(Intercept) -16.946 457457.106 0 1
Sepal.Length 11.759 130504.043 0 1
Sepal.Width 7.842 59415.386 0 1
Petal.Length -20.088 107724.596 0 1
Petal.Width -21.608 154350.619 0 1
(Dispersion parameter for binomial family taken to be 1)
Null deviance: 1.9095e+02 on 149 degrees of freedom
Residual deviance: 3.2940e-09 on 145 degrees of freedom
AIC: 10
Number of Fisher Scoring iterations: 25
It means that the coefficients are able to separate the predicted log odds so far that you have no misclassifications. Note above that after 25 fisher iterations it stopped. You can set the max iterations to be a lot lower, like 3, and then the model summary will be far more interpretable. But the standard errors of the coefficients will still be high so you can’t say much about the significance of individual predictors.
<- glm (Species== "setosa" ~ 1 + Sepal.Length + Sepal.Width+ Petal.Length+ Petal.Width, data= iris, family= "binomial" , control = list (maxit = 2 ))
Warning: glm.fit: algorithm did not converge
summary (setosa.model.maxit)
Call:
glm(formula = Species == "setosa" ~ 1 + Sepal.Length + Sepal.Width +
Petal.Length + Petal.Width, family = "binomial", data = iris,
control = list(maxit = 2))
Deviance Residuals:
Min 1Q Median 3Q Max
-0.7003 -0.2994 -0.1592 0.2188 0.7140
Coefficients:
Estimate Std. Error z value Pr(>|z|)
(Intercept) -2.4428 3.1506 -0.775 0.4381
Sepal.Length 0.4038 0.8871 0.455 0.6490
Sepal.Width 1.8157 0.8944 2.030 0.0424 *
Petal.Length -1.7249 0.8993 -1.918 0.0551 .
Petal.Width -0.2378 1.5787 -0.151 0.8803
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
(Dispersion parameter for binomial family taken to be 1)
Null deviance: 190.954 on 149 degrees of freedom
Residual deviance: 14.313 on 145 degrees of freedom
AIC: 24.313
Number of Fisher Scoring iterations: 2
How correlated are flower measurements?
Sepal.Length Sepal.Width Petal.Length Petal.Width
Sepal.Length 1.0000000 -0.1175698 0.8717538 0.8179411
Sepal.Width -0.1175698 1.0000000 -0.4284401 -0.3661259
Petal.Length 0.8717538 -0.4284401 1.0000000 0.9628654
Petal.Width 0.8179411 -0.3661259 0.9628654 1.0000000
