#################################################################################
# Bayesian model selection and statistical modeling 
# Chapman & Hall/CRC Taylor and Francis Group 
#
# Chapter5: Section5.5.2
#
# Link function selection for binomial regression
#
# Author : Tomohiro Ando
#
################################################################################

# O-ring data

y <- c(1,1,1,1,0,0,0,0,0,0,0,0,1,1,0,0,0,1,0,0,0,0,0)

x <- c(53,57,58,63,66,67,67,67,68,69,70,70,70,70,72,73,75,75,76,76,78,79,81)

n <- length(y)


# Fitting logit model 
p <- glm(y~x,family=binomial(link="logit"))
B <- p$coefficients
O <- cbind(1,x)%*%B
P1 <- exp(O)/(1+exp(O))

BIC1 <- p$aic-2+log(n)*2

summary(p)


# Fitting probit model 


p <- glm(y~x,family=binomial(link="probit"))
B <- p$coefficients
O <- cbind(1,x)%*%B
P2 <- pnorm(O)
BIC2 <- p$aic-2+log(n)*2

summary(p)


# Plot

y <- y/6

plot(x,y,xlab="Temperature",ylab="Probability of failure of the O-rings",ylim=c(0,1),xlim=c(min(x),max(x)))
lines(x,P1)
lines(x,P2,lty=2)