R help - Jan 2010 - splitting a factor in an analysis of deviance table (negative binomial model)

Dears useRs,

I have 2 factors, (for the sake of explanation - A and B), with 4 levels each.
I've already fitted a negative binomial generalized linear model to my data,
and now I need to split the factors in two distinct analysis of deviance table:
 - A within B1, A within B2, A within B3 and A within B4
 - B within A1, B within A2, B within A3 and B within A4

Here is a code that illustrates my problem:

# inputing my data

require(MASS)

my.data <-
data.frame("A"=rep(c(rep("Alevel1",4),rep("Alevel2",4),rep("Alevel3",4),rep("Alevel4",4)),4),
                     
"B"=c(rep("Blevel1",16),rep("Blevel2",16),rep("Blevel3",16),rep("Blevel4",16)),
                      "value"=rnegbin(64, 10, 10))

# fitting the model with interaction (a + b + a:b)

model <- glm.nb(value ~ A*B, data=my.data)
anova(model, test="F")

     Df Deviance Resid. Df Resid. Dev      F Pr(>F)
NULL                    63     80.639
A     3    1.374        60     79.265 0.4581 0.7115
B     3    2.285        57     76.980 0.7616 0.5155
A:B   9    9.700        48     67.280 1.0778 0.3753

#until here it's ok, now I need to redistribute the 12 degrees of freedom
(B+A:B -> 3+9) into 4 splitted factors (A within B1, A within B2...)

model2 <- glm.nb(value ~
A/((B=="Blevel1")+(B=="Blevel2")+(B=="Blevel3")+(B=="Blevel4")),
data=my.data)
anova(model2, test="F")
                 Df Deviance Resid. Df Resid. Dev      F Pr(>F)
NULL                                63     80.639
A                 3    1.374        60     79.265 0.4581 0.7115
A:B == "Blevel1"  4    3.871        56     75.394 0.9676 0.4238
A:B == "Blevel2"  4    5.236        52     70.158 1.3090 0.2639
A:B == "Blevel3"  4    2.878        48     67.280 0.7195 0.5784
A:B == "Blevel4"  0    0.000        48     67.280

However, the 12 degrees of freedom are distributed as 4,4,4,0 and not 3,3,3,3,
as I expected.
Is there a way of obtaining the split I need?

Thanks in advance, all the best!

Rafael.


     
____________________________________________________________________________________
[[elided Yahoo spam]]

	[[alternative HTML version deleted]]

Sorry by the delay. You could do:
> my.data <- expand.grid(A=factor(1:4), B=factor(1:4), rep=1:4)
> my.data$y <- rbinom(my.data$A, 10, 0.5)
> 
> model <- glm(cbind(y, 10-y)~A*B, family=binomial, data=my.data)
> anova(model, test="Chisq" )
Analysis of Deviance Table

Model: binomial, link: logit

Response: cbind(y, 10 - y)

Terms added sequentially (first to last)


     Df Deviance Resid. Df Resid. Dev P(>|Chi|)
NULL                    63     53.528          
A     3   1.6792        60     51.849    0.6416
B     3   2.9380        57     48.911    0.4013
A:B   9  13.8380        48     35.073    0.1282
> 
> X <- model.matrix(~A/B, my.data)
> 
> A <- X[,2:4]
> B.A1 <- X[,grep("A1:", colnames(X))]
> B.A2 <- X[,grep("A2:", colnames(X))]
> B.A3 <- X[,grep("A3:", colnames(X))]
> B.A4 <- X[,grep("A4:", colnames(X))]
> 
> model2 <- glm(cbind(y, 10-y)~A+B.A1+B.A2+B.A3+B.A4, family=binomial,
> data=my.data)
> anova(model2, test="Chisq")
Analysis of Deviance Table

Model: binomial, link: logit

Response: cbind(y, 10 - y)

Terms added sequentially (first to last)


     Df Deviance Resid. Df Resid. Dev P(>|Chi|)  
NULL                    63     53.528            
A     3   1.6792        60     51.849   0.64156  
B.A1  3   1.2924        57     50.556   0.73093  
B.A2  3   4.3126        54     46.244   0.22962  
B.A3  3   1.8865        51     44.357   0.59629  
B.A4  3   9.2844        48     35.073   0.02574 *
---
Signif. codes:  0 ?***? 0.001 ?**? 0.01 ?*? 0.05 ?.? 0.1 ? ? 1
> 
Walmes.

-----
..ooo0
...................................................................................................
..(....)... 0ooo...                              Walmes Zeviani
...\..(.....(.....)...     Master in Statistics and Agricultural
Experimentation
....\_)..... )../....       walmeszeviani at hotmail.com, Lavras - MG, Brasil
............
(_/............................................................................................
-- 
View this message in context:
http://r.789695.n4.nabble.com/splitting-a-factor-in-an-analysis-of-deviance-table-negative-binomial-model-tp1017744p2235345.html
Sent from the R help mailing list archive at Nabble.com.