RE: [R] multinom and contrasts

From: John Fox <jfox_at_mcmaster.ca>
Date: Fri 15 Apr 2005 - 04:09:12 EST


Dear chip,

> -----Original Message-----
> From: r-help-bounces@stat.math.ethz.ch
> [mailto:r-help-bounces@stat.math.ethz.ch] On Behalf Of array chip
> Sent: Thursday, April 14, 2005 1:00 PM
> To: John Fox
> Cc: R-help@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@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@stat.math.ethz.ch
> > > [mailto:r-help-bounces@stat.math.ethz.ch] On
> > Behalf Of array chip
> > > Sent: Wednesday, April 13, 2005 6:26 PM
> > > To: R-help@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@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@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@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 Received on Fri Apr 15 04:17:54 2005

This archive was generated by hypermail 2.1.8 : Fri 03 Mar 2006 - 03:31:11 EST