Output from mental-glm.R

## Mental health data
library(gnm)
library(vcdExtra)
data(Mental)

# display the frequency table
(Mental.tab <- xtabs(Freq ~ mental + ses, data = Mental))

##           ses
## mental       1   2   3   4   5   6
##   Well      64  57  57  72  36  21
##   Mild      94  94 105 141  97  71
##   Moderate  58  54  65  77  54  54
##   Impaired  46  40  60  94  78  71

# fit independence model Residual deviance: 47.418 on 15 degrees of freedom
indep <- glm(Freq ~ mental + ses, family = poisson, data = Mental)
deviance(indep)

## [1] 47.42


long.labels <- list(set_varnames = c(mental = "Mental Health Status", ses = "Parent SES"))
mosaic(indep, residuals_type = "rstandard", labeling_args = long.labels, labeling = labeling_residuals, 
    main = "Mental health data: Independence")

## Warning: no formula provided, assuming ~ses + mental

# as a sieve diagram
mosaic(indep, labeling_args = long.labels, panel = sieve, gp = shading_Friendly, main = "Mental health data: Independence")

## Warning: no formula provided, assuming ~ses + mental

# fit linear x linear (uniform) association.  Use integer scores for rows/cols
Cscore <- as.numeric(Mental$ses)
Rscore <- as.numeric(Mental$mental)

# column effects model (ses)
coleff <- glm(Freq ~ mental + ses + Rscore:ses, family = poisson, data = Mental)
mosaic(coleff, residuals_type = "rstandard", labeling_args = long.labels, labeling = labeling_residuals, 
    suppress = 1, gp = shading_Friendly, main = "Mental health data: Col effects (ses)")

## Warning: no formula provided, assuming ~ses + mental

# row effects model (mental)
roweff <- glm(Freq ~ mental + ses + mental:Cscore, family = poisson, data = Mental)
mosaic(roweff, residuals_type = "rstandard", labeling_args = long.labels, labeling = labeling_residuals, 
    suppress = 1, gp = shading_Friendly, main = "Mental health data: Row effects (mental)")

## Warning: no formula provided, assuming ~ses + mental

linlin <- glm(Freq ~ mental + ses + Rscore:Cscore, family = poisson, data = Mental)

# compare models
anova(indep, roweff, coleff, linlin)

## Analysis of Deviance Table
## 
## Model 1: Freq ~ mental + ses
## Model 2: Freq ~ mental + ses + mental:Cscore
## Model 3: Freq ~ mental + ses + Rscore:ses
## Model 4: Freq ~ mental + ses + Rscore:Cscore
##   Resid. Df Resid. Dev Df Deviance
## 1        15       47.4            
## 2        12        6.3  3     41.1
## 3        10        6.8  2     -0.5
## 4        14        9.9 -4     -3.1

AIC(indep, roweff, coleff, linlin)

##        df   AIC
## indep   9 209.6
## roweff 12 174.5
## coleff 14 179.0
## linlin 10 174.1

mosaic(linlin, residuals_type = "rstandard", labeling_args = long.labels, labeling = labeling_residuals, 
    suppress = 1, gp = shading_Friendly, main = "Mental health data: Linear x Linear")

## Warning: no formula provided, assuming ~ses + mental


## Goodman Row-Column association model fits well (deviance 3.57, df 8)
Mental$mental <- C(Mental$mental, treatment)
Mental$ses <- C(Mental$ses, treatment)
RC1model <- gnm(Freq ~ mental + ses + Mult(mental, ses), family = poisson, data = Mental)

## Initialising
## Running start-up iterations..
## Running main iterations........
## Done

mosaic(RC1model, residuals_type = "rstandard", labeling_args = long.labels, labeling = labeling_residuals, 
    suppress = 1, gp = shading_Friendly, main = "Mental health data: RC(1) model")

## Warning: no formula provided, assuming ~ses + mental