[R] Multinomial logistic regression under R and Stata

[Thomas W Blackwell]

Wing  -

Apologies if I'm overdoing the answer to a fundamentally simple question.

There's another experiment you might do:  Try doing the orthogonalization
between quadratic and linear terms yourself, explicitly.  (See below.)
Find out whether and how much the reported standard errors change in each
package, when you do the orthogonalization explicitly, ahead of time,
versus when you let the multinom() fitting routine use its own default
procedure.  The package whose standard error estimate doesn't change is
doing it right.

I think it's not possible to do the orthogonalization exactly, because
one doesn't have access ahead of time to the weights that will be used
in the multinomial fit.  But maybe approximate is good enough to show
whether there is a difference in the standard errors.

Suppose your "social status" covariate is a column in a data frame
 d1  named "soc.st".  The non-orthogonalized quadratic term would be
"soc.st^2".  Define a new data frame, d2, by

 d2 <- cbind(d1, soc.quad = residuals(lm(soc.st^2 ~ soc.st, data=d1)))

The column named "soc.quad" in d2 is the quadratic term, crudely
orthogonalized to the constant and linear terms in soc.st.  Use it
in your multinomial regression in place of "soc.st^2", (if that's
what you were doing before) and see whether both packages give the
same estimates and standard errors that they did before.  I'll be
most surprised if the estimates are different, but I expect that
one package or the other will change its reported standard errors.

-  tom blackwell  -  u michigan medical school  -  ann arbor  -



On Thu, 27 Mar 2003, Tak Wing Chan wrote:

> Dear Colleagues
>
> I have been fitting some multinomial logistic regression models using R
> (version 1.6.1 on a linux box) and Stata 7. Although the vast majority
> of the parameter estimates and standard errors I get from R are the same
> as those from Stata (given rounding errors and so on), there are a few
> estimates for the same model which are quite different. I would be most
> grateful if colleagues could advise me as to what might be causing this,
> and should I worry ...
>
> Anyway, with R, I have been using the function multinom under the
> package nnet. Below are two examples where the estimates for standard
> error differ substantially between R and Stata:
>
>           beta              s.e.
> R:        5.939880 2.920165
> Stata:  5.939747 5.455495
>
> R:      11.228705 2.191625
> Stata: 11.22761  4.630293
>
> The parameters concerned are the quadratic term of a quantitative
> variable (measuring social status). I notice that the s.e. for this
> quadratic term are large anyway compared to other s.e. in the model.
>
> There are other differences between R and Stata, and these concerned the
> intercept terms. Here is an example:
>
>            beta            s.e.
> R:        0.2870793 0.4512347
> Stata: -0.2109653 0.5053566
>
> Since both estimates are not significantly different from zero, I trust
> I can ignore the difference between the estimates. Or could I?
>
> Many thanks in advance for any help. Please let me know if I should
> provide further info.
>
> With best wishes.
>
> Wing
> --
> Department of Sociology, University of Oxford,
> Littlegate House, St Ebbes, Oxford OX1 1PT, UK
> tel: +44 (1865) 286176, fax: +44 (1865) 286171
> http://users.ox.ac.uk/~sfos0006