[Fred Dzanku]

Dear Stephen, Joseph and All,
I appreciate your input very much, thank you. I'll explore the options suggested. Best wishes.

[Richard Williams]

Joseph's gsem is pretty neat. I am not sure why the square root of 2 bit is in there but it seems to work. But are there any random draws in there? If so, how? I guess I should break down and read Stephen's article but if gsem somehow avoids the random component that would seem to be nice.

[Fred Dzanku]

I would like to point out that I have used David Roodman's cmp with both random and halton draws and got very similar results when compared with mvprobit and gsem. I also compared estimates from a bivariate probit between the statistical packages and programs just to see whether or not the difference really lies with the simulation. The results (below) shows that, as suggested by others, the differences in simulation approaches appear to account for the differing results. (Although the standard errors of the estimates are still different in one case such that testing at the 10% level one will reject the null that x1 = 0 in Stata but fail to reject it in Limdep).

I will contact the Limdep team for further clarification on what might be happening on their side.

************************** COMPARING BIVARIATE PROBIT ESTIMATES IN STATA & LIMDEP ********************

______________ Bivariate Probit in Stata _____________________

biprobit (y1 = x1 x2 x3) (y2 = x1 x2 x3 x4 x5),cl(vid)

______________ Results ___________________________________

Seemingly unrelated bivariate probit            Number of obs     =        650
                                                Wald chi2(8)      =      22.26
Log pseudolikelihood = -851.23085               Prob > chi2       =     0.0045

                                   (Std. Err. adjusted for 66 clusters in vid)
------------------------------------------------------------------------------
             |               Robust
             |      Coef.   Std. Err.      z    P>|z|     [95% Conf. Interval]
-------------+----------------------------------------------------------------
y1           |
          x1 |  -.0067399   .0039496    -1.71   0.088    -.0144811    .0010012
          x2 |   .4743763    .135101     3.51   0.000     .2095832    .7391695
          x3 |   .0155758   .0076112     2.05   0.041     .0006581    .0304936
       _cons |   .0349337   .2000827     0.17   0.861    -.3572213    .4270886
-------------+----------------------------------------------------------------
y2           |
          x1 |   .0058175   .0045848     1.27   0.204    -.0031685    .0148035
          x2 |   .1756506   .0819487     2.14   0.032     .0150341    .3362671
          x3 |  -.0052762   .0065862    -0.80   0.423     -.018185    .0076326
          x4 |   .0878311   .1200443     0.73   0.464    -.1474514    .3231136
          x5 |   .1879889   .1404749     1.34   0.181    -.0873368    .4633146
       _cons |  -.6281078    .224958    -2.79   0.005    -1.069017   -.1871983
-------------+----------------------------------------------------------------
     /athrho |   .0170186   .0778276     0.22   0.827    -.1355207    .1695579
-------------+----------------------------------------------------------------
         rho |    .017017   .0778051                     -.1346971    .1679515
------------------------------------------------------------------------------
Wald test of rho=0: chi2(1) = .047817                     Prob > chi2 = 0.8269

______________ Bivariate Probit in Limdep _____________________

Sample    ;all $
Skip$

        
Bivariate    ;Lhs=y1,y2
        ;Rh1=x1,x2,x3,one
        ;Rh2=x1,x2,x3,x4,x5,one
        ;cluster=vid $

______________ Results ___________________________________

Normal exit from iterations. Exit status=0.

+---------------------------------------------+
| FIML Estimates of Bivariate Probit Model    |
| Maximum Likelihood Estimates                |
| Model estimated: Mar 17, 2016 at 02:11:14PM.|
| Dependent variable                 Y1Y2     |
| Weighting variable                 None     |
| Number of observations              650     |
| Iterations completed                 16     |
| Log likelihood function       -851.2309     |
| Number of parameters                 11     |
| Info. Criterion: AIC =          2.65302     |
|   Finite Sample: AIC =          2.65365     |
| Info. Criterion: BIC =          2.72878     |
| Info. Criterion:HQIC =          2.68241     |
+---------------------------------------------+
+---------------------------------------------------------------------+
| Covariance matrix for the model is adjusted for data clustering.    |
| Sample of    650 observations contained     66 clusters defined by  |
| variable VID      which identifies by a value a cluster ID.         |
| Sample of    650 observations contained      1 strata defined by    |
|    650 observations (fixed number) in each stratum.                 |
+---------------------------------------------------------------------+
+--------+--------------+----------------+--------+--------+----------+
|Variable| Coefficient  | Standard Error |b/St.Er.|P[|Z|>z]| Mean of X|
+--------+--------------+----------------+--------+--------+----------+
---------+Index    equation for Y1
 X1      |    -.00673993       .00415199    -1.623   .1045   46.6107692
 X2      |     .47437629       .09061974     5.235   .0000    .34307692
 X3      |     .01557582       .00727657     2.141   .0323   11.2224103
 Constant|     .03493363       .20420921      .171   .8642
---------+Index    equation for Y2
 X1      |     .00581750       .00474579     1.226   .2203   46.6107692
 X2      |     .17565058       .07533074     2.332   .0197    .34307692
 X3      |    -.00527620       .00680249     -.776   .4380   11.2224103
 X4      |     .08783110       .11983375      .733   .4636    .28923077
 X5      |     .18798893       .14145261     1.329   .1839    .11538462
 Constant|    -.62810784       .23312962    -2.694   .0071
---------+Disturbance correlation
 RHO(1,2)|     .01701717       .07854286      .217   .8285

_____________________ End _______________________________

[Richard Williams]

Does anybody have a tie-breaker program? e.g. can you estimate an mprobit model in SPSS or SAS? My bias is always with Stata but it would be nice to see what another program said.

[Richard Williams]

Fred inspired me to dig out my copy of Limdep. I ran his model but WITHOUT the cluster option. The standard errors are much closer to what Stata produced. This makes me wonder if Limdep is handling clustering correctly; or perhaps something needs to be done differently with clustering in Limdep. Here is what I got:

|-> Sample    ;all $
|-> Skip $
|-> Mprobit    ;LHS=Y1, Y2, Y3
    ;eq1=X1, X2, X3,one
    ;eq2=X1, X2, X3, X4, X5,one
    ;eq3=X2, X3, X4, X6,one
    ;pts=200;$
Normal exit:  25 iterations. Status=0, F=    1242.816

-----------------------------------------------------------------------------
Multivariate Probit Model:  3 equations.
Dependent variable             MVProbit
Log likelihood function     -1242.81585
Estimation based on N =    650, K =  18
Inf.Cr.AIC  =   2521.6 AIC/N =    3.879
Model estimated: Mar 17, 2016, 11:33:39
Replications for simulated probs. = 200
--------+--------------------------------------------------------------------
        |                  Standard            Prob.      95% Confidence
MVProbit|  Coefficient       Error       z    |z|>Z*         Interval
--------+--------------------------------------------------------------------
        |Index function for Y1
      X1|    -.00703*        .00369    -1.91  .0567     -.01426    .00020
      X2|     .47399***      .06512     7.28  .0000      .34635    .60163
      X3|     .01570**       .00656     2.39  .0167      .00284    .02855
Constant|     .04693         .17507      .27  .7887     -.29620    .39005
        |Index function for Y2
      X1|     .00770**       .00370     2.08  .0371      .00046    .01495
      X2|     .17352**       .07248     2.39  .0167      .03146    .31558
      X3|    -.00602         .00679     -.89  .3757     -.01933    .00730
      X4|     .10489         .11601      .90  .3659     -.12248    .33226
      X5|     .09890         .15990      .62  .5362     -.21450    .41230
Constant|    -.69995***      .19718    -3.55  .0004    -1.08642   -.31348
        |Index function for Y3
      X2|     .17529**       .07357     2.38  .0172      .03111    .31948
      X3|     .00311         .00687      .45  .6509     -.01036    .01658
      X4|     .16664         .11144     1.50  .1348     -.05179    .38506
      X6|     .45556**       .19530     2.33  .0197      .07277    .83835
Constant|    -.55654***      .10249    -5.43  .0000     -.75741   -.35567
        |Correlation coefficients
R(01,02)|     .01254         .06507      .19  .8472     -.11499    .14007
R(01,03)|    -.05260         .06603     -.80  .4257     -.18201    .07681
R(02,03)|     .40508***      .05845     6.93  .0000      .29052    .51964
--------+--------------------------------------------------------------------
Note: ***, **, * ==>  Significance at 1%, 5%, 10% level.
-----------------------------------------------------------------------------

[Richard Williams]

For good measure, here is what mvprobit gives if you drop the cluster option. This looks very very close to Limdep without clustering. Which makes me think even more that the differences in the programs are due to differences in the way they handle clustering. The Stata results seem much more plausible to me.

. mvprobit (y1 x1 x2 x3) (y2 x1 x2 x3 x4 x5) (y3 x2 x3 x4 x6),  dr(200) nolog seed (1003)

Multivariate probit (SML, # draws = 200)          Number of obs   =        650
                                                  Wald chi2(12)   =      71.74
Log likelihood = -1242.8793                       Prob > chi2     =     0.0000

------------------------------------------------------------------------------
             |      Coef.   Std. Err.      z    P>|z|     [95% Conf. Interval]
-------------+----------------------------------------------------------------
y1           |
          x1 |   -.006999   .0036266    -1.93   0.054    -.0141071     .000109
          x2 |   .4740413    .075187     6.30   0.000     .3266776    .6214051
          x3 |   .0156802   .0066407     2.36   0.018     .0026647    .0286957
       _cons |   .0458113   .1743776     0.26   0.793    -.2959625    .3875852
-------------+----------------------------------------------------------------
y2           |
          x1 |   .0076917   .0036487     2.11   0.035     .0005403    .0148431
          x2 |   .1741295   .0687453     2.53   0.011     .0393911    .3088679
          x3 |  -.0060159   .0066718    -0.90   0.367    -.0190924    .0070606
          x4 |   .1035748   .1155102     0.90   0.370     -.122821    .3299706
          x5 |   .0982916   .1540134     0.64   0.523    -.2035692    .4001523
       _cons |  -.6988311   .1923439    -3.63   0.000    -1.075818    -.321844
-------------+----------------------------------------------------------------
y3           |
          x2 |    .175141   .0688806     2.54   0.011     .0401375    .3101444
          x3 |   .0032237   .0066113     0.49   0.626    -.0097341    .0161816
          x4 |   .1668598   .1115449     1.50   0.135    -.0517642    .3854838
          x6 |   .4567178   .1871014     2.44   0.015     .0900058    .8234298
       _cons |  -.5554719     .09958    -5.58   0.000    -.7506451   -.3602987
-------------+----------------------------------------------------------------
    /atrho21 |   .0115043    .064559     0.18   0.859    -.1150291    .1380376
-------------+----------------------------------------------------------------
    /atrho31 |  -.0480755   .0656746    -0.73   0.464    -.1767954    .0806444
-------------+----------------------------------------------------------------
    /atrho32 |   .4273434   .0693806     6.16   0.000       .29136    .5633268
-------------+----------------------------------------------------------------
       rho21 |   .0115038   .0645505     0.18   0.859    -.1145244    .1371675
-------------+----------------------------------------------------------------
       rho31 |  -.0480385   .0655231    -0.73   0.463    -.1749762    .0804701
-------------+----------------------------------------------------------------
       rho32 |   .4030988    .058107     6.94   0.000     .2833861    .5104416
------------------------------------------------------------------------------
Likelihood ratio test of  rho21 = rho31 = rho32 = 0:  
             chi2(3) =  40.9908   Prob > chi2 = 0.0000

.

[Fred Dzanku]

This is informative. Will let you know what William Greene says about this.

[Fred Dzanku]

Below is what I got from William Greene, Department of Economics, Stern School of Business, New York University, who wrote the Limdep (or NLOGIT) program. This response is in connection with the trivariate probit estimates with cluster robust SEs I posted earlier.[INDENT=2]"Fred. The estimates are not that far apart. And, note that the log[/INDENT][INDENT=2]likelihood values are essentially the same. I can tell you that nlogit is using the[/INDENT][INDENT=2]GHK simulator to do the estimation. I do not know what Stata is doing. You can be sure[/INDENT][INDENT=2]that the estimated coefficients are very sensitive to the values of the correlation[/INDENT][INDENT=2]coefficients. But, in order to compare the results, you should probably compare the[/INDENT][INDENT=2]estimated partial effects. I don't know if Stata knows how to do that. The nlogit[/INDENT][INDENT=2]command is described in the manual."[/INDENT][INDENT=2]/B. Greene[/INDENT]

[Fred Dzanku]

Apologies the previous post was not clear.

Below is what I got from William Greene, Department of Economics, Stern School of Business, New York University, who wrote the Limdep (or NLOGIT) program. This response is in connection with the trivariate probit estimates with cluster robust SEs I posted earlier.

"Fred. The estimates are not that far apart. And, note that the log
likelihood values are essentially the same. I can tell you that nlogit is using the
GHK simulator to do the estimation. I do not know what Stata is doing. You can be sure
that the estimated coefficients are very sensitive to the values of the correlation
coefficients. But, in order to compare the results, you should probably compare the
estimated partial effects. I don't know if Stata knows how to do that. The nlogit
command is described in the manual."

/B. Greene

[Richard Williams]

Perhaps Greene didn't understand the main problem. I think we've already agreed that the coefficient estimates and LLs are quite similar. What is radically different are the standard errors when clustering is used. I think Limdep is doing it wrong. Normally clustering makes the standard errors larger. In Limdep they become smaller, indeed, dramatically so, and hence the Z values become much bigger. It would be nice to get a third opinion (e.g. from SAS, SPSS, or R) but in the absence of that I would trust Stata over Limdep. At least when clustering is used. I almost wonder if Limdep is doing the opposite of what it should do, e.g. dividing when it should be multiplying.

[Fred Dzanku]

Dear All, Just to bring some closure to the above issue I can report that there was some problem with the cluster estimator in the multivariate probit routine in LIMDEP. I post below by conversation with Greene. Those using later versions of LIMDEP may not have this problem.

"Dear Fred. I may have an explanation for your finding. Are you using NLOGIT or LIMDEP?
/Bill Greene"

"Fred. I believe I have - there is indeed a problem specifically with the
cluster estimator in the multivariate probit routine. However, you
are using an old version of nlogit, and I cannot do repairs on it, so at least
for this computation (this is hard for me to say), I think you might as well
use Stata. After a patch, I'm confident they will give essentially the same
answer. But, as noted, I cannot repair your program.
Regards,
Bill"

[Richard Williams]

Perhaps Stata Corp should quote him in its promotional materials ;-)

I'm not sure if I have ever seen a maintenance release for Limdep, at least in the last few years. Maybe I just don't know how to find them. I will be on the lookout for one.