###################################################
### chunk number 1: preliminaries
###################################################
options(width = 70, show.signif.stars = FALSE)
data(Contraception, package = "mlmRev")
library(lattice)
library(Matrix)
library(lme4)
#lattice.options(default.theme = standard.theme(color = FALSE))


###################################################
### chunk number 2: Contra1
###################################################
print(xyplot(ifelse(use == "Y", 1, 0) ~ age|urban, Contraception,
             groups = livch, type = c("g", "smooth"),
             auto.key = list(space = "top", points = FALSE,
             lines = TRUE, columns = 4),
             ylab = "Proportion", xlab = "Centered age"))


###################################################
### chunk number 3: cm1
###################################################
print(cm1 <- glmer(use ~ age + I(age^2) + urban + livch + (1|district), 
                  Contraception, binomial), corr = FALSE)


###################################################
### chunk number 4: Contra2
###################################################
Contraception$ch <- factor(Contraception$livch != 0, labels = c("N","Y"))
print(xyplot(ifelse(use == "Y", 1, 0) ~ age|urban, Contraception,
             groups = ch, type = c("g", "smooth"),
             auto.key = list(space = "top", points = FALSE,
             lines = TRUE, columns = 2),
             ylab = "Proportion", xlab = "Centered age"))


###################################################
### chunk number 5: cm2
###################################################
print(cm2 <- glmer(use ~ age + I(age^2) + urban + ch + (1|district),
                  Contraception, binomial), corr = FALSE)


###################################################
### chunk number 6: anovac
###################################################
anova(cm2, cm1)


###################################################
### chunk number 7: cm3
###################################################
print(cm3 <- glmer(use ~ age*ch + I(age^2) + urban + (1|district),
              Contraception, binomial), corr = FALSE)


###################################################
### chunk number 8: ContraCat
###################################################
print(dotplot(ranef(cm3, post = TRUE),strip=FALSE)$district)


###################################################
### chunk number 9: urbanRural
###################################################
cat(head(capture.output(xtabs(~urban+district, Contraception)),7),sep='\n')


###################################################
### chunk number 10: cm4
###################################################
print(cm4 <- glmer(use ~ age*ch + I(age^2) + urban + (urban|district),
                   Contraception, binomial), corr = FALSE)


###################################################
### chunk number 11: anovacm4
###################################################
anova(cm4,cm3)


###################################################
### chunk number 12: ContraCat2
###################################################
print(dotplot(ranef(cm4, post = TRUE))[["district"]])


###################################################
### chunk number 13: ContraSc
###################################################
print(plot(ranef(cm4), type = c("g","p"), aspect = 1)[["district"]])