Output from titanic-loglin.R

# Exercises: survival on the titanic, using loglm()

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

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.3 15        0
## Pearson          643.2 15        0

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

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 105.6 10        0
## Pearson          101.4 10        0

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

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.52  9 1.288e-14
## Pearson          78.83  9 2.756e-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.32 15                                    
## Model 2     105.55 10     553.76         5          0e+00
## Model 3      85.52  9      20.04         1          1e-05
## Saturated     0.00  0      85.52         9          0e+00