[R] multinom and contrasts

[John Fox]

Dear chip,

> -----Original Message-----
> From: r-help-bounces at stat.math.ethz.ch 
> [mailto:r-help-bounces at stat.math.ethz.ch] On Behalf Of array chip
> Sent: Thursday, April 14, 2005 1:00 PM
> To: John Fox
> Cc: R-help at stat.math.ethz.ch
> Subject: RE: [R] multinom and contrasts
> 
> Dear John,
> 
> Thanks for the answer! In my own dataset, The
> multinom() did not converge even after I had tried to 
> increase the maximum number of iteration (from default 100 to 
> 1000). In this situation, there is some bigger diffrenece in 
> fitted probabilities under different contrasts (e.g. 
> 0.9687817 vs. 0.9920816). My question is whether the analysis 
> (fitted probabilities) is still valid if it does not 
> converge? and what else can I try about it?
> 

If multinom() doesn't converge to a stable solution after 1000 iterations,
it's probably safe to say that the problem is ill-conditioned in some
respect. Have you looked at the covariance matrix of the estimates?

Regards,
 John

> Thank you!
> 
> 
> 
> 
> --- John Fox <jfox at mcmaster.ca> wrote:
> > Dear chip,
> > 
> > The difference is small and is due to computational error.
> > 
> > Your example:
> > 
> > > max(abs(zz[1:10,] - yy[1:10,]))
> > [1] 2.207080e-05
> > 
> > Tightening the convergence tolerance in multinom() eliminates the
> > difference:
> > 
> > >
> > options(contrasts=c('contr.treatment','contr.poly'))
> > >
> >
> xx<-multinom(Type~Infl+Cont,data=housing[-c(1,10,11,22,25,30),],
> > reltol=1.0e-12)
> > # weights:  20 (12 variable)
> > initial  value 91.495428
> > iter  10 value 91.124526
> > final  value 91.124523
> > converged
> > > yy<-predict(xx,type='probs')
> > > options(contrasts=c('contr.helmert','contr.poly'))
> > >
> >
> xx<-multinom(Type~Infl+Cont,data=housing[-c(1,10,11,22,25,30),],
> > reltol=1.0e-12)
> > # weights:  20 (12 variable)
> > initial  value 91.495428
> > iter  10 value 91.125287
> > iter  20 value 91.124523
> > iter  20 value 91.124523
> > iter  20 value 91.124523
> > final  value 91.124523
> > converged
> > > zz<-predict(xx,type='probs')
> > > max(abs(zz[1:10,] - yy[1:10,]))
> > [1] 1.530021e-08
> > 
> > I hope this helps,
> >  John
> > 
> > --------------------------------
> > John Fox
> > Department of Sociology
> > McMaster University
> > Hamilton, Ontario
> > Canada L8S 4M4
> > 905-525-9140x23604
> > http://socserv.mcmaster.ca/jfox
> > --------------------------------
> > 
> > > -----Original Message-----
> > > From: r-help-bounces at stat.math.ethz.ch 
> > > [mailto:r-help-bounces at stat.math.ethz.ch] On
> > Behalf Of array chip
> > > Sent: Wednesday, April 13, 2005 6:26 PM
> > > To: R-help at stat.math.ethz.ch
> > > Subject: [R] multinom and contrasts
> > > 
> > > Hi,
> > > 
> > > I found that using different contrasts (e.g.
> > > contr.helmert vs. contr.treatment) will generate
> > different
> > > fitted probabilities from multinomial logistic
> > regression
> > > using multinom(); while the fitted probabilities
> > from binary
> > > logistic regression seem to be the same. Why is
> > that? and for
> > > multinomial logisitc regression, what contrast
> > should be
> > > used? I guess it's helmert?
> > > 
> > > here is an example script:
> > > 
> > > library(MASS)
> > > library(nnet)
> > > 
> > >       #### multinomial logistic
> > >
> > options(contrasts=c('contr.treatment','contr.poly'))
> > >
> >
> xx<-multinom(Type~Infl+Cont,data=housing[-c(1,10,11,22,25,30),])
> > > yy<-predict(xx,type='probs')
> > > yy[1:10,]
> > > 
> > > options(contrasts=c('contr.helmert','contr.poly'))
> > >
> >
> xx<-multinom(Type~Infl+Cont,data=housing[-c(1,10,11,22,25,30),])
> > > zz<-predict(xx,type='probs')
> > > zz[1:10,]
> > > 
> > > 
> > >       ##### binary logistic
> > >
> > options(contrasts=c('contr.treatment','contr.poly'))
> > >
> >
> obj.glm<-glm(Cont~Infl+Type,family='binomial',data=housing[-c(
> > 1,10,11,22,25,30),])
> > > yy<-predict(xx,type='response')
> > > 
> > > options(contrasts=c('contr.helmert','contr.poly'))
> > >
> >
> obj.glm<-glm(Cont~Infl+Type,family='binomial',data=housing[-c(
> > 1,10,11,22,25,30),])
> > > zz<-predict(xx,type='response')
> > > 
> > > Thanks
> > > 
> > > ______________________________________________
> > > R-help at stat.math.ethz.ch mailing list 
> > > https://stat.ethz.ch/mailman/listinfo/r-help
> > > PLEASE do read the posting guide! 
> > > http://www.R-project.org/posting-guide.html
> > 
> >
> 
> ______________________________________________
> R-help at stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide! 
> http://www.R-project.org/posting-guide.html