Re: [R] polr in MASS

About this list Date view Thread view Subject view Author view Attachment view

From: array chip (arrayprofile@yahoo.com)
Date: Tue 06 May 2003 - 04:53:30 EST


Message-id: <20030505185330.67293.qmail@web41215.mail.yahoo.com>


--- Prof Brian Ripley <ripley@stats.ox.ac.uk> wrote:
> On Mon, 5 May 2003, array chip wrote:
>
> > I misunderstand something here? Or the the way I
> used to reproduce the
> > predicted probabilities is not correct? Thanks
>
> Yes:
>
> 1) How to wrap lines in emails
>
> 2) How to post plain text to mailing lists.
>
> 3) What the model actually is.
>
> --
> Brian D. Ripley,

Sorry that I didn't realize the formatting problem.
Hopefully this time it works:

Hi, I am trying to use the proportional-odds model
using the "polr" function in the MASS library. I tried
with the dataset of "housing" contained in the MASS
book ("Sat" (factor: low, medium, high) is the
dependent variable, "Infl" (low, medium, high), "Type"
(tower, apartment, atrium, terrace) and "Cont" (low,
high) are the predictor variables (factors). And I
have some questions, hope someone could help me out.
The following commands are from the MASS book as well.

>
house.plr<-polr(Sat~Infl+Type+Cont,data=housing,weights=Freq)
> summary(house.plr)

Call: polr(formula = Sat ~ Infl + Type + Cont, data =
housing, weights = Freq)

Coefficients:
                   Value Std. Error t value
InflMedium 0.5663922 0.10465276 5.412109
InflHigh 1.2888137 0.12715609 10.135682
TypeApartment -0.5723552 0.11923800 -4.800107
TypeAtrium -0.3661907 0.15517331 -2.359882
TypeTerrace -1.0910073 0.15148595 -7.202036
ContHigh 0.3602803 0.09553577 3.771156

Intercepts:
            Value Std. Error t value
Low|Medium -0.4961 0.1248 -3.9740
Medium|High 0.6907 0.1255 5.5049

Residual Deviance: 3479.149
AIC: 3495.149

I also tried to predict the probabilities of the 3
categories of "Sat" using the predict function:

> hnames<-lapply(housing[,-5],levels)
>
house.pr1<-predict(house.plr,expand.grid(hnames[-1]),type="probs")>
cbind(expand.grid(hnames[-1]),round(house.pr1,2))

     Infl Type Cont Low Medium High
1 Low Tower Low 0.38 0.29 0.33
2 Medium Tower Low 0.26 0.27 0.47
3 High Tower Low 0.14 0.21 0.65
4 Low Apartment Low 0.52 0.26 0.22
5 Medium Apartment Low 0.38 0.29 0.33
6 High Apartment Low 0.23 0.26 0.51
7 Low Atrium Low 0.47 0.27 0.26
8 Medium Atrium Low 0.33 0.29 0.38
9 High Atrium Low 0.19 0.25 0.56
10 Low Terrace Low 0.64 0.21 0.14
11 Medium Terrace Low 0.51 0.26 0.23
12 High Terrace Low 0.33 0.29 0.38
13 Low Tower High 0.30 0.28 0.42
14 Medium Tower High 0.19 0.25 0.56
15 High Tower High 0.10 0.17 0.72
16 Low Apartment High 0.43 0.28 0.29
17 Medium Apartment High 0.30 0.28 0.42
18 High Apartment High 0.17 0.23 0.60
19 Low Atrium High 0.38 0.29 0.33
20 Medium Atrium High 0.26 0.27 0.47
21 High Atrium High 0.14 0.21 0.64
22 Low Terrace High 0.56 0.25 0.19
23 Medium Terrace High 0.42 0.28 0.30
24 High Terrace High 0.26 0.27 0.47

what I am confused is that when I tried to reproduce
these predicted probabilities using the model
coefficients, I can sometimes get different results:

for example, for low Infl, Type tower, cont low,

logit(P(low))=P(low)/(1-P(low))=-0.4961, solve for
P(low)=exp(-0.4961)/(1+exp(-0.4961))=0.38;

and

logit(P(low,
medium))=P(low,medium)/(1-P(low,medium))=0.6907, solve
for P(low,medium)=exp(0.6907)/(1+exp(0.6907))=0.67,
which is the sum of 0.38 plus 0.29.

The above 2 examples showed that I can reproduce the
predicted probabilities when using the intercept
alone. However, for other combinations of the
predictor variables, I can NOT repreduce the results.

For example, for medium Infl, Type tower, cont low,

logit(P(low))=P(low)/(1-P(low))=-0.4961+0.56639=0.07029,
solve for P(low)=exp(0.07029)/(1+exp(0.07029))=0.52;
but the probability using "predict" function is 0.26.

and

logit(P(low,
medium))=P(low,medium)/(1-P(low,medium))=0.6907+0.56639=1.25709,
solve for
P(low,medium)=exp(1.25709)/(1+exp(1.25709))=0.78,
which certainly is NOT the sum of 0.26 plus 0.27, and
is not the sum of 0.52 and 0.27.

Did I misunderstand something here? Or the the way I
used to reproduce the predicted probabilities is not
correct?

Thanks very much!

______________________________________________
R-help@stat.math.ethz.ch mailing list
https://www.stat.math.ethz.ch/mailman/listinfo/r-help


About this list Date view Thread view Subject view Author view Attachment view

This archive was generated by hypermail 2.1.3 : Tue 01 Jul 2003 - 09:11:45 EST