> data <- read.table("sex.dat", header=T)
> attach(data, 1)
> options(contrasts=c("contr.treatment"))
> uv <- sex*birth
> sex <- factor(sex); birth <- factor(birth)
> > indep.fit <- glm(count ~ sex + birth, family=poisson)
> pearson <- summary.lm(indep.fit)$residuals
> hat <- lm.influence(indep.fit)$hat
> adj.res <- pearson/sqrt(1 - hat) # standardized Pearson residuals
> cbind(sex, birth, count, fitted(indep.fit), adj.res)
sex birth count adj.res
1 1 1 81 42.41145 7.6028389
2 1 2 68 51.21382 3.0766897
3 1 3 60 86.42333 -4.1166720
4 1 4 38 66.95140 -4.8396108
5 2 1 24 15.96868 2.3282276
6 2 2 26 19.28294 1.8114763
7 2 3 29 32.53996 -0.8114805
8 2 4 14 25.20842 -2.7568090
9 3 1 18 30.04860 -2.6816389
10 3 2 41 36.28510 0.9762279
11 3 3 74 61.23110 2.2472882
12 3 4 42 47.43521 -1.0263685
13 4 1 36 70.57127 -6.0630727
14 4 2 57 85.21814 -4.6038421
15 4 3 161 143.80562 2.3845449
16 4 4 157 111.40497 6.7845091
> LL.fit <- update(indep.fit, . ~ . + uv)
> summary(LL.fit)
Coefficients:
Value Std. Error t value
(Intercept) 4.1068414 0.08950939 45.881681
sex2 -1.6459643 0.13472875 -12.216875
sex3 -1.7700229 0.16463639 -10.751104
sex4 -1.7536852 0.23430240 -7.484708
birth2 -0.4641050 0.11952096 -3.883043
birth3 -0.7245224 0.16199835 -4.472406
birth4 -1.8796640 0.24908191 -7.546369
uv 0.2858355 0.02823566 10.123209
Residual Deviance: 11.53369 on 8 degrees of freedom
Number of Fisher Scoring Iterations: 3
> u <- c(1,1,1,1,2,2,2,2,4,4,4,4,5,5,5,5)
> v <- c(1,2,4,5,1,2,4,5,1,2,4,5,1,2,4,5)
> uv2 <- u*v
> LL2.fit <- update(indep.fit, . ~ . + uv2)
> summary(LL2.fit)
Coefficients:
Value Std. Error t value
uv2 0.1460382 0.01407636 10.374715
Residual Deviance: 8.845195 on 8 degrees of freedom
> row.fit <- update(indep.fit, .~. +u:birth)
> summary(row.fit)
Value Std. Error t value
(Intercept) 4.9872158 0.14624227 34.102424
sex2 -0.6577204 0.13124054 -5.011564
sex3 0.4666360 0.16265012 2.868956
sex4 1.5019483 0.17951515 8.366694
birth2 -0.3193886 0.19821146 -1.611353
birth3 -0.7268794 0.20015708 -3.631545
birth4 -1.4903163 0.23744302 -6.276522
birth1u -0.5953263 0.06555007 -9.082009
birth2u -0.4054336 0.06067887 -6.681628
birth3u -0.1297489 0.05633659 -2.303102
birth4u NA NA NA
Residual Deviance: 8.262966 on 6 degrees of freedom
> column.fit <- update(indep.fit, . ~ . + sex:v)
> summary(column.fit)
Coefficients: (1 not defined because of singularities)
Value Std. Error t value
(Intercept) 4.9151850 0.13907409 35.342204
sex2 -1.2186815 0.25582686 -4.763696
sex3 -1.5060354 0.23760340 -6.338442
sex4 -1.3955940 0.20728940 -6.732587
birth2 0.5458980 0.11723030 4.656629
birth3 1.5926216 0.14786815 10.770553
birth4 1.5101843 0.16419368 9.197579
sex1v -0.5845440 0.05929738 -9.857840
sex2v -0.4955406 0.07989777 -6.202182
sex3v -0.2031464 0.06537858 -3.107231
sex4v NA NA NA
Residual Deviance: 7.586124 on 6 degrees of freedom