#
                                        # Binomial GLM
                                        #
                                        # Snoring data
                                        # 

snoring <- read.table('snoring.dat')
names(snoring) <- c('snoring','disease', 'total')
snoring.x <- snoring$snoring
snoring.y <- cbind(ndisease=snoring$disease,
                   nnormal=snoring$total-snoring$disease)
snoring.glm <- glm(snoring.y ~ snoring.x, family=binomial)
snoring.glm.probit <- glm(snoring.y ~ snoring.x,
                   family=binomial(link='log'))
plot(snoring$snoring, snoring$disease/snoring$total)

xx <- seq(0,5, length=100)
lines(xx,
      predict(snoring.glm, newdata=data.frame(snoring.x=xx), type='response'),
      col='red')
lines(xx,
      predict(snoring.glm.probit, newdata=data.frame(snoring.x=xx),
              type='response'),
      col='blue')


                                        #
                                        # Poisson GLM
                                        # 
                                        #
                                        # crab data
                                        # 
crab <- read.table('crab.dat')
names(crab) <- c('color','spine','width','satell','weight')
plot(crab)

crab.glm <- glm(satell ~ width, data=crab, family=poisson(link="log"))
plot(crab.glm)

                                        #
                                        # Poisson GLM
                                        # with the log link function
                                        # (the canonlical link)
                                        # 
crab.glm <- glm(satell ~ width, data=crab, family=poisson(link="log"))
                                        #
                                        # Poisson GLM
                                        # with the identity link function
                                        # 
crab.glm.id <- glm(satell ~ width, data=crab,
                   family=poisson(link="identity"), start=c(1,0))
                                        # 
                                        # not specifying the start value
                                        # produces an error
                                        # 

plot(crab$width, crab$satell)
xx <- seq(min(crab$width), max(crab$width), length=100)
lines(xx,
       predict(crab.glm, newdata=data.frame(width=xx), type='response') ,
       col='red')
lines(xx,
       predict(crab.glm.id, newdata=data.frame(width=xx), type='response') ,
       col='blue')




                                        #
crab$binary <- crab$satell>0
plot(crab$width, crab$binary)
boxplot(split(crab$width, crab$binary), col='gray')
boxplot(split(crab$width, crab$binary), horizontal=TRUE, col='gray')

plot(crab$width, crab$binary)

crab2.glm <- glm(binary  ~ width, data=crab, family=binomial(link='logit'))
plot(crab$width, crab$binary)
xx <- seq(min(crab$width), max(crab$width), n=100)
lines(xx,
      predict(crab2.glm, newdata=data.frame(width=xx),
              type='response'),
      col='blue')

library(mgcv)
crab2.gam <- gam(binary  ~ s(width), data=crab, family=binomial(link='logit'))
lines(xx,
      predict(crab2.gam, newdata=data.frame(width=xx),
              type='response'),
      col='blue')