Plot for Predicted Probability

그래프를 그리기 위해 먼저 안선생님 책 제 18장에 있는 logistic regression을 잠깐 review 하겠읍니다.

Logistic regression test 를 실행하기 위한 dataframe을 만듭니다.

dose=c(0,0,1,1,2,2)
response=c(0,1,0,1,0,1)
count=c(7,3,5,5,2,8)
toxic=data.frame(dose,response,count)
toxic

##   dose response count
## 1    0        0     7
## 2    0        1     3
## 3    1        0     5
## 4    1        1     5
## 5    2        0     2
## 6    2        1     8

glm을 이용하여 Logistic regression test 를 실행합니다. 이 내용은 안선생님 책 179 쪽에 있는 것 과 동일합니다.

out=glm(response~dose,weights=count,family=binomial,data=toxic)
summary(out)

## 
## Call:
## glm(formula = response ~ dose, family = binomial, data = toxic, 
##     weights = count)
## 
## Deviance Residuals: 
##     1      2      3      4      5      6  
## -2.14   2.76  -2.79   2.48  -2.46   1.99  
## 
## Coefficients:
##             Estimate Std. Error z value Pr(>|z|)  
## (Intercept)   -0.945      0.636   -1.49    0.138  
## dose           1.105      0.519    2.13    0.033 *
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 41.455  on 5  degrees of freedom
## Residual deviance: 36.196  on 4  degrees of freedom
## AIC: 40.2
## 
## Number of Fisher Scoring iterations: 4

exp(coef(out))

## (Intercept)        dose 
##      0.3889      3.0194

exp(confint(out))

## Waiting for profiling to be done...

##              2.5 % 97.5 %
## (Intercept) 0.0985  1.271
## dose        1.1660  9.339

dose의 odd ratio는 3.02입니다. 즉 dose가 1 증가하면 죽을 가능성(response가 1일 가능성)이 3.02배 증가합니다. 이 모델을 시각적으로 나타내기 위해 predicted probability 그래프를 그려봅니다.

plot(response~dose,data=toxic)
curve(predict(out,data.frame(dose=x),type="resp"),add=TRUE)

여기서 문제 입니다. 다음과 같이 predicted probability 의 95% 신뢰구간을 같이 표시하고 싶습니다. 제일 먼저 정답을 맞추신 분께 치맥 시원하게 쏩니다.