## ----echo=FALSE,eval=FALSE----------------------------------------------------
#  options(width=80)

## ----results='hide',eval=FALSE------------------------------------------------
#  library(catdata)
#  data(addiction)
#  attach(addiction)

## ----eval=FALSE---------------------------------------------------------------
#  library(mgcv)

## ----eval=FALSE---------------------------------------------------------------
#  minage <- min(na.omit(age))
#  maxage <- max(na.omit(age))
#  ageindex <- seq(minage, maxage, 0.1)
#  n <- length(ageindex)

## ----eval=FALSE---------------------------------------------------------------
#  gender1 <- rep(1, n)
#  gender0 <- rep(0, n)
#  university1 <- rep(1, n)
#  university0 <- rep(0, n)
#  
#  datamale <- as.data.frame(cbind(gender=gender0,age=ageindex,university=university1))
#  datafemale <- as.data.frame(cbind(gender=gender1,age=ageindex,
#                                    university=university1))

## ----eval=FALSE---------------------------------------------------------------
#  ill01 <-ill
#  ill01[ill==1] <- 0
#  ill01[ill==2] <- 1

## ----eval=FALSE---------------------------------------------------------------
#  gam1 <- gam(as.factor(ill01) ~ s(age) + gender + university, family=binomial())
#  gam2 <- gam(as.factor(ill) ~ s(age) + gender + university, family=binomial(),
#              data=addiction[addiction$ill!=2,])

## ----eval=FALSE---------------------------------------------------------------
#  summary(gam1)
#  summary(gam2)

## ----fig.width=9,eval=FALSE---------------------------------------------------
#  plot(gam1)

## ----fig.width=9,eval=FALSE---------------------------------------------------
#  plot(gam2)

## ----eval=FALSE---------------------------------------------------------------
#  predmale1 <- predict(gam1, newdata=datamale, type="response")
#  predmale2 <- predict(gam2, newdata=datamale, type="response")

## ----eval=FALSE---------------------------------------------------------------
#  predfemale1 <- predict(gam1, newdata=datafemale, type="response")
#  predfemale2 <- predict(gam2, newdata=datafemale, type="response")

## ----eval=FALSE---------------------------------------------------------------
#  p2 <- predmale1
#  p1 <- predmale2 * (1-predmale1)
#  p0 <- (1-predmale2) * (1-predmale1)

## ----eval=FALSE---------------------------------------------------------------
#  pf2 <- predfemale1
#  pf1 <- predfemale2 * (1-predfemale1)
#  pf0 <- (1-predfemale2) * (1-predfemale1)

## ----fig.width=16, fig.height=8.5,eval=FALSE----------------------------------
#  par(mfrow=c(1,2), cex=1.8)
#  plot(ageindex, p0, type="l", lty=1, ylim=c(0,1), main="men with university degree",
#       ylab="probabilities")
#  lines(ageindex, p1, lty="dotted")
#  lines(ageindex, p2, lty="dashed")
#  legend("topright", legend=c("Weak-willed", "diseased", "both"),
#         lty=c("solid", "dotted", "dashed"))
#  
#  
#  plot(ageindex, pf0, type="l", lty=1, ylim=c(0,1),
#       main="women with university degree", ylab="probabilities")
#  lines(ageindex, pf1, lty="dotted")
#  lines(ageindex, pf2, lty="dashed")
#  legend("topright", legend=c("Weak-willed", "diseased", "both"),
#         lty=c("solid", "dotted", "dashed"))

## ----eval=FALSE---------------------------------------------------------------
#  gam3 <- gam(as.factor(ill)~ s(age) + gender + university,
#              data=addiction[addiction$ill!=2,], family=binomial())
#  gam4 <- gam(as.factor(ill)~ s(age) + gender + university,
#              data=addiction[addiction$ill!=1,], family=binomial())

## ----eval=FALSE---------------------------------------------------------------
#  summary(gam3)
#  summary(gam4)