## Chapter 8: {G}eneralizedLinear Models, and Survival Analysis 

## Sec 8.1: Generalized Linear Models 
## ss 8.1.1: Transformation of the expected value on the left 
## Footnote Code
## Simplified plot showing the logit link function 
p <- (1:999)/1000 
gitp <- log(p/(1 - p)) 
par(pty="s")
plot(p, gitp, xlab = "Proportion", ylab = "", type = "l", pch = 1)

plot of chunk unnamed-chunk-1

par(pty="m")

## ss 8.1.2: Noise terms need not be normal 
## ss 8.1.3: Log odds in contingency tables 
## ss 8.1.4: Logistic regression with a continuous explanatory variable 
library(DAAG) 
## Loading required package: lattice
anestot <- aggregate(anesthetic[, c("move","nomove")],  
                     by=list(conc=anesthetic$conc), FUN=sum) 
## The column 'conc', because from the 'by' list, is then a factor. 
## The next line recovers the numeric values 
anestot$conc <- as.numeric(as.character(anestot$conc)) 
anestot$total <- apply(anestot[, c("move","nomove")], 1 , sum) 
anestot$prop <- anestot$nomove/anestot$total 

## Footnote Code
## Plot proportion moving vs conc: data frame anesthetic (DAAG) 
plot(prop ~ conc, data=anestot, xlab = "Concentration",  
     ylab = "Proportion", xlim = c(.5,2.5), ylim = c(0, 1), pch = 16) 
with(anestot, 
     {text(conc, prop, paste(total), pos=2) 
      abline(h=sum(nomove)/sum(total), lty=2)}) 

plot of chunk unnamed-chunk-1

## Fit model directly to the 0/1 data in nomove 
anes.logit <- glm(nomove ~ conc, family=binomial(link="logit"), 
                  data=anesthetic) 
## Fit model to the proportions; supply total numbers as weights 
anes1.logit <- glm(prop ~ conc, family=binomial(link="logit"), 
                   weights=total, data=anestot) 

## Footnote Code
## Graphical summary of logistic regression results 
anestot$emplogit <- with(anestot, log((nomove+0.5)/(move+0.5))) 
plot(emplogit ~ conc, data=anestot, 
     xlab = "Concentration", xlim = c(0, 2.75), xaxs="i", 
     ylab = "Empirical logit", ylim=c(-2, 2.4), cex=1.5, pch = 16) 
with(anestot, text(conc, emplogit, paste(round(prop,2)), pos=c(2,4,2,4,4,4))) 
abline(anes.logit) 

plot of chunk unnamed-chunk-1

summary(anes.logit) 
## 
## Call:
## glm(formula = nomove ~ conc, family = binomial(link = "logit"), 
##     data = anesthetic)
## 
## Deviance Residuals: 
##     Min       1Q   Median       3Q      Max  
## -1.7667  -0.7441   0.0341   0.6867   2.0690  
## 
## Coefficients:
##             Estimate Std. Error z value Pr(>|z|)
## (Intercept)    -6.47       2.42   -2.67   0.0075
## conc            5.57       2.04    2.72   0.0064
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 41.455  on 29  degrees of freedom
## Residual deviance: 27.754  on 28  degrees of freedom
## AIC: 31.75
## 
## Number of Fisher Scoring iterations: 5
## Sec 8.2: Logistic Multiple Regression 
## Footnote Code
## Presence/absence information: data frame frogs (DAAGS) 
plot(northing ~ easting, data=frogs, pch=c(1,16)[frogs$pres.abs+1], 
     xlab="Meters east of reference point", ylab="Meters north", asp=1) 

plot of chunk unnamed-chunk-1

## Footnote Code
## Pairs plot; frogs data 
pairs(frogs[, c(5:10,4)], oma=c(2,2,2,2), cex=0.5)