Survival on the titanic, using loglm()

Exercises: survival on the titanic, using loglm()

library(MASS)     # for loglm()
library(vcd)      # for mosaic, aka plot.loglm()

data(Titanic, package="datasets")  # effects::Titanic gives a case-form version

Baseline model

Titanic <- Titanic + 0.5   # adjust for 0 cells
titanic.mod1 <- loglm(~ (Class * Age * Sex) + Survived, data=Titanic)
titanic.mod1

## Call:
## loglm(formula = ~(Class * Age * Sex) + Survived, data = Titanic)
## 
## Statistics:
##                  X^2 df P(> X^2)
## Likelihood Ratio 659 15        0
## Pearson          643 15        0

plot(titanic.mod1, main="Model [AGC][S]")

Associations with survival

titanic.mod2 <- loglm(~ (Class * Age * Sex) + Survived*(Class + Age + Sex), data=Titanic)
titanic.mod2

## Call:
## loglm(formula = ~(Class * Age * Sex) + Survived * (Class + Age + 
##     Sex), data = Titanic)
## 
## Statistics:
##                  X^2 df P(> X^2)
## Likelihood Ratio 106 10        0
## Pearson          101 10        0

plot(titanic.mod2,  main="Model [AGC][AS][GS][CS]")

update() is a simpler way:

update(titanic.mod1, . ~ .+ Survived*(Class+Age+Sex))

## Call:
## loglm(formula = . ~ Class + Age + Sex + Survived + Class:Age + 
##     Class:Sex + Age:Sex + Class:Survived + Age:Survived + Sex:Survived + 
##     Class:Age:Sex, data = Titanic)
## 
## Statistics:
##                  X^2 df P(> X^2)
## Likelihood Ratio 106 10        0
## Pearson          101 10        0

titanic.mod3 <- loglm(~ (Class * Age * Sex) + Survived*(Class + Age * Sex), data=Titanic)
titanic.mod3

## Call:
## loglm(formula = ~(Class * Age * Sex) + Survived * (Class + Age * 
##     Sex), data = Titanic)
## 
## Statistics:
##                   X^2 df P(> X^2)
## Likelihood Ratio 85.5  9 1.29e-14
## Pearson          78.8  9 2.76e-13

plot(titanic.mod3, main="Model [AGC][AS][GS][CS][AGS]")

compare models

anova(titanic.mod1, titanic.mod2, titanic.mod3, test="chisq")

## LR tests for hierarchical log-linear models
## 
## Model 1:
##  ~(Class * Age * Sex) + Survived 
## Model 2:
##  ~(Class * Age * Sex) + Survived * (Class + Age + Sex) 
## Model 3:
##  ~(Class * Age * Sex) + Survived * (Class + Age * Sex) 
## 
##           Deviance df Delta(Dev) Delta(df) P(> Delta(Dev)
## Model 1      659.3 15                                    
## Model 2      105.6 10      553.8         5          0e+00
## Model 3       85.5  9       20.0         1          1e-05
## Saturated      0.0  0       85.5         9          0e+00