In Venables \& Ripley 3rd edition (p. 231) the proportional odds model
is described as:
logit(p<=k) = zeta_k + eta
but polr apparently thinks there is a minus in front of eta,
as is apprent below.
Is this a bug og a feature I have overlooked?
Here is the naked code for reproduction, below the results.
------------------------------------------------------------------------
---
version
library( MASS )
data( housing )
hnames <- lapply( housing[,-5], levels )
house.plr <- polr( Sat ~ Infl + Type + Cont, data=housing, weights=Freq
)
summary( house.plr )
newdat <- expand.grid( hnames[-1] )[1:5,]
cbind( newdat, predict( house.plr, newdat, type="probs" ) )
# Baseline probs:
diff( c(0,tigol( c(-0.4961,0.6907) ), 1) )
# First level of Infl:
diff( c(0,tigol( c(-0.4961,0.6907) + 0.5663922 ), 1) )
# But the change of sign for eta is needed to reproduce the fitted
values:
# Line 2:
diff( c(0,tigol( c(-0.4961,0.6907) - 0.5663922 ), 1) )
# Line 5:
diff( c(0,tigol( c(-0.4961,0.6907) - 0.5663922 + 0.5723552 ), 1) )
------------------------------------------------------------------------
---
Here is the resulting output:
------------------------------------------------------------------------
---
version
_
platform i386-pc-mingw32
arch i386
os mingw32
system i386, mingw32
status
major 1
minor 8.1
year 2003
month 11
day 21
language R
library( MASS )
data( housing )
hnames <- lapply( housing[,-5], levels )
house.plr <- polr( Sat ~ Infl + Type + Cont, data=housing,
weights=Freq )
summary( house.plr )
Re-fitting to get Hessian
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.3661908 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