[R] polr probit versus stata oprobit

[Thomas Lumley]

On Wed, 10 Nov 2004, Jean Eid wrote:

> Dear Thomas,
>
> Where you also able to replicate the second example?  (the exaample
> that I turned the housing data into numerical variables) That is the one
> that my estimates differ.
>

I don't have your second example, but I get the same results from
    polr(formula = Sat ~ as.numeric(Infl) + as.numeric(Type) +
      as.numeric(Cont), data = housing, weights = Freq, method = "probit")
and
    oprobit Sat Infl Type Cont [fw=Freq]
for example.

 	-thomas


>
> On Wed, 10 Nov 2004, Thomas Lumley wrote:
>
>> On Wed, 10 Nov 2004, Jean Eid wrote:
>>
>>> Dear All,
>>> I have been struggling to understand why for the housing data in MASS
>>> library R and stata give coef. estimates that are really different. I also
>>> tried to come up with many many examples myself (see below, of course I
>>> did not have the set.seed command included) and all of my
>>> `random' examples seem to give verry similar output. For the housing data,
>>> I have changed the data into numeric vectors instead of factors/ordered
>>> factors. I did so to try and get the same results as in STATA and to have
>>> the housing example as close as possible to the one I constructed.
>>>
>>> I run a debian sid, kernel 2.4, R 2.0.0, and STATA version 8.2, MASS
>>> version  7.2-8.
>>>
>>>
>>> here's the example ( I assume that you have STATA installed and can run in
>>> batch mode, if not the output is also given below)
>>>
>>
>> That example shows the same results with Stata and polr() from MASS.
>>
>> For the housing data, I also get the same coefficients in Stata as with
>> polr():
>>
>> In R:
>> library(MASS)
>> library(foreign)
>> write.dta(housing, file="housing.dta")
>> house.probit<-polr(Sat ~ Infl + Type + Cont, data = housing, weights =
>> Freq, method = "probit")
>> summary(house.probit)
>> -------------------------
>> Re-fitting to get Hessian
>>
>> Call:
>> polr(formula = Sat ~ Infl + Type + Cont, data = housing, weights = Freq,
>>      method = "probit")
>>
>> Coefficients:
>>                     Value Std. Error   t value
>> InflMedium     0.3464233 0.06413706  5.401297
>> InflHigh       0.7829149 0.07642620 10.244063
>> TypeApartment -0.3475372 0.07229093 -4.807480
>> TypeAtrium    -0.2178874 0.09476607 -2.299213
>> TypeTerrace   -0.6641737 0.09180004 -7.235005
>> ContHigh       0.2223862 0.05812267  3.826153
>>
>> Intercepts:
>>              Value   Std. Error t value
>> Low|Medium  -0.2998  0.0762    -3.9371
>> Medium|High  0.4267  0.0764     5.5850
>>
>> Residual Deviance: 3479.689
>> AIC: 3495.689
>> ------------------------
>>
>>
>> In Stata
>> -----------------
>> . use housing.dta
>> . xi: oprobit Sat i.Infl i.Type i.Cont [fw=Freq]
>> i.Infl            _IInfl_1-3          (naturally coded; _IInfl_1 omitted)
>> i.Type            _IType_1-4          (naturally coded; _IType_1 omitted)
>> i.Cont            _ICont_1-2          (naturally coded; _ICont_1 omitted)
>>
>> Iteration 0:   log likelihood = -1824.4388
>> Iteration 1:   log likelihood = -1739.9254
>> Iteration 2:   log likelihood = -1739.8444
>>
>> Ordered probit estimates                          Number of obs   =    1681
>>                                                    LR chi2(6)      =   169.19
>>                                                    Prob > chi2     =   0.0000
>> Log likelihood = -1739.8444                       Pseudo R2       =   0.0464
>>
>> ------------------------------------------------------------------------------
>>           Sat |      Coef.   Std. Err.      z    P>|z|     [95% Conf.
>> Interval]
>> -------------+----------------------------------------------------------------
>>      _IInfl_2 |   .3464228    .064137     5.40   0.000     .2207165     .472129
>>      _IInfl_3 |   .7829146    .076426    10.24   0.000     .6331224    .9327069
>>      _IType_2 |  -.3475367   .0722908    -4.81   0.000    -.4892241   -.2058493
>>      _IType_3 |  -.2178875    .094766    -2.30   0.021    -.4036254   -.0321497
>>      _IType_4 |  -.6641735   .0917999    -7.24   0.000     -.844098    -.484249
>>      _ICont_2 |   .2223858   .0581226     3.83   0.000     .1084676     .336304
>> -------------+----------------------------------------------------------------
>>         _cut1 |  -.2998279   .0761537          (Ancillary parameters)
>>         _cut2 |   .4267208   .0764043
>> ------------------------------------------------------------------------------
>>
>>
>>
>>  	-thomas
>>
>
> ______________________________________________
> 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
>

Thomas Lumley			Assoc. Professor, Biostatistics
tlumley at u.washington.edu	University of Washington, Seattle