```postscript("l09.ps");options(width=68)
#****************************************************/
# Housing data from Cox and Snell (1980), Example W,*/
# on satisfaction with housing in Copenhagen.  Data */
# was extracted from R package MASS.  Data may also */
# be found at                                       */
# www.stat.ucla.edu/data/cox-and-snell/exampleW.data*/
#****************************************************/
housing<-as.data.frame(scan("housing.dat",
what=list(nn=0,sat="",infl="",type="",contact="",freq=0)))
housing\$nsat<-(housing\$sat=="Medium")+2*(housing\$sat=="High")
housing\$ninf<-(housing\$infl=="Medium")+2*(housing\$infl=="High")
housing\$ncont<-(housing\$contact=="Medium")+
2*(housing\$contact=="High")
fullhouse<-housing[rep(1:72,housing\$freq),]
housing\$rsat<-as.character(housing\$sat)
housing\$rsat[housing\$sat=="Medium"]<-"Intermediate"
housing\$rsat<-as.factor(housing\$rsat)
#*************************************************************/
# Mark C: Polytomous Regression.  Baseline logit model first.*/
#*************************************************************/
cat('\n Baseline logit model  \n')
#To make the numbers come close to those we get from SAS,
#rename the categories to change the ordering. Resulting
#effect coefficients will be double those of SAS, since
#proc logistic by default uses the -1,1 parameterization.
housing\$rsat<-as.character(housing\$sat)
housing\$rsat[housing\$sat=="High"]<-"Top"
housing\$rsat<-as.factor(housing\$rsat)
library(nnet)
multinom(rsat~contact,data=housing,weight=freq)
#*************************************************************/
# Mark D: Cumulative logit model.  First reuse a known tech- */
# nique.  Simple model exhibits different parameterizations. */
#*************************************************************/
cat('\n Binary response for comparison purpose  \n')
housing\$oksat<-housing\$nsat>1
glm(oksat~factor(contact),family=binomial,data=housing,
weight=freq)
#*************************************************************/
# The following data set will not give us what we want, be-  */
# cause it will order sat incorrectly.  Recode data to put   */
# medium between high and low, and add descending to make    */
# low the baseline.  First fit the simplest model.  This     */
# should give logits of raw proportions.  Check this.        */
#*************************************************************/
cat('\n Ordinal Logistic Regression with wrong order  \n')
cat('\n Ordinal Logistic Reg. with right order, empty model  \n')
library(MASS)
# Venables and Ripley, Modern Applied Statistics with S
# polr does regression with ordered response.   Allows var-
polr(factor(nsat)~1,weights=freq,data=housing)
cat('\n Ordinal Logist. Reg. with added covariates  \n')
polr(factor(nsat)~infl+type+contact,weights=freq,data=housing)
#*************************************************************/
# Mark E: Cf. the complimentary log log.  Odds ratios are no */
# longer given, since they don't come naturally out of the   */
# model.  Intercepts have a completely different meaning and */
# won't be compared.  Linear parameters are multiplied by a  */
# constant, approx., and z statistics are similar.           */
#*************************************************************/
cat('\n Complimentary log log  \n')
polr(factor(nsat)~infl+type+contact,weights=freq,data=housing,
method="cloglog")
```