[Rd] weird behavior of drop1() for polr models (MASS)

[Peter Dalgaard]

Jeroen Ooms wrote:
> I would like to do a SS type III analysis on a proportional odds logistic
> regression model. I use drop1(), but dropterm() shows the same behaviour. It
> works as expected for regular main effects models, however when the model
> includes an interaction effect it seems to have problems with matching the
> parameters to the predictor terms. An example:
>
> library("MASS");
> options(contrasts = c("contr.treatment", "contr.poly"));
>
> house.plr1 <- polr(Sat ~ Infl + Type + Cont, weights = Freq, data =
> housing);
> drop1(house.plr1,attributes(house.plr1$terms)$term.labels,test="Chisq");
>
> house.plr2 <- polr(Sat ~ Infl * Type + Cont, weights = Freq, data =
> housing); 
> drop1(house.plr2,attributes(house.plr2$terms)$term.labels,test="Chisq");
>
> Notice that model 2 has a * instead of a + between predictors Infl and Type.
> In model 1, estimated parameters are nicely attributed to the right term,
> however in house.plr2, only 2 of the 4 terms are evaluated.
>
> I am using R version 2.7.2 (2008-08-25) i386-pc-mingw32, and MASS_7.2-44.
>   

(It would have been better form to actually print the result so that
people can see what you're on about:
> drop1(house.plr2,attributes(house.plr2$terms)$term.labels,test="Chisq");
Single term deletions

Model:
Sat ~ Infl * Type + Cont
          Df    AIC        LRT   Pr(Chi)
<none>       3484.6
Infl       0 3484.6  2.086e-08
Type       0 3484.6 -1.301e-10
Cont       1 3497.8       15.1 0.0001009 ***
Infl:Type  6 3495.1       22.5 0.0009786 ***

I take it that the perceived problem is the two 0-Df lines.)



This is slightly peculiar, yes, but notice that in lm() or glm() you'd
normally just not get entries for main effects in the presence of their
interaction. This is deliberate -- type III tests can be seriously
misleading and/or uninterpretable in that case, so you really are better
off without them.

The direct cause of the behaviour above is that in R,

Infl+Type+Infl:Type

and

Type+Infl:Type

are just two parametrizations of the same model. Try


summary(house.plr2)
summary(update(house.plr2,~.-Type))

-- 
   O__  ---- Peter Dalgaard             Øster Farimagsgade 5, Entr.B
  c/ /'_ --- Dept. of Biostatistics     PO Box 2099, 1014 Cph. K
 (*) \(*) -- University of Copenhagen   Denmark      Ph:  (+45) 35327918
~~~~~~~~~~ - (p.dalgaard at biostat.ku.dk)              FAX: (+45) 35327907