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

From: Peter Dalgaard <P.Dalgaard_at_biostat.ku.dk>
Date: Tue, 30 Sep 2008 16:12:01 +0200

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

______________________________________________
R-devel_at_r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-devel
Received on Tue 30 Sep 2008 - 14:16:13 GMT

Archive maintained by Robert King, hosted by the discipline of statistics at the University of Newcastle, Australia.
Archive generated by hypermail 2.2.0, at Tue 30 Sep 2008 - 19:30:14 GMT.

Mailing list information is available at https://stat.ethz.ch/mailman/listinfo/r-devel. Please read the posting guide before posting to the list.

list of date sections of archive